{ "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.905229355379 NaN 1.95e-01 4.1 \n", " 2 -7.909828405416 -4.60e-03 2.98e-02 1.0 \n", " 3 -7.910018524177 -1.90e-04 2.99e-03 1.1 \n", " 4 -7.910051933691 -3.34e-05 1.35e-03 2.4 \n", " 5 -7.910052567947 -6.34e-07 7.95e-04 1.0 \n", " 6 -7.910052843377 -2.75e-07 2.43e-05 1.0 \n", " 7 -7.910052853975 -1.06e-08 1.17e-05 2.7 \n", " 8 -7.910052854039 -6.44e-11 3.73e-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.0795109 \n AtomicLocal -2.1806274\n AtomicNonlocal 1.7339358 \n Ewald -8.3979253\n PspCorrection -0.2946254\n Hartree 0.5417537 \n Xc -2.3920752\n\n total -7.910052854039\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.0883273 … -0.0562598 -0.114939 -0.0700162\n 0.201344 0.0909044 0.0122929 0.0111153 0.0420624 0.0176508\n 0.249296 0.174774 0.176138 0.132965 0.220115 0.11233\n 0.249296 0.231428 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.460337\n 0.369971 0.401677 0.412475 0.565555 0.388191 0.462728" }, "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+gvaeTAAAgAElEQVR4nOzdd2AUdfo/8Oczu+mF9GyyCQQChh6QlgQISEJyeCKIIFKEU8926Cn6RSxnwwZ6KAqeh/hT1DsLYD+FQEAgkFBClxoglPQe0rOZ+fz+WIwRCCRhd6e9X3+R3WV5Mswzz8ynMs45AQAA6JUgdwAAAAByQiEEAABdQyEEAABdQyEEAABdQyEEAABdQyEEAABdQyEEAABdQyEEAABdQyEEAABdQyEEAABds2UhPH78eH19fds/L4qiDf91wPG0Ic45Vh+0IUmS5A5BU5DstmXLQnjHHXccP3687Z+vra214b8OOJ421NTU1NjYKHcU2oGT07ZwPG0LTaMAAKBrKIQAAKBrKIQAAKBrKIQAAKBrKIQAAKBrKIQAAKBrKIQAAKBrKIQAKnbmzJnS0lK5owBQNxRCALX67sf/RY//yw0xCadPn5Y7FgAVQyEEUBlOlFHEH80QZ3+XUx3Yu9I1eMyq4hf3iscrsSYcQEcY5Q4AANrqcDlfnS19foobGd3RjaUvuifjxy9MQbd0jRu6Olsat050M9CUbuyu7kKkN5M7WADVQCEEULqz1fy7M/zTk1JxHd0WwT4aaRhhstY5Q5+7Z1s/08fX8PxASi/kq7OluB+bwj3YXd2FqZGCyU3GwAHUAYUQQKFya/iabL46WzpWwW/vKrwTYxhuYld50BMYjTCxESbDWzGGX/L4p1nSi3stfXzZlK7CjO5CgKvjIgdQF+O3337bv3//yMhI68+NjY0ZGRkXLlwYMGBAeHg4EeXm5hYUFDT/hYEDBwoCehYB7KW8gX48J63OljIK+c3hwvxo4U9hglN7cs7AKNHMEs2GetGwIVdafZov2GeJC2ZTugqTIgRPJ7uFDqBOxlWrVj344IOvvPLKfffdV1BQMGLEiM6dO4eFhd19992vvPLKgw8++O6773766adms9n6F7Zu3eru7i5v0ADaU9dE/zsvfZolbc3n8SFsSlfhqzGC+/U12bgaaHxnYXxnqmw0fH9WWp0tPZIuxoewWT2ECV0EZ9zQAhARkfGLL77Yvn37hAkTZs6c+emnn4aHh2/atImIvv3224cffvjBBx8kor/85S+vv/663KECaFCDSOtzpdWn+Q/npCGB7K7uwuc3CV62fmjr5EyzegizeghlDfS/c9IHx6QHt4l/DhemdGPjwgQjKiLom5GIhg8f7uHhkZaW5uHhwX7rgxAEwdPT0/rnsrKy9PR065OibJECaIjIKaOQf3ZSWp0t9fZhU7oKi2OcAu3fjefncrEi5tTwr7P5ogPSQ9ul2yPYlK7C1TsgATTsYstLeHh4Tk7Ovffem5mZmZSUZDabf/3115UrV1rf3bp167Fjx/bv3z9u3LjPPvvMyenK96sWi+XHH3/cv39/yxcHDx7cs2fPK35eFEVRFG32q+gejqcNiaIoSZINj2dTU9O/Pvh/Li7Og26d/d/TbHU293NhMyPp4EQhxJ2IOJFD//dCXOnhXvRwL3amiq/K5vdtE2ubaEJnmmyq3rzqw149uk2aOMGG/xxOTtvC8Ww7QRDYte7xLhZCJyenhoaGkydPbtu2bfr06Waz+ejRo2vXro2NjX3xxRcXLVpERCUlJbGxscuXL3/44Yev+F0Wi2X79u0+Pj4tX/Tz8+vatesVP9/Y2NjQ0NDuXwtageNpQxaLRZIkG37h6tWrn/p2v6W+pvOZoHsmJW8cK0ZcbHAhef/TQpzp0Sh6NIoOlLNVZ9gtT71T3ci9Pl96Q/fIHj162OpfaWhoaO0GGjoAyd52rq6ubS2ExcXFoaGhr7322sSJE1966SUiuvnmm7t06TJnzpzg4GDrZwICAiZNmrRr167Wvsvd3X3hwoXR0dFtjE8URYy7sSEcTxuyFkIXFxdbfaG5e++m7KU+TtLPM7v37m2zr7WhWHeKNdOQ3D4PLHi3rvqCOSzchqeTJEk4OW0IyW5bAhHl5+efPHly6NChTU1NRuPF0mg0GiVJuuSm+MiRIyEhITKECaBy39CN9y1fdyZzc+/eveWO5Wqm3zE5K3XVsLe2fVfic+1PA2iCceXKlStWrLjrrrtCQkJmzZo1ffp0Pz+/kJCQ5cuXJycnh4SETJ48ediwYd7e3ps3b05PT1+2bJncMQOozOFy/vUZ6cjkAG8lPgpeymQyvTuK37yuaXJXwcdZ7mgA7E/YsWPH7Nmz33//fSK65ZZbUlNTy8vLd+3aNXv27O+//56I7rzzzqKiooMHDw4aNOj48eNdunSRO2YAlXk4XVwwyOCvhipoNdCf3dpFWLAXwzFAFxjnNluxPjo6+tNPP217H2FVVZWXl5et/nXA8bQhG/YRfnVaWnhAypxoNKhqdkJZA/VaY9l4s7Gvrw3irq6ubp6OBdcPyW5bmEkLYEd1TfTUbmlJjEFdVZCI/FzomWjDYxl4KATtQyEEsKOFB8TYIDYqRG1lkIiI5vQWiurph7O2nEYCoEAohAD2cr6G/+uo9PoQtWaZUaAlMYa5O6R6PBaCpqk1RQGUb+4O6dE+hi6eqnwctBoTyvr7sbd/xUMhaBkKIYBd/JLPdxfzx/upPsUWxwj/PCieq7bZqDoApVF9lgIokMjpsQxxScz17qOkBN282EO9hGcz8VAImiVbISwpKfn4k09PnjwpVwAA9vP+USnAlW6L0MiN5jMDDGkFPK0AD4WgTbIl6m133fd4hmXYzZPRDw8aU9pAL+8T34k1yB2Izbgb6bUhwtwdooRSCFokWyEMCgxwL/i1XnAP/dwya7P4v3O8EU0voAnP7xGndhNsMg9dOaZFCm5G+ugEshQ0SLYejNUrl2/dunXYsGG1gtNP56Tlx8S7t/JxYcKUbiw5THDWSJMS6M7hcr4mWzoyWWtbDjGid2MNWIAUNEm2giMIwqBBg9zc3PxdaFYP4cck49HJTiNMbNEBKeS/llmbxR/PSRbcfYLaqG5Z0bYb6M/GYwFS0CIFPXkFuNL9PYVt440HJhkHBfyhIjahIoIarDotVTTSX6MUlFa29foQw39PSb+Wo6sQNEWJGRvmwR7tK2wbb9zXXBE/R0UEpatrovnqXFa07fxd6Jlow9wdeCgETVFiIWwW/ltFzJz4e0V8YJu4rcB2W2YA2MjCA2KMapcVbbs5vYXCOixACpqi6ELYrIvnxYq4e4Kxtw97arfY+YumRzNQEUEprMuKLlTtsqJtZxTobSxACtqisryN8LpYEbfeYgh1Z/dvE7v8sSJKEm5UQQaPq39Z0bZLwAKkIB97XORVVgibdfVi86OFI5ONPyUbvJ3p3jSxx6qmYfcvCOoz7J5H/k/u6EBfNufzPSX8CfUvK9p2i2OEtw6JebVokQHHkSRpePIEc9+h//1qjW2/WfWp28+PvTzIcHyKcXWC4cTWn0pnrPzf+o1yBwU6InL6e4b45lDBTf3LirZdNy/2QE9h/i48FILjlJeXHzxXXJD0wuof19n2m1VfCJsN9Ger/7Wod9rLXncsEHGfCo7y/lEp0JVu76qdVGojLEAKDsY8/SlmetyZVYuem2fbb9ZU9iaOGf1rypedY8atOIYbVXCEsgZ6ZZ+4JEY7y4q2nXUB0kczsAApOMhTu8W/PvDQ9u8+i4qKsu03a6oQEhEjWhZneGGvWFwvdyigA//IFKd2E/r56WKMzOWmRQpeTvQxFiAF+8ss4f87J71wo11uOrVWCImojy+bHik8sxuDu8G+Dpfzr89Iz9snM1WBES2JNfwjU6xolDsU0DSJ05zt4qKhBjutc6vBQkhECwYZ1uXwHUVosgE70vCyom2HBUjBAT44JjkbaGZ3exUsbRZCLyd6fYgwJ13EqBmwE80vK9p2WIAU7KqsgV7cKy6Ls+PahZpN4xndBW8nwqgZsAc9LCvadliAFOxq/i5xWqQQbc+eeM0WQoyaAfvRybKibYcFSMFOMkv4T+ftNUammWYLIWHUDNiHfpYVbTssQAr2YB0j84bdxsg003gyY9QM2NzjO6S/62ZZ0bbDAqRgc9YxMjPsNkammcYLIUbNgG39ks93FetrWdG2sy5AmluDZAMbKGugl+w8RqaZ9vN5ZnfB24k+wKgZuG4ip8cyxLdjBHc9LSvadt282P09had2I9fABubvEu+08xiZZtovhES0LM7wIkbNwHV7/6gU4EqTInSRNR3zLBYgBVtwzBiZZrpI6T6+bAZGzcD10fOyom2HBUjh+jlsjEwzXRRCInppkGFdDs/AqBnoqH9kinfoeFnRtsMCpHCdljtqjEwzvRRC66iZhzFqBjrEuqyowxpqVA0LkML1cOQYmWZ6KYRENLO70AmjZqBDsKxou1gXIH15HzojoN3m7xKnO2qMTDMdFUIiWopRM9B+WFa0A14fYvjPSSxACu3j4DEyzfSV2xg1A+2FZUU7xt+FnsYCpNAezWNkOjlqjEwzfRVCwqgZaKdFB7GsaAc9jAVIoT0cP0amme4KIUbNQNudr+HvHcGyoh2EBUih7WQZI9NMjxmOUTPQRlhW9DolhLJ+fmwJFiCFa5FljEwzPRZCwqgZaIP0IoZlRa/fWzHC4kNiXq3ccYCCyTVGpplOk9w6auZpjJqBK2lqapr/4qtTHn1+Yf86LCt6nbp5sTt9c5P+Ov+9Dz6SOxZQIhnHyDTTaSEkohexQxO0IiUl5V878ssahboda+SORQuKvlt8NjjmH0tX5uXlyR0LKM4HxyRXmcbINNNvIfR2okXYoQmupF+/fsLZPT5Z6wcPHiR3LFowPiHeffNSF6MhICBA7lhAWUob6MW94rsyjZFpput2nxndhf93XPrgmPRQL/3eEMDlzOGdXZ/blvlnsVuAh9yxaMHMO6c490t8L9vL2VnXFxy4nLxjZJoJf/rTn959911JkojIYrG8/fbb48ePv/XWW99///3mF1977bXk5OTZs2dnZWXJG67NYdQMXG5XMQ/3NJi9cNW2mdFmp32lvNoidxygJJklfO15roQlfIUnn3xy+fLl77zzDhG99NJLn3zyyRNPPPHII48sWbJkyZIlRPTcc8/98MMPzzzzTERExE033VRfr6migVEzcLnUXJ4QKncQ2uJupEEBbCv2KYTf/DZGRpBxjEwzYcyYMW+++eY777zDOd++ffv9998/evTosWPHzp49Oy0tra6ubvny5cuWLRs1atRLL70UHBy8Zo3Whg+8NMiQgrVmoIWNedKYELmD0JxEs7AxDxMK4SLrOjLTZR0j00wgori4uLNnz5aWlo4dO/aHH34oLS0tLCz8+eefk5OTT58+XVdXN2jQxSEDsbGxe/bskTVg2/NyooVYawZ+U9tEe0v4iGC549CcxFCWmoscAyKiUlnXkbmckYg6depkMBgKCgrmzp2bmpoaEhLCOb/55pvvvffetLQ0X19fxi5G6+/vf+LEida+q7q6evr06W5ubi1ffOSRR26//fYrfr6mpqb5m+U1wUQfHHF6d7/lvh4qbiNVzvFUtfV5wgBfA6+vrZYkiwWdWrZRU1PT053OVzufLK4xuaEcXi+1J/vju4xTOlOkc0N1td3/LXd3d0G4xnOnkYgaGhpEUfTy8po1a1ZERMTatWtFUZw9e/acOXP++te/1tXVNX+6trbW29u7te9yc3N76qmnevTo0fLFiIgIT0/PK36ec97aW473r5F8zM9NM3u5BbrKHUpHKep4qtf2cjEpnHl6ekqS5OKCHQhtxtPTc3SouLPSOCNQEa1hqqbqZM8s4akF4pHJRk8F9A5aGYno9OnTzs7OJpNp7dq1GzZssCb/Pffc88ADDyxYsKCqqqqkpMQ6Aej06dM33nhja99lMBj69esXHR3tsOhtqI8vm9ldeHq3+OFI+YcwgYxSc/nyEbhS20VCKNuYx2d0lzsOkI+ixsg0E4hoxYoVt912m4uLS2Rk5JYtW6xvbNmyJTIy0mQyxcfHr1ixgojOnDmTmpp6xx13yBmvPb14I0bN6F1RHZ2r5oMDVNzopGSJZrYe3YT69u+jChoj08w4YMCAurq6tWvXEtF77703bdq01atXNzU11dTUWAeIvv3227feeuuaNWvOnj375JNPXtLyqSXWUTNztou7JxqV0ocLjpWaJ40OEYwCWVTcWaxcUZ2YkdHxSh7VCQmmR6UNtGCfmDLOqLT/fuPnn38eFRVlMBiIaMSIEadPnz579qzBYAgPDzcajUQ0YMCAU6dOHTt2LDQ01N/fX+6A7cu61szyo9LfeivrhgUcY2MuTzQrLUk1ZUwoS81FIdQphawjczmhd+/e1ipo5eTk1L17965du1qrYPOL/fr103wVtFoaZ3hpn1hUd+1PgvZsykchtK9EMyZR6NTuYqWsI3M5PPdcqnnUjNyBgKMdr+QWifCwYleJocLmfKkJE+t1RuL0cLrixsg0QyG8ghdvNKzPxagZ3UnN5Ul4HLSzIDfq7MkyS5Bc+qLMMTLNFBqWvJpHzWCtGV3ZmMcTQlEI7Q6to3pjHSOjnHVkLodCeGUzugs+zrT8KFpw9ELktCVfGhOKjLC7hFAsOqov83eJM7orcYxMM6R9q5YNx6gZHcks5mYPFuIudxw6MCqE7SnhNU1yxwEOYR0j8/xAJY6RaYZC2KrePhg1oyOpeXwsOggdwsNINwawNGzJpAMKHyPTDIXwal680fDzvtOPvbYsJydH7ljAvlJzpQS0izpKQqiQmovWUY07ePDgnc//S6gtVewYmWZKj09eXk4kfTDz3ZNu4+64S+5YwI5qmyizhMeb8EToIGPNbAPGy2haU1NT4qQZq0/WG1c9ofy8QiG8hghTADu5zdsvSO5AwI7SCvhAf+bpJHccujE4gJ2r5oXogNcug8HAnF2dz+7qHWGWO5ZrM177I/q2PeX7Oz49OHpIf7kDATtKzZUSzbgpdByjQKNChE150rRIHHZtYowlvbM5oibr5YkD5I7l2nAWXoPRaJx+04D1eXLHAfaEkTKOl2hmG/PQOqpZnGhTscvdCerYlQ+F8NoSzUJaAW/A6FGNKqmn7CpsveRoiaHoJtSyQ2Xcw4m6eakjrVAIr83HmXr7snSsuKZRG3Kl0SGCE1LBsXr6MM7pRCXSSpvWq2rBQmR/mySHsZQcjPbWJqysJpcErLWmXSk5UnKYatIKhbBNksOElBxkrDZtzMPWS/JICEU3oTbVNdGuIj4qRDX1RTWBymtoIDtXzfNr5Y4DbC2rkjeK1MsHhVAGY83CL/kSlrbXns35/MYA5q2e+UgohG1iYHQTVgrWIowXlVGwG5nd2R5syaQ563OlpDA1FRc1xSqvZDND66j2pOaiXVROY9FNqEUpOTxZVWmFQthW48JZSo4kIWc1ROS0OR9LjMopwYxFR7Ump4aX1PMB/iiEWhTmwfxd2YEyVELt2FPCQ92x9ZKcRoew3diSSVvW5fCkMEFQUx1EIWyP5DC0jmoK2kVl52Gkgf5sG7Zk0pD1OWqaQWiFQtgOSWZhPWYTasjGPLSLyg8b1muJyGlTnjRWbSv3qixceY0OYZklvNoidxxgC/Ui7S7m8SEqu3XVHoyX0ZLdxTzMQ33dDSiE7eBupKGBbHM+klYLthbwAf5qmuqkVUMDWXYVL66XOw6whZQcrqIFZZqhELZPcpiAtda0AVvSK4RRoPgQYRNaRzUhJUdKVtUMQiv1RSyvpDC2Hs04moCRMsqREIrWUS2oaKTD5Xx4sPrSCoWwffr7sSoLP12FpFW30gY6XcWHBqovYzUp0YwtmbRgY640wsRcDHLH0X4ohO3DiMaaBSSt2qXmSvEmbL2kFL19mMjp5AWklbql5HI1tosSCmEHYDahBmDHCaUZg9ZR9UvNVeVIGUIh7IAks/BLnmRB176aoYNQaRLN2JJJ3Y5VcIlTVCdVphUKYbsFuFKkN9uJDetV69QF3oCtlxRmrFnYmIctmVRMpRMnrFAIOyLJzNZjpWDV2pDLE81MrSmrUSY3CnVne7Elk2qtz5VUt7JaMxTCjsCG9aq2MY8nhKo1YzUs0cxS0TqqTo0SbS/kN6l2Yq5a45ZXXDA7UYm1MFRJ4rQ5XxqDQqg8CaFsIxpa1CmtgPfxZX4ucsfRUSiEHeGEtTBUa28pD3JlYR4ohIozOkTYVczrsCWTCq1X54IyzVQcurySzGw9WkdVKDWXj1Vtl762eTlRtD/bVoi0Up8UFW691BIKYQclh7GUXI6UVZ3UXAkdhIqVEMqwYb3qFNTR+Ro+RM3rNKEQdlB3b+ZqoCPlKIVqUi/SrmI+KgSnvUIlhgoYL6M663OkhFDBoOI6iEJ4HZLMWGJGZdIKeH8/bL2kXDFB7PQFDENTGVXPILRCIey45DDMJlSZjXkSFpRRMqNAI0zsFwxDUw9OlJonjVV5WqEQdtyYUCGjEIPc1CQ1lyeqdqqTTiSGClhrTUX2lfAAF9bZE4VQr6yD3NIwyE0lShvo5AU+LEjdGat5iWZs+akmKbk8SeXtooRCeJ2wYb2KbMqTRpoYtl5SuD6+zCIRtvxUC5VuSX8J1f8C8sKWTCqSmssT0C6qBmNCsE+vOtQ00b4SHm9S/xPhsWPHKioq5A5DrW70Z8X1/HwNklYFsPWSWiSY2UYUQjXYmCsNC2LuRrnjuG7C1KlTu3Xr9vbbbxPRF1984fcbb29vxtiZM2fmz5/v5ubW/HpNTY3cMSuIwCghVMCGosp3porXNPE+viiEKjDWzDblSRKySvHW5/Ik9beLEpFw4MCBjIyM55577vz589OmTSv7zYIFC2JiYiIiIojosccea37dw8ND7piVBbMJVWF9Lh9rFlAGVSHUnQW5sX2lSCulS8nhyZpoZRGIKCoqKj4+ftWqVS3fWLly5T333GP9M+e8vLxchujU4E/hwoZcbCiqdBvz0C6qJmPNDA0tCpddxastvK+fFtLq4lNtZGTkuXPnml/dvXv3yZMnp06dav3x3XffveGGG7y9vV966aWrfJcoiocPH979R0VFRfaLXglMbmT2YHuwoaiCSZx+ycMSo2qSEMpSMa1e2dbl8D+Fa6SV5WIvp7u7e8uK9dFHH02ZMsXb25uI5syZ8/LLLzs7Ox86dCghIaFHjx7Tp0+/4nfV1dW9/PLL7u7uLV985JFHJk+efMXP19TUME3sEz4myPjj6abebjJPrdfM8bS5fWXMz9nJh9dUV7f1r1gsFkmSLBaLPePSkfaOLRjaie0odCquqHZT/0AMe1BCsv98xmlSZ7G6Wun3K+7u7oJwjY7Mi2dZWVmZyWSy/rmuru7LL7/84YcfrD927tzZ+od+/frdddddqamprRVCT0/PTz/9NDo6uo3xcc49PT3b+GEl+3NXvmCf+LKnq7xhaOZ42lz6KSkpvH0Hx1oIXVxUu9Oo8rTr+HsS9fdvOlDjgQbtK5I92Zsk2l5i+X+jXTzdZIzCZi7WyYyMjIEDB1r/vGbNmsDAwBEjRlz+6eLiYutjIrQ00sQOlvLKRrnjgFZsRLuoCiWGso1oHVWqjCIe6cWCNFEFiciYkZGxbt26qqqqiRMnWl/66KOP7rnnnubn7qeeeio2NtbPz2/jxo2rV69OT0+XL1qFcjVQXDDblCfdFqGFkcQaUy/SjiK+OgH/NSqTYBae2CG+PkTuOOBKUnIkte840ZLxmWeeCQ8P37hxo7OzMxHV1tZGRUXNnj27+RO+vr4rVqyoqamJjIxMT09vfnCElpLDhPW5/LYIueOAy2wv5H19WSdnueOAdooNYlkXeEk9Bcjc5wBXsD6XLx5mkDsKmzH+8ssvLX92d3f/97//3fKV+fPnz58/37FRqU9SGHvnMJpxlGhjrur3iNEnJ4FGBLPN+dLkrniaV5bSBsqq5DEaWr8eZ5ht9PZhEqesSkyiUJwNuTzBjPNclRLMWLZJiTbkSKNCBC2tX6+hX0VuY80sBUmrMOUNdKKSDwvUzq2rriSGslTsTag8Kbmq35L+EiiENpMcxrAlk9JszJNGmJiLdvoy9KWvH6tt4tnYkklhUnN5kra6G1AIbSYxVNiazxtEueOAFjbmYeslFWNEY7CovcIcKuOuBor0RiGEK/F1od6+LKMISasgqbkcI2VULTGUbUTrqJKs11y7KKEQ2laSGa2jCnK2mldpZVFg3RprZhuxJZOSpORIGmsXJRRC20oOE7Alk3JsyOWJ2HpJ5cweLMCVHShDWilCXRPtLOKjNdfdoLXfR17DgtjZap5fK3ccQEREqbkcK6tpQGIo24BuQmXYUsAHBjBvJ7njsDUUQlsyMLopVMACiUrAiTblSViyWQMSzGxjLnJKEdbnSMma2JL+Ehr8leSVZGbr0TqqAPtLuZ8LC/dAIVS9MaFCRhGvk3mXMyAiSsnR2sQJKxRCG0sOY+tz0bUvv9RcbEmvEd5O1AfjsRUgp4YX1/OB/hpMKxRCG+viyXyc2YFSJK3MNqJdVEOsY0fljkLvUnL4WLOgyeFnKIS2lxzGMHZUXo0SZRTyUSac3hqRECpgvIzsUnI0OIPQClcK20sOE9ajb19W2wp4b1/mi+3ltSI2mB2v4GUNcsehYyKnTXnSWI2uX6/N30peo0PY7mJebZE7Dh1Du6jGOAs03MR+QeuofDKLeZgHC3GXOw77QCG0PXcjDQlkWwrQkiOb1FwsMao1CaEC1lqTUUouT9JouyihENpJcpiAtdbkUtFIxyp4rIZ2DQUiGmvGlkxyStHoDEIrzf5i8koOw2xC2WzKk4Zj6yXN6efHqhr5GWzJJIcLFjpczkcEa/bmEoXQLvr7sQsW7KMmD7SLapJ1Sya0jsoiNVcaHqzlm0tcL+yCEY01Y8C3PFLzMJVemxLQOiqTlByepNHxolZa/t3klWTGbEIZnKvmlY28P7Ze0qIkM0vNxZZMMtigxT0IW0IhtJfkMGFTnmTBiBnH2pDLE0I1ufYFUJgH83NhB7Elk2Mdr+QWiXr6aDmrUAjtJTXJoIMAACAASURBVMCVunmzXcVIWofaiHZRTUtE66jDpeTwceEazykUQjtKNrP1mEThQJzolzxpTIjGk1bPEkOxJZOjrc+RkrV+c4lCaEdJYUIKxss40MEy7u3MIrw0nrR6dlOokF7IG0S549CNRom2FfKbtD4MW+O/nryGB7MTlbwUCyQ6ygZsvaR1Ps7Uy5elY0smR9lWwHv7MD+tL9uLQmhHTgKNCBbQkuMwG3OlhFAUQo1D66gjaXtBmWba/w3llRzG1qN11CEaJUov5KNDcEprXIJZwHgZh1mv9YkTVrhq2Bf2JnSY9ELeUwdtODA8mB0t5+XocbC/wjo6V82HBKIQwvXp7s2cBTpSgVpodxvzpLHoINQBZ4HigtnmfLSO2t36XCkhVDDoIKtQCO0uCQ+FDpGayxM0vQoUNEswY9FRR0jJ0fLWSy3hwmF3yWEMWzLZW0UjHS7ncdh6SR8SQxkW8rU3TrQxV0rSRysLCqHdJYQKGYW8rknuODTtlzwpTtOr40NL0f6sElsy2dm+Eu7rwjp7ohCCLXg5UbQ/21aIpLWjjXnYeklHGNFNocKmfOSUHelkvKgVrh2OkGTGhvX2lZrLMVJGVxJD2Ua0jtpTSo6k7a2XWtLL7ykvTKKwH4vFMuX+x06/NcO/Lk/uWMBxggr3rH5i4stvLpE7EG2qaaI9JTxeN8v2ohA6wqAAVljHc2tQC21v586da7MqLTeM+fg/n8sdCzjOkrfestz8zFsf/qe2tlbuWDTolzw+LIh5GOWOw1FQCB1BYJSADevto3///q4Xcv33fDJhXJLcsYDj/HXGFPev/h5+Q193d3e5Y9Egnays1kw3FV9uSWaWksv/coPccWiOt7d38NM//+cmQ7S/XppxgIimTZlU03/iVoyXsY+UXL4mQUeFUEe/qrz+FCak5koS0tbWKhrpfA3v54sqqDvDg7ENhV2cqeJVjbyfn45yCoXQQULcyeTG9pQgb20svZAPC2JGnMj609OHVTTwfHQR2tq6HJ4cJuioDKIQOhLGjtpDeqGEBWX0iRHFBLGMIkxMsjFdzSC0QiF0nOQwYT32UbO19EIeG4zTWKdig4UMLFVhU00Sbc6X9LY8hb5+W3mNNLEDpbyyUe44NKRJoj0lfJgOtomBK0I3oc3tKOKRXizITe44HAuF0HFcDRQbzH7B9jG2c6CMd/FivtiDUK+GBbKDZbxelDsODUnJlXSy40RLxjlz5gwdOvSuu+4SBOHnn38+dOhQ83uurq6PPvooEX3zzTebNm0KCwt74IEHfH195YtW9ZLDhJQcPrGL3HFoxfZC7Diha25G6tmJ7Snhw4NxGthGSg7/5zDdrV4vDBw48J133nnyySeJqKampvw3a9as+eabb4hoyZIl8+bN69+//8GDB0ePHi2KuPvquOQwtg7jZWwno4jH4Qqob3HBLB3dhDZS3kBZlTxWhzeXnPMTJ064u7uXl5fzFnr16vXJJ59YLJbQ0NDU1FTOuSiKPXr0+P7773kr+vfvv3///tbevdyFCxfa/mHN6PyFJatSssc36/B42u9gNjY21tfX2+Ob9amqqspO3/zlKfG2DU12+nLFslOyf3lKnLBedweTcy4QUY8ePQICAvbs2dNcHbdt25abmzt58uRTp04VFxePGjWKiARBuOmmm9LS0mQr2pow1oxJFLaRW8PrRd7dW393r9DCiGC2vRD97raRksN1shPvJS4usRYUFJSfn9/86kcffXTnnXe6u7vn5+f7+voajb9/LDs7u7Xvqq2tfeaZZ3x8fFq+OG3atDFjxlzx83V1dQaD7hqjRwWwL88Y7o6w/Ua9ejueG8+xYf4GO625bLFYJElCR4Ct1NbWCoJdhub5MnJmTocKayO97PH1CmWnZN+Q4/R4lEVjy5i7urpe89y7WOEsFouzs7P1z9XV1WvWrNmwYQMROTs7NzX9fsm2WCwuLq0O0XNycho2bFh4eHjLFyMiIlr7K42NjVf5Nq0a14Ue3iWR0cXm26nr7XjuqeDDTWSnX1kQBEmSdHU87erql47rFBfM91Q49w7Q0aOMPZL9cAU5GXivAK2d84xd+8QwEhHnPC8vz2w2W1/68ssvO3fuPGzYMCIym80VFRXV1dWenp5ElJOT07Vr19a+y8nJacKECdHR0W2Mz2Aw6OoJxirAnXr58l2lwmhb7/Wlt+OZXti0JNZgMNjl8idJEmNMV8fTrux6cg43STuK+V+idPSfZY/juSFPGhdO+jznBSLatGmT0WgcOnSo9aWPPvronnvusf65S5cu/fv3X716NRFVVFSkpKRMmDBBrlg1Y1DjsVcWLy0sLJQ7EBWra6JjlXyQnh4CoDVxmFZ/3Y4dO7biX8uGOuv0oiTMnDlz+vTpCxcudHJyIqLjx4/v2bNn5syZzZ9YuHDhvHnzpk+fHhsbO27cuMGDB8sXrUasenrGpmLnSbPulzsQFdtVzPv6Mlc93rzCpQb4sTNVWLPpuiTcNv14nevSJ++TOxB5GCdOnLhgwYJu3bpZf/b399+3b19QUFDzJ5KTkw8cOLBz585HHnkkNjZWpjg1JdC3U+nJ7cE9QuUORMW2F2IONVxkFGhQANtRpLulom3IxcPLeGpbl646vSgZJ0+e3PLngICAgICASz5kNpsnTZrkwKg0bv+21FErfr37lgFyB6Ji6YXSPVFYIBAuGh7M0gul5DA0EXTQzPdSck8d/XCaTi9KuJTIwNnZ+ea4AduwHEZHcaKdxTw2CGcvXBQbLKCb8HqklxpvHzmgLQMsNQmXEnnEm9jWAuRtBx2r4N5OLMRd7jhAMeKC2M4i3oSJ9R1ikWhXsa5XK0QhlEdMEPu1nFdb5I5DndILdZ20cDlfFwrzYL+W4+ayIzJLeI9OzMdZ7jjkg0IoDxcDDfRnO9CY0yEohHC54Vh9u6O25vN4k64TCoVQNvEmllaAppyOSMemE3CZWMwm7Ki0AmkkCiHIYqRJQDdhB5Q2UF4N7+ur67yFy8UF4YmwIyRO6UV8RLCua4Guf3l5DQ9mmSXYXLvdMgr5sCBmn4XVQMWifFi1hefWoBa2z4EyHuLGgtzkjkNWKISy8XSiXj4ssxh52z7phRLaReFyjCgmSEC/e3ttzefxtl73WHVQCOWESRQdkF7E4/TdjAOtQTdhB6QVcp13EBIKobxGYrxMO1kk2lvChwXqPW/hijBwtL04RsoQEQqhvEaahPRCzAJuh/2lvJsX66TjCU9wFUMC2KEyXmf7Ta8161gF9zSycA8UQpCPnwt18WL7y3AP21aYQQhX4WakPr5sTwkSqq22FqCDkAiFUHbxJrY1H3nbVulFPBaFEFoXF8y2o3W0zdIK0EFIhEIou5EmlobxMm2Wjt2X4KriglkGxsu0WVqB3teUsUIhlNmoEGFrgSQhc9vgXDW3SLybF/IWWjU8mG0vRD61SXYVt0i8uzcSCoVQbiY38ndlRyqQude2vZAPx8QJuKpQd+ZhZFmVSKhr21rAR4UgoYhQCJUA3YRtlIEOQmiDOEyiaJu0Aj4SCUVEKIRKMNLE0pC3bYAOQmgLdBO2EYaMNkMhlF+8iW3Jx1zCa6hpomMV/EZ/5C1cQ1wQBo5eW0Edldbz3j5IKCIUQiXo6sWcBHbyAlL3anYW8QH+zMUgdxygeNH+7Hw1L2uQOw5l25IvxZsEAXWQiFAIFWIkFh29FkylhzYyMBocyHaidfSqMIOwJRRCRcBswmvKKMKmE9BWccEsowjdDVeDDsKWUAgVAQNHr44T7SjiMUE4XaFN4oIEdBNeRVkDna3iA/xQCC/ClUURevqw6iZ+HnuKtuJIOfdzYSZ97x0KbRcbzHYXYzn7VqUVSHHBzIjL/29wJBSBEY00CWgdbQ06CKFdfJypsyc7iOXsW5FWwEeacPH/HY6FUowMRjdhqzKKUAihfYZjk97WbcUSo3+EQqgU8SHoJmzV9kIeF4S8hXaIDcL6MldWbaGjFXwwdrduAYVQKaL9WH4dL6qTOw7lKamnwjre2xd5C+2AhdZas72QDw5grpiS2wIKoVIIjOKC2LZC9O9fKr1QigliBtRBaI8enVidyHMwAO0yaQUS2kUvgUKoIBgvc0UZRTwOm05AOzGi2CABi45ebitGylwGh0NB4kOwvswVoIMQOgbdhJdrEGl/KXZxuRQKoYIMDmBZlbyiUe44lMQi0b5SPhSFENoP3YSX21HE+/gyD6PccSgMCqGCOAk0NBCp+wd7S3h3b+btJHccoEKDA9jhcl7bJHccSoKJE1eEQqgsI01CWgHGy/wuvQh7EEIHuRmpnx/LLMGd5e/SCiR0EF4OR0RZ0E14ifRCHot2UeiouGDsTfi7Jol2FfMReCK8DAqhssQEsoNlaMz5HdaUgesRG8QyUAh/s6eEd/ViPs5yx6E8KITK4makaD+2A2O+iYjoTBUXJd7VC4UQOmiESUgvlJBOVuggbA0KoeLEm9hWdBMSEVF6ER+B/gy4DiY38nZmJypRComwGW/rcJVRHEyrb4YOQrh+6Ca0kjilF2KkzJXhoCjOCBPbXcwb8UxIlF6IIaNwveLQTUhERIfKeaAbC8amnleCQqg4Xk50QyeWWaz31K1poqwLfGAACiFcF0yrt9qajw7CVqEQKlG8CZMoaEcRH+jPnHGGwvXp58fyanlZg9xxyA0dhFeBy4wSjTQxTKvfjl3pwRYMjAYHMqy+va0Qm060CoVQiUaFCOmFXNR35mYUSlhrG2wiLohl6HuDs+OV3MXAOnsioa4MhVCJ/FzI7MEOlOq3EnKiXcU8FrsvgS3EBQs67yZEB+HVCe7u7jExMVlZWdafDxw4MHLkSGdnZ29v74ULFxLRokWLIluora2VNWC90Hk34eFyHuDKAl3ljgM0ITaYZZZwi46fCdFBeHVCdXV1cnLyPffcQ0RFRUVJSUnTp0+vrKzMy8u75ZZbiKisrCwpKWnDb1xdcXFyhJEmpufZhOggBBvydqIIL3agTL8JhTVlrk4QBGHevHm7d+/OyspasWLF4MGDH3roITc3N09Pz759+1o/5OPj0+03goDWKkcYHSKkFeh3aagMbMYLNhWn4016z9fwBpHf0AkJ1SqBiDw9PTt37pyVlXXw4MEuXbqMGzcuLCxswoQJ2dnZ1g99+OGHwcHBQ4cO/fLLL6/yXZzzmpqaqj+yWCyO+D00J8SdvJ3Z0Qqdpm461toGm4oL1u/A0c35PD4EDzBXc3GjYm9v77KysoKCgp9//nnt2rU33njj008/PXny5D179kyZMuWee+7x8/PbvHnz3Xff7e/vP3bs2Ct+V1VVVWJi4iWPjK+88spf//rXK36+pqaGMVzsWhXrb9xwxtK5u9jGz2vmeJY2sOI6p3BjfXW1bDFYLBZJknAbZys1NTXyBjDAiz2d71xdXSdvGLbSrmTfdN441JdXV7f1SqIx7u7u12zIvFgIKyoqAgMDAwIC/vznP48YMYKIXnjhBX9//7y8vMGDB1s/M2XKlM2bN3/99detFUJvb+/vvvsuOjq6jfFxzj09Pdv4YR0aEy6l5vJHPQ1t/LxmjmdqqRQbLHl7yfm7WAuhi4uLjDFojLwnZ39PkrilnHmEe2jhZrFdyZ5R0jQ32uCJuROtE4iovLz8/PnzUVFRvXv35vxi64H1D5fcdHDOtfHMoQrxJrZFl+NlMgp5XBBacsDGYoL0OImiuJ6K63lfX1y3r0aorKx87rnnRo8eHRERce+996akpKSlpdXU1Lz44osxMTEhISEffvjhiRMnioqKvvzyy08++WTSpElyx6wXkd5MIMqu0l3qYsgo2IM+Fx3dki+NMAkC8umqhKioqPPnz3/88cdEFBER8dlnnz366KNRUVHFxcWrV68moh07dvz5z3/u37//smXLPvnkk9baRcEeRuhvNmGDSPtL+dBAJC7YmD4LIWYQtoWxoKCg5c/jx48fP358y1c+/PBDx4YEv7POJpzdQ+44HGhvKY/yYZ5OcscBmjM4gB2r5NUW0tXZtbWALx+BjoZrwAFSNB2uL4M9CMFOXAzU349llugooSob6fQFPtAfCXUNKISK1tuXVTbyvFodpS52pQf70du0+rQCHhPEnHCZvxYcIUVjRMODBV2ttZZRJGGkDNhJbDBL19M2FGkF0kgTLvLXhmOkdLpadPR0FWfEumDCE9jH8GAhvYjrZ+nCrQU8PgTZdG0ohEqnq27C9EI+AiPcwG6C3cjPhR2v1EVC1TbR4XIMwG4TFEKlG+DPcmp4aYPccTgEOgjB3vTTTZheyAf6M9e2rkylayiESmdgFBPEthXoomMjHVPpwc5idTObMK1AQrtoG6EQqsBIky7Gy1RZ6FQVH4Ch3mBPw4NZuj62odhawDFSpo1wmFRAJ92EO4r4oADmjFMS7KmvLyuo5SX1csdhZw0i7SlBR0Nb4aqjAkMC2fEKXqX1HYHSsRkv2J/AaGgg26H1h8LdxbyXD/PS0xo61wOFUAWcBRoUoP2OjfRCKRYdhGB/ephNuLWAx2MAdpuhEKpDfAhL0/R4GYnTrmIeg92XwP6sswnljsK+0gokrLXddrjuqMNIk6DtbsJfy7nJnQW6yh0H6EBMENtTwhu1e2Mpcsoo4sMxUqbNcKTUITaI7S/ldU1yx2E36CAEh/FyokgvdqBUs3eW+0p5F0/m7yJ3HOqBQqgO7kbq68t2FWs2dTGDEBwpLpht126n+9Z8dBC2DwqhasSHaHkSRXoRCiE4Tlwwy9BuNyE2420vFELVGGkStDpeprCOyhp4VCekLjhIXDDbptHbSk60rVDCmr3tgkKoGiOC2Y4ibfbwby+U4oKYgMwFR+nmxTjxs9UarIWHy7mfCwt1Rzq1AwqhanRypu7ebJ8W99fOKORxwTgVwaFigwRNzs1FB2EH4OqjJlrtJkQHITieVrsJ0wrRQdhuKIRqMjJYg5v0Noh0sIwPDkDqgkNpdeAongg7AIVQTeJDhLQCSdRW8maW8J6dmCcWRQTHutGfHa/g1dpawvfkBS4wivBCIWwfFEI1CXSlEHf2a7mmKiFmEIIsXAw0wF9rc3O3FvDR2IOw/VAIVSbexLbmayp1M9BBCDIZrrlNejGDsGNQCFVmpElr3YQZhRIKIcgiNpilF2lqQhI6CDsGhVBlRoWwLQWSZirhyQvcILBwD6QuyGBEsLCjiGsmnXJreJWFR/kgm9oNhVBlwjyYh5GdqNRI7qYX8hF4HASZBLiSvws7WqGRbNpSwEeFYF2KjkAhVB8tdROigxDkpaVuQnQQdhgKofpoqZtwO4aMgqxig5lmNulFB2GHoRCqT7yJbdbEE+EFC52p4tF+SF2QjWaeCEvqKa+W90M2dQgKofr06MREThpYLzijkA8OYE44B0E+fXxZST0vrpc7juu2tUAaYWIG1MEOwUVIlUaatLDoaEYRJk6AzBjR0ECWUaj6SRRpBXykCdfzDsKBUyVtdBNuL+Sx2HQC5BYbLGigm3BrAToIOw6XIVXSwMBRkVNmMY8NQuqCzDTQTXjBQlmV/EasXN9RKISq1NePldTz/Fq547gOh8p4qDvzc5E7DtC9mCC2r1TdW15vK+DDgpgzLucdhSOnSoxohEnYpuaODay1DQrhYaQeKt/yOq1AQgfh9cCxUyu1dxOmF/FYFEJQBrXvTYgOwuuEQqhWau8mTC/kw1EIQRlUvVt9XRMdLONDA5FNHYdCqFYD/dmZal7WIHccHVJYRxca+Q2dkLqgCHFBKh4vk1HEB/gzd6PccagZCqFaGQWKCWLb1dlNuK1AigvG6sCgFBFejDE6U6XKWri1QEK76HVCIVSxkSZBpd2E6CAEpYkNUms3IabSXz8cPhWLV+36MuggBKVRaTehRaLdmI973VAIVWxoIDtczqsscsfRTg0i/VrOB2PyLyiJSrsJdxfzqE6sk7PccagcCqGKuRhoUADbobbb2N3FvLcP+vZBWW4MYFkX+AW13VZuLeDxIbinvF4ohOoWb2JpBSobL4M9CEGBnAQa6M92qe22Mq1Awma81094//33z5492/yzJEnr1q1btmzZ119/XVlZaX3x0KFD77333rffftvU1CRTnHBlI02C6roJsSs9KNNwtW3SK3LKKOIjsHL9dRMOHz48cODA3bt3E1F9fX1SUtL8+fOPHz/++eefr1+/nohWr16dkJBw6tSpN954Y8KECXIHDH8QF8z2lPB6Ue442owTZeQ3oG8fFCg2iKXnqmlm7oFSbnZnAa5yx6F+xmXLlplMpldfffW7775btGhRU1NTZmamk5OT9W3O+fPPP79s2bI77rijtra2e/fuaWlpI0eOlDdoaOZhpD6+bHcxV0XzSE1NzYDRN5eW1hwMfi0sOUnucAD+IOOLpes/+O+fVt6wbs1/5I6lTdBBaCsCEU2cOHHdunWc89WrV8+ZM2fv3r0//fRTaWkpEWVnZ2dlZY0fP56I3N3dk5KS1q5dK3PI8EcqmkRx/vz5PNFTHPXgzxs2yx0LwKXWr1vLpy7O2LOfc3UkVFoBH4leBlswElFoaGhDQ0NxcXF2dvbbb78dEBDg6up69913//TTT/X19b6+vm5ubtZPh4aG5uXltfZdDQ0NK1asMJlMLV9MTEwcNGjQFT9vsVgsFrUN0lKeuEB6/xg92UdU/vEMi4g09EkYa9nz5N+fVnioFotFkiRBQO+LbSj/5CSif//zlfsWvCfOfEn5gyEsFkujxbKtkN4eShaLerpG5GA0Ghm7xu2CkYisHxJF0WKx9OvXb/ny5UT0wgsvPPPMM88991zLrxAEQZJaHaPIOa+urr5w4ULLFxsaGlr7K5IkXeXboI1i/Gl2kdDYJCn/eH6cxeJnPPL1TRIRKTxUSVLB8VQRVRzM/v36bfvq332/FzKK+LAART8USpJ0pFzyMgrBLso/ripgJKKCggInJ6fg4GCTyRQbG2t9Iy4u7qOPPgoJCSkvL29sbHR2diai/Pz80NDQ1r7L1dV17ty50dHRbfy3GxsbXVywMev1MrlQV6+mI9XOPV1dlHw8RU7vHm1aOcrg4qKCxhzrPZ+Sj6e6WCwWVRxMF6L/6y+9fYR/O9YgdyxX09jYuCPHaVQId3FxkjsWLRCI6Oeffx4zZowgCElJSUePHrW+ceTIkYiIiMjIyPDw8JSUFCJqbGxMTU1NSEiQM164kvgQFXQTfn5KCvcgrKwGCvfXKGFXMf+1XOkJlVagjiFyqmB87rnnli1b9v333xPRvHnz4uPjRVF0c3N7//33v/zyS0EQnnvuuQceeGD//v1paWldunRJTEyUO2a41EgT++9Jfn+E3HG0jhO9cUB6K0bRd9kARORqoEf6CIsOSJ+NVvTpmlbAXx6MPmzbENzd3bdv3x4fH09EUVFRmZmZZrPZ398/PT3dWvPuvvvuNWvWMMamTZu2fv36a/Y6guONMglpBZKk4FvY785ITgIlmnHygAr8rbewLkc6eUG5GXWmhklE3byQULZhfPrpp1v+HB4ePnfu3Es+FBcXFxcX58CooH2C3CjQlR2tZDHecofSijcPSv8YiA0IQR28neiBnsLiQ9L7wxX6ULi9SBiFdlHbwZO1RsSbKOVEmdxRXFlqLi9vpIldcLKBasztZ1h9WsqrVehD4aasUnQQ2hC2ANACzvnPfx9TYnFqmv6nl5+ZJ3c4l3r9gPh0tIDnQVARfxea2UNY8qv0xlDFPRSOmzJzw4my8t6hD37xodyxaARu0rXAYrHUVlU0xP5l6659csdyqZ1F/HQVTYvEmQYq83/9hI+OS2XKW3x03+Hj4k1/yz5+RO5AtAOXJy1wdnZetXzJIMvRwBmvyx3LpV47IM3rJzjhRAO1CfNgEyOEpYeVNV9d5OR9778Tq7d+9+kHcseiHWga1YjEhDH/u3Ho4HWuu4r50ECltEIeqeC7iqQvb8KcX1Clp6OF2B+aHu8neCnmFF5xTDL3jP7m/hu8vLzkjkU7cKOuHR5G/upg4dEMUTn9+6/ukx7ra3DD7RaoU6Q3uylU+OCYUh4KL1jo5X3SP4cprttS7VAINWVWD0Hi9MUpReTt6SqekiM92AvnGKjYMwOExYdEhWz5uWCveEtnNihAKU0+moGLlKYwoiWxhqd2STUKWD3/jQPS33oLnZzljgPgOkT7sYH+7NMs+W8uT13gn2RJLw3C46DtoRBqTWwQG25i/zwoc94W1NGabOmRPkhaUL3nBhoWHZCa5C6Fc3dIT/Y3mNxkDkOTUAg16M2hwtLD4rlqOfsK/3lQnNVDCHSVMQQA24gJYmEetCpbzkq4KY8fLud/74Mrtl3gsGpQmAf7W2/hmUzZ8rasgT4+IT3WF2cXaMTTAwyv7pNtOV+R02M7xLdiBBe0sNgHLlXa9HS0Ia2Ab5Npb6Z3D4uTIoTOnujSB434UxhzM9L/zslzc7n8qBTgQhOwSKHd4Mhqk5uRXhks/N8uGaZS1DTR+0elef1xaoGmzI8WXt0vQyGsaKSX94lLYvEwaEe4WmnWzO4C5/Tfk45O3fePSmNChRs64XEQNOX2CKHaQr/kO/recsFecWKE0N8PCWVHKISaxYjeiTU8vduhUykaRFryqzSvH84r0BqB0f/1F17f79AZhScv8M9OSi/ciMdB+8IFS8tigli8ib1xwHGp+/EJaaA/3YgJv6BFM7sLJyopo8hxD4WPZYhPRWPKhN2hEGrcoqHCv45KZx0ylULktPiQ9FQ07l5Bm5wEeqKf8KajJuluzOPHKunh3rhK2x0OscaFebA5vYWndzsidb84JZk9aHgwHgdBs+7rKews4r+W2/3OUuQ0d4f4NqZMOAQKofbN729IL7T7VApOtOiA9DQeB0HTXA00p7fwxgG731m+f1QKdKXxnXGJdgQcZe1zM9Krg4VHd4h2nQ78w1nJKFBSGB4HQeMe7iP8fF46dcGOu+oNAQAAEVhJREFU6VTeQK/sE5fE4LbSQVAIdWF6d8HNQP+x51SKNw5K/xggoAyC5nk70QO9hLd+tWM2vbRPvL2r0A9TJhwFhVAXGNHiYYZnM6Vqi12+f2MeL2ug2yJwOoEuPN7P8NUpKa/WLg+Fxyr45yellzBlwoFw5dKLYUFsdAhbdNAuUyle3y8+FS3geRB0wt+FZnQX3rHPQ+HjO8VnBhgCsGC9A6EQ6sjCIcK/j0pnqmx8G7urmJ+opGmROJdAR+b1F1Ycl0obbPy1qbk8q5L+hikTjoXDrSNmD/ZIH8NTtp5K8dp+6clowRmnEuhJmAeb2EV474gts6lJork7xCWxBmSTg+F468u8fsKOIp5mu6kURyv4ziLp3htwIoHuPDNAeO+IaMN+938dlYLd6M/h6GNwNFy/9MXNSK8PER7NsNlUitf2S4/2NbgZbfNtACrS3ZuNMgkrjtvmobC8gV7bL76NKRNyQCHUnTsjBU8n+jTLBtmbXcXXnpce7IWzCHTq2YHCmwfFelsMQXthrzgFUyZkgkuY7jCiJTGGZzLFC9fdpPPGQemh3oKPsy3CAlChaD82wJ99dt23lUcr+FenscuEbFAI9ejGADbWLCy6vl0pCupo9Wnp732QuqBrT0cbFh6Qmq6vFD6+Q3wWUybkg0KoU28MNXx4/LqmUiw+KM7sIQQidUHfRppYqAetzu54JfzpPD9TTQ+hi0E+OPQ6FexGj/Q2PLmrg9lb1kAfnZDm9sX5A3DxobBjN5VNEs3fJb41zOCEZJIPjr1+/V9/YXcJ39qhqRRLD0u3RQhdPNGxD0A3hzMngX4615FUWnZECnGncZgyISsUQv1yNdDCIcJj7Z9KUdNE7x8Vn+yPkwfgoif7C6/sb3ene1kDLTyAKRPyw7VM16Z2E7ycaGU7x7z9+6g0OkS4oRPuYQEumtxVqGigzfntu6l8fo84tZvQ1xepJDMUQr1bEmv4R3umUjSI9PavEh4HAVoSGD0ZLbzenofCoxV8dbb0PKZMKAAuZ3o30J8lh7UjgVdmSdF+dGMA7mEB/uCu7sLxStpd3NaHwrk7xOcGGvxd7BoUtAkKIdDCIYYPj0sn27Djtshp8SHp6WjcwwJcykmguX2FhQfa1NHw4znpXDU90BNXYEXAfwNQsBvN7dumXSm+PCWFuNEIEx4HAa7gvp7C9kLp1/Jr3FNaJJq3U3orBlMmlAL/D0BE9EQ/4UApT829WgJzokUHpacH4HEQ4MrcjfRIH8ObB69xT7n0sNTdm/4UhhtKpUAhBCIiFwO9NkSYt0sUWy+FP56VDIySkb0ArXukj/Dz+aut2VTWQIsOiouG4oZSQVAI4aIpXYVOzvTxiVZvZhcdlJ4dIKAMAlyFtxPdFyX881CrefSPTHF6pNAHUyaUBIUQfrckxvD8HrGy8QpvbcrjpfU0KQInDMA1PNbX8PkpKb/2Cm8dqeBfn5H+MRCPg8qC6xr8boA/GxcmvHalqRSvHxDnRwt4HgS4piA3mhEpvHP4Cnk0N0N8AVMmlAeFEP7g9SGGj05IWZV/6OHYXcyPV9CM7jhbANrkyWjhw2NSecMfXvz+rHS+hu7HlAnlMfbq1Ss8PPzVV18dMmRIXV3drbfe2vzelClT7r///hUrVqxatar5xR9++MHNzU2OUMERgtzoiX6G+bulbxJ/b715bb80r7/gjPwFaJtwDza+i/DeEekfAy+mTaNET+6S3o01GJFHymP83//+l5qaOm7cuOzsbMZYamrq2rVrjUYjEXXp0oWITp48GRAQcO+991r/grMz9iPXuMf7CX2/btqQy8eaGREdreA7iqT/3uQkd1wAavJUtDDyx6bH+gqeTkRE7x6WojoxDLpWJmNkZGRkZOTKlStXrVo1depUIhozZswl1S4iIiIxMVGmCMHRnAV6fbAwd4e4/zajUaDX90t/72twN8odFoCqRHVio0KED49Lj/UViuvpjQPitvHIIoW6+JQ+ePDgQ4cOWf986623jhs3btGiRfX19dZXvvvuu9GjR8+aNSszM1OeMMGxbu8qhLjTRyek7Cr+03kJe2cDdMCzA4TFh6RGiZ7LFO/qgQ1blOviHYqPj8+pU6ecnJzeeuutIUOGlJaWLliwIDMzc/Xq1aNGjYqNjQ0MDNy4cePIkSPT09MHDhx4xe+6cOFCYmKik9Mf2tCeffbZWbNmXfHz1dXVtv1ldM62x/PFPuy2b/L/U3N49ogxhob6qoZr/xUtsVgskiQ1Nl5pKgm0X01NDecd28JdxSKdqYexfsLCzXv8Ruy9w6uq9Vn27YWLZ9u5u7sbDNeYr3KxEFZXV/v5+bm4uMydO9f6Sp8+fXr06FFeXn7zzTdbXxk+fPjJkydXrlzZWiH08vL64IMPevfu3fJFHx8fDw+P1v55Ly+vNv4y0BY2PJ79qap6ye1p3UZEle/xGvWyrb5WLayF0MUF49xtgzHm6ekpdxQyKHxv2pFOg4Kzl4Q9uN2234yLpw1dLITHjx9PSEho+UZQUBARVVdX+/r6Nr8YGBhYXl7e2ncxxoKCgsxms31CBUfzdGYikzyc0J4D0EE+LgYnxn1c0LmgaAIR7d+/f+vWrXfeeefx48dzc3OJqKGh4dlnn+3Zs2dYWNiWLVtEUSSiAwcOfPbZZ2PHjpU5ZHAILy+vzNQfv5s38Z8vPy93LABqlfL1f9fcOyRj/Q9yBwJXYwwODm5qalq6dKnZbP7222/vvfdeznljY+PgwYO//vprxtgLL7ywY8cOd3d3xtjjjz8+bdo0uWMGB+ncuXPnzp3ljgJAxTw9PW8dP17uKOAajFlZWe7u7taJg7fddtttt91WWVnp7u7ePOZl8+bNjY2N9fX13t7esoYKAABge8bLy1unTp0uecXZ2Rnz6AEAQJPk7MKtrq7W4YhqO5EkqaamRu4otMPaCiJ3FNpRW1trHWoANlFVVSV3CJoiZyEcNGhQXl6ejAFoyeHDhy8Z9wvXY+nSpS+99JLcUWjHHXfcsXnzZrmj0AiLxWJd/xJsBYN6AQBA11AIAQBA12y5CKwoiidOnGj75y0Wy5EjR0pKSmwYg25lZWXV19cfOHBA7kA0Ij8/v6ysDMfTVqqrq0+dOhUQECB3IFrQ1NTEOcfJ2UZRUVGurq5X/wyz4XCVxYsXf/zxx9dc1a3ZhQsXvLy8GMPCJTYgimJtbS1WXbKVhoYGzvk18wfaqKamxsXFxTpNC65fZWXl5cP74YpWrVoVFRV19c/YshACAACoDvoIAQBA11AIAQBA11AIAQBA11AIAQBA1wwvvviiY/6lhoaGbdu2nTt3LjQ09JKRpZzzvXv3Hjp0KDAwEOP02qiurm7Lli0FBQWhoaGC8IcbmnPnzu3atausrMxkMl3yFrTm8OHDe/fu9fLyuuL+sZIkZWdnM8ZwfraFKIo7d+48duxYSEhI8/L9zWpra7du3Zqdne3j4+Pm5iZLhOpSU1OzZcuW4uLi0NDQS4bZ19XV7dy5Mysry9vb293dXa4IVY87RHFxca9evUaMGBEbGxsdHV1eXt78liRJkydPvuGGG2655ZagoKC9e/c6JiRVO3v2bOfOnRMTEwcOHDhy5Mi6urrmtx588MGgoKDk5ORevXr16dMnPz9fxjjVYt68eWFhYRMmTAgICPjpp58u/8DSpUsFQXj11VcdH5vq1NfXjxo1Kjo6OikpKTw8PDs7u+W7mzdvDg4Ojo2NTUpKiouLkylGNTl58qTZbE5KSoqOjh4zZox1Yo9VVlZWWFjYmDFjbr/9dj8/vx9//FHGOFXNQYXw+eefv/XWWznnkiQlJye//vrrzW9t2LAhPDz8woULnPNXXnnl5ptvdkxIqvbQQw9Zd460WCxDhgz58MMPm9/KzMxsbGzknIuiOHbs2CeeeEK2KFXi5MmTHh4eOTk5nPPPP/+8Z8+ekiS1/MCZM2f69euXmJiIQtgWK1euHDBggPUkfOihh+67777mtyorK4OCgr766ivrj9aJ4XB1d99995w5czjnDQ0N0dHR//nPf5rfevTRR6dOnWr98zvvvDNs2DB5QlQ/B7WbffPNNzNnziQixtj06dO/+eablm9NmDDBOhN8xowZ69atq62tdUxU6tV8PI1G49SpU1sez0GDBlkbowRB6NmzZ2VlpWxRqsS3334bHx9vNpuJaNKkSWfPnj1y5Ejzu5zzBx54YPHixVdsMoXLffPNN1OnTrWehDNmzGh5cv70008mk2nSpEm//vpraWlp2xff0LPmZHd2dp4yZUrL4+nq6uri4mL9s7OzM9rtO8xBhfD8+fPNe5137tw5Jyen+a2cnJzmt8LDw4koNzfXMVGpVGNjY1FRUWvHs1leXt5XX3115513OjY69cnJyWley9/FxcVkMp0/f7753RUrVoSEhIwdO1am6NTnkmQvLS1t3iDs5MmTLi4uAwYM+Pvf/x4VFbVo0SL5wlSHCxcuVFZWtpbs8+bNKysru+OOO+6///5PPvlk6dKlMoWpeg5a8aihoaG5z9zZ2bnlTm8t3xIEwWAwNDQ0OCYqlbI2OjUfNBcXl8t3zquqqpo0adKsWbOwN9M1NTY2tlz6y9nZufkMzMvLe/PNN9PT02UKTZUaGxtbJjsRNTQ0eHh4EFFVVdX+/fsPHjzYs2fP48ePDxgw4LbbbrvhhhvkDFfZrKdia8l+4sSJY8eOzZw509vbOz09fceOHf369ZMnUJVzUCEMCQlpXly7pKQkJCSk+S2TydT8VkVFhcViafkuXM7T09PT07OkpMT6AF1cXHzJEaupqbnllluio6PfeOMNmWJUE5PJdPTo0eYfWx7PxYsXBwUFLV68mIgOHz5cVlYWEhJy9913yxOoSrTM6OLiYldXV19fX+uPISEh3bt379mzJxFF/f/27h8knTCMA/jriZHL4XkUYRCVxlFecEmWJ0E0iENgRRJRFDg6BOEoBAoGDkFBYlAY0RC5NgQHQULkEiQ4tEhDBUERBC1RR70NwiG/H/2hQZP7ftb3hoeXe3ngfZ/3eQXB4XAUCgUkwi/wPG8ymR4eHpqamsh/iz0ajYbD4UgkQgjxeDw+ny8UCqGh6y9UaWtUlmXtWc7j42Ov16sNeb1ebSiXywmCYLVaqxNV/aqcz1wuVzmfz8/PgUCgo6NjY2MDDc1/Qpblk5OT8vvpxWLx7e3N6XSWh/x+fyAQ4DiO4ziTyWQ2m3FS+K1/VrQsy9p/ODQ0dH9///r6SghRVfXu7q65ublWcdYFhmE8Hs9ni/3l5UU7F2xsbFRVlaJ39O9Upybn7OyMZdn19fXV1VWWZYvFIqV0eHh4ZWXl6emptbV1cXFxd3e3ra2tsgASPqMoCsdxW1tbiUTCYrFcXV1RSnt6enZ2dubn53meX15eTiaTyWRyf3+/1sH+de/v7y6Xa2ZmZm9vT5KkaDRKKd3e3hZFsfKz8fFxVI3+xM3NDcdx8Xg8k8nwPH94eEgpDYVC4XCYUurz+crlXdPT0263G4Wj3zo4OOB5PpPJxGIxq9V6e3tLKe3s7Mxms6lUqqWlZXNzM5vNulyuubm5Wgdbr6p0od5ms42MjBwdHT0+Pq6trUmSRAhhGKavr89utweDwUKhcHFxsbCwMDs7W4V46p3dbne73YqiqKqaTqcdDgchxGg0DgwMsCwriqJ2j95isYiiWNNg/zqDwTA1NVUqlc7Pz4PBYCQSMRgMDMPYbLb+/n7ts4aGht7e3vJ2NHyBZdmxsbHT09PLy8tYLOb3+wkhDMN0d3d3dXVNTk5eX1/n83lJklKplFb0CJ8RBEGSJEVRCCHpdLq9vZ0QYjQaBwcHR0dHnU5nPp8vlUoTExNLS0soxP0dPMMEAAC6hv5bAACga0iEAACga0iEAACga0iEAACga0iEAACga0iEAACga0iEAACgax+Oj2GtzbP87gAAAABJRU5ErkJggg==", "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: 6%|β–‰ | ETA: 0:00:02\u001b[K\rDiagonalising Hamiltonian kblocks: 7%|β–ˆβ– | ETA: 0:00:03\u001b[K\rDiagonalising Hamiltonian kblocks: 19%|β–ˆβ–ˆβ–ˆ | ETA: 0:00:02\u001b[K\rDiagonalising Hamiltonian kblocks: 28%|β–ˆβ–ˆβ–ˆβ–ˆβ–Œ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 39%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–Ž | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 48%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–Š | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 59%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–Œ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 72%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–Œ | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 83%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ– | ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 94%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–| ETA: 0:00:00\u001b[K\rDiagonalising Hamiltonian kblocks: 100%|β–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆβ–ˆ| Time: 0:00:01\u001b[K\n" ] }, { "output_type": "execute_result", "data": { "text/plain": "Plot{Plots.GRBackend() n=210}", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOydd0BN/R/H3/feJkVJCWWEQsPITLZC9iibIiPy2DtEHjt7ZpOVbCFZj9DzICMtsiKziLRv935/f9x+LdVdZ9xb9/VX59zv+Xw+1T3nc77f72dwCCFQoUKFChUqyitctg1QoUKFChUq2ETlCFWoUKFCRblG5QhVqFChQkW5RuUIVahQoUJFuUblCFWoUKFCRblG5QhVqFChQkW5RlJH+OXLl/nz59+5c4dWa1SoUKFChQqGUZNw3JQpU+7cuVO1atUOHTrQapAKFSpUqFDBJBLNCAMCAjQ1NVu1akW3NSpUqFChQgXDiHeE379/9/b23rhxIwPWqFChQoUKFQwjfml06tSp8+bNq1atmtiRJ06c8Pf353ILOdexY8f26tVLcoMSEjj6+qRiRQCIjeUIhWjcmAiFePmS07ChclSDEwgEPB6PLe2JiYiN5dSuLdTV5ejrS315VFSCu3vVJk2+/v23sYGB5p8D+HwSF4fGjTkU2Col8fGcatWIllbu4bdv4PFgYECLLqEQL15wGjWS9CtHCCGEFPnyl8S3byAEEtxSioVQKORwOByO1P/6hAROpUqkUqXcw69fweOhalWKzUtNRWKisG7d4v8FkZHfpk7VrFkzY8eOCpXyTJGY9+85uro5+vq8Fy84Fhbk3TtO1apER0duo6UnJQW/fnFMTaV4GMr8RIqN5VhYkLx/+MuXnHr1COXPtvh4TpUqAl1duiI3eTye2BtTjCO8cOFCUlLS6NGjJdG3fPlyd3f3unXrFjxpZWUl+f/g61fY23NPnhTa2YEQdOrEbdAA9+8LCUHHjtzoaKGRkYSS2CQ9PV1XV5ct7X//zTl/Hk5OxNJSbepUqV8d/Pz2PHmi9ezZ5RMnzpqa6jRtSpo3R7NmaNaMVKkCANeuvfXwELx7V59608Uxdy6nRw+4u+f+UocPc96+xa5dtLweEYL27bnv3gkrV5ZoPJ/PFwgE6urqkgzeuJFTpQrmz1eOF7s8+Hw+j8dTU5M0sCCPceM406ejd+/c33fpUk7TppgyheJf/7//ODNnnnnwoHuFChUAvH2Lp085T57gyRM8ecJJSbmekfGvltbHKVNm29vbSyt88mSOu3tWnz4a3btzL14U+vtzHjzgXL4sZN4XhoZy/Pxw6ZIUfz3Znkhv38LRkZuQIBQ5wowM2Nlx4+OF0r9FiGH69G8NG4avXduTYrkAAKFQKBAIxL+hklJxd3c3NDQ0MzMzMzOrUKGCgYHBpEmTShpsbW399OnT0gWWglBIevcmS5bkHj5/TipXJitX5h4OGECOHpVZNqOkpKSwqL1ZM9Kzp3Dz5ix3d1ku//3794YN22/f/ofPJ8+fk0OHyPTppGNHUqkSqV2btG9/j8OpB9jPm+dFteHiuXaNNG2af/jxI6lShaSl0aWuY0cSHCzp4Ozs7IyMDAkHd+1KrlyR0SoWycjIyM7OlvaqX7+Iri75/Tv/TK1a5OVLKg0TYWLSCujP4zXt0oXo65OaNXOfJ2fPknfvCJ/P3737QGDgWdmEV69OoqN/E0JWriSjRxOhkEyaROztSWoqpb+DBHz+TAwMiFAoxSWyPZG2bCFjx+YfXrlCOnaUQYwY3r9/z+HU4HDsBw1yo146IQKBQJIvrRhHmJaW9uP/ODg4LF++/HfBb3Rh5HSEmzeTli1Jns0bNhB9fRIennu4axcZPVpm2YzCriPU1SUbNwquX89o3ZpKsUIhiYsj06Yd43DsgMkDBoygUrrENpibk//+yz/Tsyfx96dL3YIFxNtb0sFSOcKqVcnnzzJaxSKyOcLAQNKjR/5hVBSpU4dKq/LQ0WkA7OVw6l+9Sr5+pVLy9++kcmXy61cKIeTHD6KvT96/J0IhcXcnDg5E4n87ZZiaktevpRgv2xOpe3dy+nT+4fTp+dMSCnn06BGHYwksadasG/XSJXaEYiaMFSpU0P8/6urq2traOvSsBURFYcUKHD2KvLWloCAIBGjaNPfQyQnBwVD1jCqdjAykpmLwYGJlRaKiIBRSJpnDQf362LRp2IIF3YYOzTh8eBdloqWxYfx47NyZf8bNDQcO0KWubVuEhVEvNj4eGhowNqZesmISFISCQQLBwejRgxZFoaEBjo5B589v6N4d1O6hRETAxgaiFUJ9fYwZg23bwOFg1y4YGqJ/f2RlUalOLK1a4cEDelWkpeH+fXTtmn8mOBjdu1OvyNbWdsmSyZqaBlevHqVeuuRI7lo/f/788+fPUgbIPCPMzCRNmpBDh/LP8PlEW5s4Oxca1rAhefxYBvFMw+KM8OhRoqVFcnJy0tLSatcmr15Rr4LP56enp1MvVzKSkoi+Pvn+PfcwK4sYGkr3diw5iYmkcmUiEEg0WPIZ4ZkzpFcvuQxjCxlmhEIhqVGj0Pewe3dyVsblSfHQdOtt3kwmT84X/u4dqVKFiJ6FOTnExYUMGED4fDo0F8/q1WTmTCnGy/BnOX+edO2af/jhAzE0lPRekBaBQFixIqHpqUnNjLAgxsbGlSWMHJCSOXPQoAEKRuSEhUFLC337FhrWvTuCg+nQX3YICkKtWrk/W1vj+XNWraEBAwP07o3Dh3MPNTQwZAj8/WnRVbUqjIwQHU2x2CdP0Lw5xTIVlvBwVK6MevVyDzMzERaGzp1ZtUl6nj+HtXX+Ye3acHTE3r0AwOPhyBHw+Rg2DDk5DNnDwIywyDz+6lU4OkKymGip4XBQp47w1StahEsI+7VGg4Nx7hx27y50MiQEGRmFJuZQOUIJePQILVrk/mxjUwYdIYBJk7B9e/4iuWh1lMJF4ILY2VG/OvrkCZo1o1imwlLkeXrnDpo0AT2v0zRSxBECmD0bmzaBzwcADQ2cPo20NLi70/U9LEKLFnjyJFc7TVy9WnRBm4510Tzq1SvfjjAxEePG4dAhiELz8zh3DqamqF690MlOnfD4MVJSmDRQyXj/Pv/rWyZnhADs7KCri1u3cg+bN0flyqCpCC4d24SPH5dfR0j385QOhEJER8PSstBJW1s0aICAgNxDDQ0EBuL16591645zdp6UlpZGq0m6uqhTB1FRdMl/9gxqajA3zz0UCHDrFrp1o0sdgHr1SFwcjfLFwqYjJATjxsHNrehSye/fePEC/fsXHa+tjTZt8p+AKorw9SuysvLXk62tERHBqkG0MWECdhUI1nF1pStkhnJH+O0bMjJQuzaVMhWWb98QF4d27fLPKKMjfPsWVapAT6/o+dmzsXZt/spEhQoYOfLChw9mFy7oX716lW6raF0dDQpCnz75hw8eFDMtoRYzs3I8I9y2DZ8/Y8mSoudv3oS2NootR6NaHS2FEyego4O8qF4LC3z4gIwMVm2ih5EjceMGPn3KPRw1Chcu4Ncv6hX9+nX/5cv2PXuOFQgElAgUTQelr82ilFy+DAeH/Djwjx/x9avy7Y/+uS4qomdPCAS4cSP/TPfuHUxMrhByv2XLtnRbRbcjZHgeX78+KY+OMDs7OyTknY9PoXyJPC5fRnY22hb3XereHfS/bCkrISGoX6Dei5oaGjRATAx7BtGGjg5cXLB/f+6hgQE6d0ZgIPWKdu48LhQuv3//W0JCAiUC//svtTyvizo40BVwQR/Pn8PGppjzHA5mzoSvb/6ZOnXqvH9/t3Xrfx49qkG3VfQ5wh8/EBWFgk2GGHCE5XFGmJOTY25u17PnlDZt1uYtQxckKAitW0NDo5iPLC0hEIDd1WSF5dkz2NkVOlOGV0c9PODnh7x5Gk0JhbNmuZmY/K2jU9/U1FR+ab6+u1at6nHqVFdSDvJh+XzcuFEoZVAZ10VR8owQwIgRiIjA06eFTk6Zgm3baLfK2hpv3+L3b+olX7mCzp2h+f8yw8nJiI4u+mChnBo1SEoKLb+OhLDgCLOyspKSBEJhD03N139++vEjkpLg7Fzi5Y6Oqklh8Xz5ggEDCp0pq/EyAGxsYGqKy5dzD3v2xJs3iI2lWIutbfNXr67n5Gx684aCO+X69bCsrOkpKcnp6enyS1Nw7t5Fgwb5hcVFq4gODqzaJBMRESU6Qk1NeHpi06ZCJwcNQlwc7S+g6uqwscHjx9RLLjKPv34dHTrk+0Wa4HBgZobXxTgEhmDBEWppVdTT2zZhAn/XrpV/fnrtGrjc0t4cVduExRIVBYEAnToVOlmGHSGASZPyQ2bU1DBiBA4dol6LpiZGjqRmutmr17Lq1cP8/JZUFHVXKdMUeZ4+eoSaNVGD9iVDisnMxPv3aNCgxAGTJuHCBXz4kH9GXR3u7oXqH9EEHaujAgFCQtCzQPlrxubx9euDxdVRFhxhQADq12+3a9dMg+I66Jw5A13d/AzcP+nWDaGhTNc0UnxOnoSBAYr0+SjbjtDFBeHhePs299DdHUeOgKKglkJMmID9+ylI2zpzps7Gjb7Ozn/EQ5dFykDiBICoKJibF79NI0JfH6NHF10L9fDAyZNITqbXtpYt8fAhxTLDwmBqipo188+EhDD0j2vQgM09LxYc4fr1mDu3+I8IwZ07YkoR6unBygp379JhmhJz+zYaNSp60sQEfD6+fWPDIPrR1MSoUfDzyz20sICJCa5do16RhQUsLHDhglxCXr9GVFQxSUFlkrdv8fNnoQBRJXWEpWwQ5jF9OvbuLRS0bGSEnj1pWZ8oCB0zwiKvL9HR4HJLmxBTSPmaEQYHg88vNPUuSGQk+PzSNghFqFZH/yQmpui6qAhLS0RGMm0MY0yahP3785cH6KvBPWFCvseVjZ074eZG+16LgnDpEpyc8rNEUlIQFQXpmwCyjySOsE4dODpi375CJz09sW0bvYVm6tVDaio+f6ZSJmMV0v+kfDnCtWsxd26JeVSXLxez0fUnPXqo4mUKIRTi+3cMHlzMR2W10JqIevXQpAnOns09HDoUISFISqJe0eDBePpU9hs1OxtHjmDcOEptUlQSExMPHw5xcMjfvQgJQbt2SvkSIIkjBDB7NjZuLLR43rYt9PXpfV/ncGBrm3XzJmVf96CgRwkJz1u2zD/D5Dy+HC2NPnqEuDgMGVLigNOnYW4OsY2eWrTAly+FNqjLObdugccrPtupbG8TApg0KT8wQVs7m8frXa9ei6Agip9AGhoYOTI/c1FaAgPRpAlDS0zsQghp1qx7ePiV48en5p1U0nVRlJxEWARbW9Svn19xTcTkyfTmUSQlJYWFtZkwof/p0+fll3b58hUXl8Xp6VPCw3PXWxmukF6zJljMoGDUEa5Zg1mzismgF8Hn49kzDBwoXg6Xi65dcf06tdYpMWfPltjfrgynEoro2xdv3uQu/3769Cknh5+SsvDkySuUK5owAQcOyBgys3s3Jk6k2iBFJTlZyONVBPL/UtevK6UjTEpCVlahyJFSKFJxDcDw4QgPx8uXNFmHb9++cThGGRnd79+noGrGly85mZmampqaOf9vosFwhXQOB3Xr4s0bhtQVgTlH+OYN/vmntNWhsDBwuRI5Qqi2CQtz/36J761WVoiJYagoPiuoqWHcuNwNvDp16kyd2k1f/1zLllPFXSc1MofMxMYiLq5oT7Gyyo4dnNq1gw4fbnb8+HbRmdhYCIWwsGDXLlkoJYPwT5ycIBDg5s38M5qacHMrVBSXWho3brxjx7h27bSeP58iZ4UGPh/79/eZOHFGcLCP3f+T55mfx7O4OsqcI1y7Fh4epS17nj5d4vren/TogZAQWmLllZFXr+DoWPxHurowNGTtPYsZJk7E8eMQVfz38Zlz7NjhXbvq0eH7ZQuZ2bULY8eWuBBSlnjwAD4+uHix5rBhAytUqCA6WebXRUWIKq6tXZudU6AtoYcHDh8Gfb0oRoxwuXVrTmqqrpxrsNOno2pVbN/esU2bNnknmf/HsRgvw5Aj/PYNgYHw9CxtzKVLaNdO0nrExsYwNaU+jUYZSU9HamrxkTIiyvzqaPXqsLfH8eO5hz16QFcXZ85Qr0iGkJmMDBw9Cnd36o1RNH78wNCh2LmzaBKwUjtCyWeEAJo2jbxxw65mzTav/18ipVYtdOhAV+NoEWpq8PeHj0/RSm+Sc/w4QkJw6FChZy8rFdLLviPctAnDhsHQsMQBv3/j/XsMHy6FTNXqqIhz56CpWdpORpmPlwHg4YHt2/MPlyyBtzf1C8IyhMycOoWWLVGnDsWWKBpCIUaMwLBhRYv8ZWXh/n3la0kvQqqlUQDR0c84nM7fv7eKLJCx5OmJLVtAa3FZMzP4+mLECMhQuS8yEn/9hYCAonuBrFRIL+NLo79/w88P06eXNiYkBBxOifmFxaJyhCKCglCrVmkDyoMjdHBAenr+CoGTEypWxLlz1CuSNmSmnITJeHsjKwvLlxc9HxoKKyvo67Nhk3wIhYiJKdqPt3QGDx7k6aljYGCSnu6Ud7JLF/B4+Ocf6i0syKhRaN4cM2dKd1VqKlxcsH49mjYt+hEr8/gyPiP084OjY2lV0wAcO4Zq1WBkJIVYe3tER+PHDzmtU3oePUKLFqUNKNuphCI4nKIFHhcvxtKl1E8KpQqZiYnB27dwchI/UqkJCsK+fTh6tGiFPyjzuujr1zA0RKVKUlyipaW1cePSixcXzp6tXrC+2qRJTPSj2LkTt24VTeEoBULg5oYuXTBmTNGP2KqQbmKCnz/ZyaCg3RHy+diyBbNnixl2+7bUJQw0NNChgyqJAu/fF9/EOI8GDZCQQOOOvYIwdizOnct/MerdG9ra8tZFKxbJQ2Z27MD48WU8TCY+Hu7uOHmy+A7myusIpd0gzKNVKwwciPnz88+MHo1btxAfT5VpxaOjg6NHMXUq3r+XaPz69Xj/vlA/xTzYqpAu6kHBSmQf7Y7Q3x+NGonZdP34Eb9+YdQoqYWrVke/fkVWFvr1K22MmhrMzctmh96CGBigVy8cOZJ/xssLy5ZRvz0zeDCePRO/hpORgRMnMHYsxdql5eHDhwXjGKklMxODBmHhwuLLp335go8fYWtLk3J6kdkRAli5Eleu5KdS6OgUKopLHy1aYPp0jBwpPpz+/n34+uLkyeLL/bD4+sLW6ii9jpAQ+PqWWGI7j7NnweEU35K+dJycEBxM70a0gnPiBHR0ILarT3lYHcX/GzPlfR/69oWaGi5epFiLhCEzx4+jbVvUrk2x9jzOnDlzroRd0J8/8fQpLl5ErVrOrVvP1tGhq6TN1KmoXx9TS0javHoVDg7FrJcqBfI4Ql1d7NwJDw9kZuae8fTE3r35h/Qxbx40NLB6dWljvn7FkCE4dKjEGC4WHSFb8TL0OsKLF6GtjS5dxAw7cQLW1rKsINWtCy0tREXJZl1Z4No11K8vflh5iJcB0K4dUlMP9Ogx89OnT6IzixbB25v6V6Xx48WHzNAaJrN169bBg9cMGLDW1PRM/fqoXh16etDWhpoaOBzo68PWFoMHIyEhnZCJWVkVzcwE69eD2pmhvz9CQ0ub6Fy7pqzropAyifBPevWClVW+Q6pfH82bS7GBJzNcLg4dwrZtuH+/+AE5OXBxwcSJJf5r2K2QXjZnhGvXYt48MWMIwePHGDRIRhXlvAB3RAT+XwiiNMp8KqGI9+/fJycfuXbN2strvehMv37gchEURLEiCws0bIjzJZd4jIjA5880Vu7/998cQoyASpUq/erUCW5u8PHB0aMIC8PXryAEAgGysvDff96tWx/x9p7RrBlvyRJoa6NjR2oiGCMiMGsWzpwpMZxEKMSNG+jWjQJdzJORgYQEeWvDbtuGnTsRHZ176OmJzZvlN008NWti3z4MG1Z8Q8T586Gri4ULS7z8+nU2K6SzFjhKqMPa2vrp06d5h6GhxMyM5OSIuerZM8LjkTdvZFR64QLp1k3Ga2kiJSWFMV1qaiQkpOjJnJyctLS0gmc+fiSGhtRo5PP56enp1MiimoyMjIYN2/N4rRcvPp138vRp0rw5EQop1nX0KHFwyD/Mzs7OyMjIO5w0iSxfTrHGPG7fJjwead9+38aNG6W68MIF0qYN4XJJhQpk5Ejy5Yt0ejMyMrKzswkhycmkXj1y8mRpgx8+JFZW0smXH6puvYcPSdOmFAjfvp20aUMEAkIIEQqJuTn5918q7JOAyZOJs3Puz3mWnztH6tQhSUmlXThhAtm8mWbjCiMUClNTU0U/v39PatSgUrhAIBB9aUuHRkfYpw/ZtUv8VbNmET092ZWmppJKlcjv37JLoBzGHGFkJOFwinnV+NMREkIMDaV+8BWLIjtCEf7+6ZaWhM/PPRQKSZMmJCiIYi1ZWcTIiMTF5R4WdIS/f5MqVUhCAsUaRURFEQ0NMmKE7BKyssiaNcTEhHA4xMSErFlDbt68PW/evOTk5NIvFDlCoZAMGEBmzBCjZcUKMmuW7EbKBlW33r59ZPRoCoQLBKRdu/zH4Pr1ZORIuY2TjIwMYmNDDh4k5P+Wv3xJjIzIo0diLqxbl8TG0m9fAQo6QqGQVKhA/n9EASw7wpgYYmxMJHlgNmxInJzk0tu5M/WPOXlgzBEuXkyqVi3mfLGOsFOnYuaOMqD4jpAQ0qlToZewU6eIrS31k8JZs8j8+bk/F3SEfn5kwACKdYl484ZoaRWaicrDs2ekTx+iqZkOWAKz9fV73btX2niRI/z7b9KmDcnKEiO8fXty7Ro1dkoOVbfe9Olk3TpqhMfEECOj3Lei5GSir0/N+6gkREYSQ0MSG0tSUlLS00nTpmTPHjGXxMYSU1NGjCtAQUdICLG0JM+eUSZcQkdI1x7h6tX46y9oa4sZxufj1SuMHi2XrnKbRHH7Nho1knRwOYmXEbFtG5Ysyc8pHDQIfD71X5Lx43HwYDEhMzSFySQloUkTWFnh2jVqBNrY4MIFJCUJuNzfwM/fv2t26AAuF1WqoEULTJqEM2eKBjrevs3Zvh2BgdDQKE3y79949kwpW9KLkCdktAgNG2LiRMyYAQB6enB2xp491EgWi6Ulli7FiBHIzoaHB6ysxNe8DQ6WrrwXHbASOEqLI/z4EZcuwcND/Mhbt0AIeveWS125jZeJiUGnTpIOLleO0NISAwZgxYrcQw4HixZh2TKKtRQbMvP0KRITqQ8SSU1Fo0YwNsZ//1EsWUdH58mTi+vXN/z9e3NODiIjMWcOjIxw5QqGDYO2NnR00LAhLC1XaWubOzq2HzfuhtgWfTduoG1b8e/BCguFjhDAokWIjMz9nkyZgt27KQ7fLYUpU/DsWfOqVVsdPjxh927x44ODS+xjwxisxMvQ4gg3bICrK/T0xI88cACmpuLT4ErHxgapqWW809CfCIX48QPOzpKOLyephHmsXImjR/NTawYPRloaZXOpPP6sMrN9OyZNojh5LicHlpZQV0dEBC11kG1sbGbNmqWlpQWgcWMsWIDLlxEfj6wsfP2KjRvRrBni4gIAD6DJunW3unfHzJnYswf37hUfmqi8BWUAJCYiJ4fKuiqamti1C3/9hd+/YWODatWejRmz6evXr5QpKJmICOTkfAEuEhLcti2GDsXy5Th1ClFRyM4uOjgrC/fuoWtXBuwqjTLiCJOTcfAg/vpLosG3b1Nww3A4cHCg/hmn4Ny6BQ5HivdWUYfe8tPBsUoVzJ+fX+qdy8XChfD2pljLoEF49ix/Jef3b5w+XUzxRjlp2RK/fyM6GlpaFEsWi5ERxo/H8eO4fn0zj7dbXT04NPSv6dNRowb++w+zZ6NuXVSvjq5d4emZW+5y7twN+/d309G5KV66QvLsmVwZhMXSoQMcHLB4MQghr16NOn68orPzFIp1FObZMwwciB490LChHYfTrWlT84MH0a8fsrNx4gQGDYKeHho1wuDB8PLCiRO4dSuxZUsXDmckj8dGrc8CsJNTT9mm5P+DZVasIG5uEo1PSSFcrvgoJkk4doz070+BHEpgJlhmyhRiYlL8R8UGyxBCzMzIixfy6lWKYBkRfD6xsiKXLuUeCgTEyoqaiKGCiEJmRMEyO3aQwYMplt+uHalYkXz6RLFYGchLnyjChw/k2jWyaROZMIG0by/gcJoAcS1a9GLeQkpuvQ0byNSp1Av/+ZOYmJD794mFRTtNTZfq1Wf++iWPvBJ5/Jj0709q1CAbN+aGKxZreVYWef6cBAQQb2/i4kJq1twL+GpoLDh37hwtZpVMkWCZ+HhSsyZlwtmJGn3w4Fn16iQyUqLxu3cTTU1qYvmSkoieHpHg92UCZhxhs2YlRtuW5Aj79yeBgfLqVSJHSAgJCSH16pHMzNzDo0eJnR3FKl68IMbGJC0tOyMjo1kzih2tszNRVyfPn1MpU2ZKcoRFcHObZWLSyt9f7q+a9FBy67m5ET8/WoQfO0asrUlaWnZUVPTMmaR5c/L1q5wiCxERQZydSbVqZPXqQhH7klgeHx/fuHHnpk27JyYmUmmTBBRxhAIBlRkU7ESNnj9fpWVLSZt4HT8OGxtJW9KXjoEBzM1LrCpUJnn1Supt7XIVLyOiWzc0aoStW3MPhwzBjx/51ZApwdwcDRvi4kXuw4fclBTxBQUlx8MDZ8/i9m1YWVEmkwH271//4cN/I0bIWiyKbeQsrlYKw4ahTh1s3areuHEjX1+4uKBtW2r2wyIi4OICBwfY2uLtW8ybJ3WkUq1ataKibj55crVq1aoUGCQHXC7q1sXr18wqpVAWIdwjRwzFltjO49Ej2Sur/Um5SqLIyMhtqikV5aTQWhE2bcKaNfj8GQB4PCxahCVLKFYxYQL27uXu3cubOJGyYJbFi+Hnh/PnJSqhp4IqBALExKBxY7rkb9uGNWsez5zp+/nz53nzMHs2OnXCs2eyiOLz+QMGuDdq1KNnzxeOjrK7QAWE+XgZKh1hYqKlUHigXTuJBr97h7Q0KpvUlCtHePYsNDSK7wBXCuVwRgigXj24uuY7v2HDkJSE27epVNG9e+qtW938/ad06SJZL7hSmTp1vp5e15Urnx86VPab+ioar1/D2Bi6unTJNzUlfL7rpvsaR9YAACAASURBVE36Q4dOBeDhgQ0b0KMH7t2TWtSFC88vX+bHxrpnZx+Pjy8jLlAE8/Ey1C6NNk9O3pWRkSHJ0F27UKUKDA0p092mDeLjc1/8yzxBQbL092nQAJ8/IzWVBoMUmyVLcPkyHj4EAB4PCxdSnFO4dauvQFBPKDT19fWSU1R2dvb27Sd//VqooeE6ciQl1qmQgogIKjMI/4TD4ZiYVObxrqir1xWdcXHBsWMYPFiK0vAfP2LiREyaZGlqmmluvm3jxoFs1cimCeWeEaqrX7K0rKkt2WvJpUsUr/nweOjUCTduUClTYXn0CC1aSH0VjwcLi/xy+OUHXV0sW4bp03P7MY0YgU+fhAcOPBcKhZTI79evn7p6MI+3Z9gwifM6S0BNTQNoyOFM79xZ4qJBKqiD2lT6ElTcOnz477i4denpuWc6d8bFixg/XnyTy9RUrFmDJk2gpYW4OM1Xr06+eHHbhqYtTfZQbkeor//dz2+lJCN9fbdHR3s4OLykUDuApk1frlu3+uVLisUqIPHx6NVLlgvL5+oogLFjwefj5EkAomz3QRMnbujSRV6/JaJp06Zpaa9+/ozs06ePnKI8PFChwpWkpLDLl/0psU2FVDDgCNXU1IYNM2/TBuvX559s0QKhoVi5EuvWFX9Vdjb8/NCgAcLD8egRNm+WqGKJkqLsS6MSkZiYOGfONkIc160bRa3kvXtdIyLqOzm5UitW0fj8GdnZ6NdPlmvLrSPkcrFpE+bORVoaAAiF3/j8rjExVFb3UFNTk1NCXBz27oWfH6pU0aHEJBXSwoAjFOHri61b8e5d/pl69RAaCn9/TJtWqJW0UIhTp9CoEU6dQnAwAgJK7CxfZjAxQXJy7q3KDCw4wooVKwJVgOWmptWoldygQS01tVNVq9aiVqyiERAAHR0Z69KVt0JrBbGzQ7t2uW/c1675z5+fBBy5dYttswrQsyeaN8fw4WzbUV5JS8Pnz/L245UQExNMnVq0b3n16rh9G+HhcHPD16/fBQLB9eto0QK+vti/HyEhdOV1KBpcLurUYbRqprzvsDJw4UIFDuefXbsuubpSHBJ37drxsWNfm5vXo1asonHtGurXl/Ha8plBkce6dWjWDG5uqFu37qpV03v1grMz7t2DmRnblgFbtiA+XpboQRVUERUFCwuK68SWwty5aNwYt28XKp2vr4/gYLRqdbBWrQM8Xmbt2v+sWqXVvz9DJikOotVRZmbnYH5GSAg8PeHiojZhQn+N0lu5SA+Hw3F2rh8SQkWKvgITEYG2bWW81tgYHA6+fKHUIEUiOjr6zJmz2X9WFAYAmJjA0xPz5+ce2ttj3jwMHIi8sAW2SEnBnDmYNw/VKF4lUSEFjK2LitDSwrp18PQs2oyiYkXY2j7j84fweLyQkJ/l0AuC8XgZph3hypX4/ZvGjlydOuHRozKeIfDlCwYOlP3yMjwpTExMbN9+1NCht93d15QUEDp3Lv79F3fu5B5Onw5bW4wfz5iNxdO7N4yM8vtGqWAFhh0hgEGDYGKCPxsk+founDEj8dCh+SYmxowapDCUZUeYnY2//8a0adChLRSgQgW0aoV//qFLPus8fw6BAB07yi6hDMfLREer/fwpAJKDg7WMjDBkCPbuRXx8oTHa2li9GtOn5zfi2L4dcXGFQvgY5upV3LuHc+dYM0CFiIgIFjbhNmzA8uVISip00tDQ0Nd36cCBfZm2RmGoX5/RwFFGHeHkyVBXx6pV9Gop2y2ZAgJgYAB54hPLqiO8cwcuLvo7dgTdvDnh8+dZ0dEYPBgPH8LeHvXqYeJEnDqFnz8BYMgQCAQ3mzYdExJyC4CWFs6dw+bNuHKFBbOFQgwfjiFDYGvLgnYVBYmMZHpGCKBxYwwdiqVLmdar4DRowGwqITUlvgkh/2/DVNKnyclETY3s3k2hwuJ58oQ0bEi7llKgtftE+/akfXsxY0rqPiFi4MAF6uq2S5f6ymaAYnafuHiRGBqS4ODiP339muzeTZydiZ4esbUl8+aRihUtgcjKlS3zxty7R4yNyatXMhogasMkw4WjRpGKFQmfL6NeJpGw+wSLyHzr5eTk2NkN5nJbnD59kXLhYvnxgxgbk5KfnfLCTD8cmSnSfUKEQEC0tSnoQcFO94lSGDoUhoaYMIF2RU2a4NevQgk6ZYmYmEIxZjJw//4VPv/M4cMB1BikABw+jAkTcPlyie04zMwwYQICAvDpE1atglCIjIxswC0lxSFvs9DODl5eGDiQ0eylyEgcPYp9++Sa4quQn+/fv0dFJQqFK/39LzCvXV8fy5ZhypRC6YPlHIYzKBhyhHFxCAnBvn1M6OJw0LUrrl9nQhfDCAT4/h2DB8slZP16Ly2taXPn+lBkFMts3IglS3DrlkQ157S14eCAtWuxdOlkAwOhtnb7gQNx4kTup1OmoFUrjBrF3POod2+0bo0hQxhSp6IkjIyMWrbsXbPmnhUrprNigLs7srMRUHbeTimAydVRhhyhiwssLdGzJzPayuw24T//gMuVdz9/xIhBvXuf1dNzoMgoNlmzBnv2IDQUFhbSXbhkyfSkpEc7dw7U0MCsWfl1rXbsQGIi1qyh3NJiWLkSHz+qYmQUherVZ/v4BDSmrwNTqYgqH82Zw+iChILDZLwME47wzh1EROD4cQZU5eLoiBs38sMCywynT0vdeqlYykAGhUCASZMQGIh//oGpqYxCRo/GwoXQ1MTBg5g2DUIh1NVx6hR27JCiFYBs/PgBb2+sWAEjI3oVqZCQZ89YrttiZ4f27Rl6CVMKmMygYMIRjhyJrl0lbVtPCcbGMDXFo0fMaWSGsDA0aUKBnCZNlNsRZmdj+HC8fIkbN+Tt5OXpiZEjweXi4UMMGYLMTBgb4+RJjBtH79uokxNMTIpW2FLBFnw+4uJo7McrIWvXYudOvH3LshkKQplaGj1yBJ8+wZ/xSvqOjmWwT++rV3CgYkXTxkaJHWFaGvr2RVYWLl9GpUoUCFy+HI6OEAohFMLJCb9+oW1brF6Nvn3x6xcF8v/kzBk8fIhLl2gRrkIGYmJQty77jW1r1sSMGZg9m2UzFISyszRKCP76C6NHs7D+4+iIkBCmldJKejpSU+HiQoGoOnWQkoLv3ykQxTDJyXBwgLExAgOhpUWZ2PXrYWmJlBRYWMDeHgkJcHVFx44YM4b6wJmcHLi6YswY9ucfKvJgJZW+WGbPxvPnZfAlXgZq1cL375Cs0bu80OsIvb2RkYGdO2lVUjz29oiIoOuNnhXOnYOGBjV7hBwOrKyULK3+58+fb96kd+yIDh1w4ADF+QYcDvz8oKeHb98wdizatkVEBLZuxadPry0sxixbtoFCXcOHg8vF3r0UilQhL4rjCDU0sG4dZswAn8+2KWwjyqB4/ZoRXfSJzsrCmjVYsACamvQpKREtLdjZQaGa7MhJUBBq16ZMmnJtE967d69eve7m5u369Hm9ejU4NJRV5/Fw5Ah+/sTLl1i5Eo6OCAtDrVqb4+L6rV9/5f379/KriI2Nbd16SGDgP0ePgstCAzQVJcJ6pExB+vVD7drCOXPufSnD1fElg7HVURpvR3d3aGtjyRL6NIjB0bHsJFHk5OTcvfuwefPimyrIgHJtE8bGxv38aaeu3rBLl3f0aRHVWnv8GDExOHYMzs6oW7eHoeEqgUDDx6e6/G/onTqNevCgIyETW7dOEj9aBYNERFAThkYVtWuv2bJlt41N95SUFLZtYRPG4mXocoSJiTh+HBs30vLyLiGOjrh6lTXt1GJh0f79+73BwfZUCVQuR6ijM7xatTq+vg5dunShVZGuLoKDceUKHj7E9es4csQAqODgUC8xUa1rV3z7Jrvky5eRmOgEnFVTS61ESZCPCopITER2NkxM2LajAFpaGTyeYXo6L6dIi6ZyBmMZFHQ5wiFDYGoKV1eaxEuEpSVychhaYqabDx++ArYpKZQtldjYIDpaOVItBQIsW6Zx8OC0yZPHcuh/sdLTQ1AQ/Pxw/z5at/ZPTFwWHPxmwYL3XbqgZUs8fiyLzPnz0acP+vVbFhq69Nevl5R34lQhDwq1LipizZpFPj5t9fWP6ulVYdsWNlFuRxgZidu3cegQHbKlo1u3MhJ/1b37Bg7nrbc3ZfWfKlaEsTGjjU5kZv9+GBmVWEeUDmrUQEgIVqxAq1ZjGjZc3qBBHScnUwArV6JnT5w9K4WozEy0bQtfX+zdizNnYG9vX6FCBbrsViETCugINTU1588fbGjYqCxFOciAqE89A9DiCIcOha0tOnSgQ7Z0lJkkivT0/q1arfLymkmhTKWIl8nMhI8PVq9mWq+ZGa5cwZYttRs27ODl1eW//7gPHsDHB8uXY+ZMzJ+Pkhr/FuThQxgZ4e1bvHwJNzf6jVYhE8+fK5wjFOHuTmMPc6WgVi0kJTGRQUG9I7x+HTExOHqUcsGy4OCA27fLQiDy8+fo3JlimUpRaG3LFrRpgzZtWFBtZYXWrVefO2c4duyG1NRPly9j3TqsXo2WLXHnDoYOFVMW0tcXbduieXMkJKBuXaaMViE9CjgjFDFiBK5dQ2Ii23awB5eL2rWZ6EFBvSMcPRpOTjA3p1ywLBgYoH59/Psv23bIR04OkpIwbBjFYhV/RpicjPXrsWIFawb06dNUX/+kmhq/Sxf9wYOhoYHnz9G4MV69wvv3aNu2+G5fogo1c+di9Wrcvq1qsaTQ5OTg5UtGC0BKTuXK6NsXR46wbQerMLM6SrEjDAw0SEzE4cPUSpWLMrA6GhwMHo/6l1YbGzx7RrFMalm1CoMGsflSNX78qFevzv34cT8+XrtvX6xfj4YNkZmJgABUroyfP9GiBW7fLnTJy5cwMkJYGMLDVbWylIDYWNSqBYXdt3V3h59fue5TyEy8DMWOcOPG6uPHQ1+fWqlyUQayCQMDaYntNjPDz59ITqZeMiV8+oT9+7FoEctmVKlShcfj6epi9GiEhOR+l5ydIRRi+HDweOjZM38L89AhWFrCxASfP6Np06Kinj2L6Nx56Lp1bFRaUlECCrsuKsLeHjwe7t9n2w72UD5H+Pp15/T0uVu2UCiSAuzs8PKlUtbVzCMsDK1aUS9WwQutLVmCiRMVK7sLQOPGWL0aHz5gwgSEh0NNDZaW8PKK0dCor61tN3bsBy8vPH0KLS38/o3oaFy9in374O2NsWPRpcua27cn+/gcSE1NZfv3UJGLwkbK5DFuXLkOmVG+pVGhcAwhV5KTFWtvV10d9va4cYNtO+Tg3Tv060eLZIVNq4+JwcWLmDOHbTtKQEsLzs4ICcHly2jdGsAQYCAwlscb899/sLKCnh5q1ICzMzZuxP374HBgb49x4/rp68/g8+v5+emw/RuoyEXBZ4QAxozB+fOKu3JDN69eXbp71+ngwZO0aqF2H38HkP79e2U5W8RRjmibkJK+DcwTH4/sbBod4ZMntEiWk0WLMHcu9PTYtkMcTZpg+3ZoaTlt2HAU4NSr12XKFNSuDROTYo13WbvWJSEBvXohPh4bNoDHY95kFYVQtOJqf2JggO7dcfw4Jk9m2xQ2WLvWJzv7woIFPVxdh9CnhcoZYe3a9zw9j3TqpFEkfIB1lLrWmr8/KldGxYq0CFfMGeGDB3j4UJlue1/f1Rs2zFy0yDUqan+vXrkzwpIwMcG9e4iLw6BBSE9n0EoVf/D9O9LTYWrKth3iGD8eu3ezbQRLjB07XF29d9euA2nVQqUj1NBQd3ev5O+PIUMUoqxMHubmUFdHTAzbdsjEtWto1Igu4TY2iIqSKDecSebPh48P+11SpcLT09PLayFXsqYSOjq4cAHVqqFzZ7nql6qQE9G6KIv1kCWkSxekp+PRI7btYIPFi6cNHfqwa9fFtGqhPo+wWzfcvYuVKzFtmgI9YR0clDV2NCoK3brRJVxXF4aGDFXzk5CgIHz+jJEj2baDZtTUsHs3Bg5E27Z48YJta8orir9BKILDwdix5Tdkxtqa9pg+WkqsNWiA+/fx7BlcXBRl8UdJswn5fPz4QX0qfUEUKq1eKISXF9auLS9J6PPmYfFidOyI0FC2TSmXKE4/XrG4ueHUKfz+zbYdbGBjo5yOEICBAYKDUbEi2rXDhw80KZGCrl0RGoqsLLbtkJJLl6CmRuPSKBRsm/DoUWhqondvtu1gEFdXHD2KwYNx4gTbppQ/FD9SJg9jY3TsiIAAtu1gA2tr2kt/0NiYV1MThw7B1RV2dggPp0+PROjpwdIS9+6xbIa0nD2LWrXoVaE4jjA7G97eoKkBvSLTtStu3MC8efD2ZtuU8kRODmJjFbS4WrGMH19OV0dr1AAh9O6m0+gIRUybhk2b0LMnzp+nW5UYlHF1NCxMlKZGI4pTaG3HDlhZoVMntu1gAysrhIXh/HnF2lkv27x8CRMTuuKx6aBHD3z5gqdP2baDDays6H1fp90RAhg0CNeuYepUll94HR2Vrzfh+/cYMIBeFfXq4ft3/PxJrxaxpKZizRosX86yGSxSowbu3MHLl2jT5ryzs2dMTCzbFpVxlCVSJg8uF66u2L+fbTvYgO54GSYcIYCmTREWhosX4e7OWlOk1q0RH48vlPV4p53Xr8Hn075hxuXC0hKRkfRqEcvatejRQ2k2bGhCVxdnzggiIhYHBvYdN24h2+aUcRS/uNqfuLvj2DFFiUBkkjLiCAHUrIl//kFSErp0QVISY2rz4fHQubMy1Vrz94e+PrS0aFfE+jbht2/YtQvLlrFpg4Kgrc1r0cJMS8v79+9O5bnnAAMo3YwQgIkJWrXC2bNs28E4dAeOMucIAejo4PRptGtH6tQZU71668DAC0xqh7JlE4aEMLSTz7ojXLYMo0fTHhakLNy9e+7LlxtVq/7l5cW2KWUaJQoZLUj5bFtvZYWYGAgEdMkX7wiTk5OPHz++dOnSrVu3fpF7YZHHw4IFKRxO3Jcvm/fsCZRTmrR0745r15SmuVd0NI2p9AVhMV4mLS1tx47AEyfiFyxgxwDFpHJl7cBABATAz49tU8ooP34gJQW1a7Nth/T07Yu4OGWtkyUzFSvC2BivX9MlX7wjHDp06MmTJ3k83oMHDywsLKKjo+VUWbly5WnT+uvp/W1tPV1OUdJSpw50dRW38VBBsrLw8ye9qfR52NggMpKdYMXRo6f/9dczPr+vnh5tL3vKiYEBLl7E4sW4eZNtU8oiylJc7U/U1DBqFA4cYNsOxrG2pnHhSrwjDAwMPHfu3JIlS44cOeLg4LBv3z75ta5YMTc09KK/f/OUFPmFSYdoUqj4nD8PdXU0aMCErsqVYWCAN2+Y0FWE7GwtoTBZV5fDUcZnEs00bIiAAIwYgZcv2TalzKFENWX+ZMIEHDqkfOVB5ITWeBnxjlBXVzfvZw6Ho6mpSYliKys4OmLjRkqESYGybBOeO8foug1bhda0tDZMnDgwIuKGhOWqyxsdO+Lvv+HkhETF6vKp9Ci1IzQzg6UlLl5k2w5modURSlHSMTQ0NCQkZP369SUNyMjIWLVqlWHhboT9+vWzt7cvdvyiRRw7O42xY7MNDZnbtWvbFiNHav34kVmhAl0qMjMz1dXV5RTy77+abdsKMzNlyTURCARZWVlSuZbGjdUeP4aTU47YkTk5OXw+n5IJXHQ05+5djagouwoVkJmZKb9AFuHz+QJ6tvKHD0dkpNqAAdygoGyK3kJlJDMzk8fj0fRrUoLkt97TpxpjxuRkZkqxH0DJfU0Vrq683bt5vXtnSzJYoSz/E0KI6KtV+jALC05EhEZmpnQTYaFQKFYyJHeEkZGRLi4uBw8erF3yPIXH49WsWbNGjRoFT+rp6ZVkh5kZhg0jvr7q69Yxd2vp6aFZMxIWpuboSJf35fF4kvzpSychgTNwIJFZjrQ2NGnCOXGCI8klhBAJv1ti+ftv3qxZQl3dstCdVigUAqDkz/Inq1aRIUMwdarGvn1sOiHe/2HRhtKR0DyBAC9ecK2tJfrCSyucGQYNwpw53Ph4NTMz8c8xhbL8Twghklhobo4vXzgZGTwdHSmEcziSbbsQCYiNja1Ro8axY8dKH2Ztbf306VNJBObx7RsxMCBv3kh1kbysWEFmzqRRfkpKipwSoqMJh0OysmS8PCcnJy0tTapLYmOJmZlEI/l8fnp6uixmFSYyktSoQaQ0U3HJzs7OyMigT356OmnViqxcSZ8G8WRkZGRnZ7NpgTgkvPViYkj9+nQJZ4zp04mXl0QjFc3yIgiFwtTUVElG2tqSf/+VTrhAIJDkSyt+9ezVq1fdunXz9vYeRkMIo6EhJk+Gjw/lgktD8eNljh5FlSrQ0GBOY4MG+PYNTMYuLV2KOXNA3wJ1GUNbG+fOYdcuHD/OtiniCAg45+zsERkZxbYhJaKMqfR/MnEi9u1DjvjdDIlIT0+fMmXR7NnL+WyV/pIA+rYJxTtCNze39PT0gIAABwcHBwcHb6oLhs6Zg6tXEcXgXdO8Ob5+RUICcxql5cYNWFszqpHLRePGzBVai4xEWBgmTGBIXdmgenVcvowZM/Dvv2ybUjJCoXDSpCWBgc5OTgvDwxU0Z1epI2XyqFs3KyNjWN26nZ9T4RxOngzYvRubN3+eOzf461f55dECfY5Q/B6hr69vSoGZgpGREbUW6Opi1iwsXowzZ6gVXCJcLrp2xfXrcHVlSKO0xMRg3jymlYrS6u3smNC1dCnmzlVNB6XG0hIHDmDwYNy7p6DJ4FwuV0vLXEtrYe3aI0aNwo8f6NEDPXvCwQFVqrBt3P+JiMDYsWwbITdRUVE5ORoJCZ4HD5729ZX3xdnCogUh2ytVQlzc/EaNYGaGHj3g5ITWraE424s2Nrh0iRbJ4h1hq1ataNFcgClTsGULwsLQti3dqnIRJVEopiNMT0dKCoYOZVovA22gRURG4t9/4e/PhK6yR8+emDULTk64fx+VK7NtzR+8egU+P/Ddu+xq1TQAvHuHa9dw+jQ8PGBmhm7d0Ls37OzAbrKMkhZXK4K1tXWXLmpXr/q1bbtFfmlBQVZubmF+flwulysQ4OlTXL+OuXMREYGOHdGnD3r1Qs2a8uuRCxo79Eq381gqMgTL5LF3L+nQgUJbxJCQQKpWJQIBLcLl3Jo+coRoasplgAzBMoSQ27eJnZ34YfIHy/TvTzZvlkeAIkJ3sEwRxo//aWg4xslpzK9fvxhTKkmwzMiRxMenmPN8PgkNJfPmEVtbYmhInJ3J8uVPGjd2cHObKRQKqbJQklsvOZno6hIZdCpmyMm2baRfPzFjxFqemFhi0OK3byQggIwaRfT1SePGZN48MmqUr4VFx/Pnr8hqclEkD5YhhBgZkU+fpBAuYbCMojjCnBzSuDG5epVCc8TQoMGLAwfC6JAs5w3j7EwaNpTLANkcoYQPCDkd4ZMnpEYNQkXYqWLBsCM8fPgol7uMx1vq73+cMaViHWF0NKlWjYj9+r99S3bsIDVregI3K1VyevfuHVUWSnLrSfjCJ5tw5snMJDVrkvDw0saItXz2bDJlihhFfD65c4csWEDU1GyAD1ZWDlJaWiJSOcKuXaVzE5RFjTIDj4flyzFvHkMVL1+8eJGQMMrDY83hw8eY0CcNDx8yt0RcED096Ovj7Vt6tXh7Y8ECaGvTq6XM07lzhwYNbmpr/3PoUPuMDLat+T/e3pg1CwVKURVPnTrw8MCJE8Nq1lycnW2kq2vCiHW5lI110Tw0NTFrFlavll3Cly/Yvx9iq96rqaF9e6xciYkT+1asOLhWrfGyq5QDmuJlFMURAhg0CNraOHWKCV05OTkaGhrZ2RWysxUuVjghAc7O7Kimux/Tkyd4+BDjxtGoopxgYmISG3s7OfmWqWnNHj3w6xfbBgGRkfjnH3h4SDre3t4uIeHusGEHVq5kNBgjIoLpkGy6mTgRoaGyh3yvWgVXVyn2/7Zt84mP/zc83JmVljVl3xECWL0aixcz0cLe0tLyxo3N1auPsbEZTbsyaXj2DAIBQ92X/oTufkze3li4UDUdpAw1Nezdi1atYG+Pjx9ZNsbLCwsWQKqqHwDWrcPRo4x2AStjM0IAFSpg2jQZJ4WfPuHYMcydK91VBgbw9oanJwvpMTTF9CmWI+zYEXXqYP9+JnTZ2jb39HTcv1+xmh4cP46qVcFWXUBaZ4SPHyM8vCyErSsUHA7WrcPIkWjfHnFxrJkRHo6HD2VJDDUwwLJlmDKFoUeqUIjoaFhZMaGLSTw9cf26LF1Kli+HuzuqVZP6wgkTIBDg0CGpL5QTS0u8eEFZGYE8FMsRAli1Cj4+SE9nQperKwIDkZbGhC4JuXmTzVRfWh3h0qWq6SBdzJuHJUvQuTOePmXHgMWL4eUl4z/X3R1CIQ4epNikYomLg5ERKlViQheT6Ohg8mSsWSPdVfHxOHUKM2fKopHLxa5dmD8fSUmyXC4z2towMaH+nU/hHKGtLdq1w9atTOgyNka7dggMZEKXhLx4gR49WNNubo7Pn5GaSr3kx48REaHaHaQRV1ds346ePREayrTqsDBER8s+1+dysXs3Fixg4pFaNmrKFMu0abh4Ubpgt2XLMHUqCrcLkgIbG7i4wMtLxstlho4OvQrnCAGsWAFfX/z4wYSuceNARadhakhNxe/fGDKENQN4PDRqREuhtSVLsGAB2O0iVObp1w/HjmHQIOaKNIlYvBhLlsj1z7W2xpAhWLSIOptKoOxtEOZRuTImTMDatZKOf/UKFy9i2jS5lK5YgaAghIXJJURa6IiXUURH2KABBgxAyX0PqcTJCa9fIyaGCV1iCQyEpiZMTdm0gY54mfBwPH8ONzeKxar4k86dceUKPD2Ze727exfx8Rgtd8yZjw8uX6b9kVqGZ4QAZsxAQADev5dosLc3ZsyAvr5cGitVwurVmDSJ+k27UigvjhDA0qXw82OiLraaGkaPZmh/QiwXLsDMjGUb6AjKWrIECxeqpoMMYWuL0FCszZfaCQAAIABJREFUWiXF5EAevLywdCnUpOjwXTyVKmHNGkycSO8jtWw7QgMDjBsHX1/xI6OjERKCqVMpUDpiBKpWxc6dFIiSEDqeUQrqCGvUwLhx+PtvJnSNH4/Dh5nI2RDLo0ewt2fZBspnhOHhiIpSTQcZpV49hIbi6FFMm4aPHz8JaatSce0aPn+mrC7u8OEwMsKOHdRI+5Nfv5CUxP67Jq3MmgV/f3z6JGaYtzfmzBFf+kBCdu6Ej494pVRhZobERIp7ximoIwSwYAFOn0ZsLO2KzMxgbo6gINoVlQ4h+PQJgwaxbEaTJoiIoDKW3csLCxcy2ltRBYDq1XHzJgICvOvXH9ekSVeatCxdCh8fCqaDeezciRUr6HqkRkTAyorlet90U60aRo7Exo2ljXn+HPfuYfJkypSam2P8eKmTEWVG1DOO2s59ivul0NPD9OlYtowJXYoQMiNq3talC+tm3ExLa9269UBK+nM+eoSYGAXt8lHmMTBA/fqRmZmjY2KS7t+nflJ46RJSUzF4MJUyGzTAhAmYPZtKmSJiY2MHDmz38mXX79+/Uy9dkZg/HwcOIDGxxAFeXpg3j+ImaF5eCAvDjRtUyiwFygNHFdcRApgxA/fu4fFj2hU5OyMsDB8+0K6oFE6cgKEhlS/XsnHkyAWBYElsbGYCFTu0Xl7w8lJNB1nj2LGNs2a9XLFi59Ch3IkT8fs3ZZIJwfLl8PGhfoK1aBH++4/6R+q1aze/f3dJT2/86NEjikUrGNWrw9kZmzcX/2l4OMLDMZ7qQqHa2ti2DR4eyMqiWHKxUB4vo9COUFsbc+bw3dyOX7sWQrciZ2ccOUKrEjHcvq0Qgd2LFnk0bnxAV7dNnTp15BQVFobYWAriCVXIjKmp6fr1S+fPtxe9PjdqRFlmxdmzEArRrx810gqirY3t26l/pHbv7qymFt67Nzp16kSlXIVkwQLs3o3k5GI+Er2b0lHXomdPNGokUaiO/FAeL6PQjhCAUHjg+fM7gwevffLkCa2Kxo3D3r0slM7L48UL9OzJmvY8LCwsIiIC1dWXPH0qV/G5YcM8O3du0abNHtV0UBHQ08Pu3Th+HIsWoU8feQuTCoW500EOPQUKe/RA48YUJ1CdPWs4ZszhU6e2apaD8OVatdC3bzFlSe7fR0wMjWUOt27Fxo20d7DB/5dGKXxcK7ojNDExqlgxLiPje6VK9HbjbtEClSvjn39oVVIiv34hLY3NVPqC8HgYNw579sgl5MqVO1lZR6OjAygySgUFtG+PJ09ga4vmzeHnJ/tzJCAAGhr0lkDasgWbN+PNG2qkCYXw88OkSdRIUwoWLsSWLUXbkixejKVLadyqqFULM2diyhS65OdRtSq0taksNK/ojnDQoP7h4TubNbvy4AHtUc9jx7IWMnPyJLS0UL06O9r/ZPx4nDwpe4Bydja0tWdZWMzbunUppXapkBctLXh7IyQE+/ahUye8eCG1BIEAy5dj1Sq6poMiqH2kXr0KQ0PY2lIjTSmoVw89ehRK77t7F+/fY9QoevXOmoV373DhAr1aQPU2oaI7QgDm5g02b662cCEyM+lVNGoULl8ufmGdbi5dQv36LOgtCWNjdOmCY7I2Ld60CU2ajImNPdexI9t5kSqKw8YGYWEYMQLt28PbG9nZSJO49ry/PwwM0JWujIx8Zs7E+/c4d44CUTt3lq/poIhFi7BxY37pYC8veHvTHo6noYFduzB1Ku3NDKgNHFUCRwigbVs0bYpt2+jVoqeH7t1lf/rLw+PHaN+eBb2lMGmSjNUivnzB+vXYsoVqg1RQCpeLCRPw4AHCwoR6el1NTBw3bBC/Gs7nw8cHK1YwYGDuI3XaNHmrwL9/j7AwRdl3YJJGjdCxI/z8ACA4mMrSB6XToQPat6e9HEq5mxGKWLMGa9eC7hSgceNyvzdMQgi+fGGtK31JdOmCrCzcvy/1hdOnY9IkmJvTYJMKqqlTB4GBaTzez58/Z86de8/ODn/9hSNHEBuLYsvRHD7MNTNDx44Mmde+PTp0kPeRuns3Ro2iOG1OWVi6FL6+yMjA0qVYsQI8HkN6fX2xbx+9NZypDRxVGkdobg4XF9pfRbt0QVoaaA5QLUpYGAhRuBkhh4OJE6WeFIaGIiwM8+bRY5MKGtDV1fXzWzhyZGhEhPfatahXDyEhGDQI+vqwt8e0aTh8GFFRyMrKnjp13oIF55cvZ9S8NWsEGzdOtrBwjJSpKwqfj4MHqU+bUxYsLWFu/qR16/nfvr1gsmpVtWpwdX3WsqXD8OFTaarw17gxXr1CdjZF4gh1WFtbP336lEKBRfj2jVStSl6+pE8DIYT4+JApU+SSkJKSItX4v/4i1avLpbEIOTk5aWlp8stJTib6+uTbt0In+Xx+enp6seP5fGJtTc6elV+z8pGdnZ2RkcG2FVTy/Tu5epWsWEH69SM1ahBt7cnASKDntm3bmDQjOjpaR8cZuDRu3Fyxg/+89U6cIF26UGOJtPe1gqCpWRfYw+PVZ1jv6NEzgGAdnf5xcXGljxQKhampqTKoaNSIPH8uZoxAIMjOzhYrSmlmhAAMDTFjBu1Ny9zccPw40tPp1VKQO3fQrBlz6iRHTw/9+uHQIUnHb94MY2P070+nTSqYokoVdO+ORYtw7hw+fsT06bU5nEguN65WrVpMmlG/fv127TR0ddd++uQiQ77Hrl3w8KDBLOVBR6cCl3tTU1PAsF4Pj6G1a/vw+ZVev65Dkwoq42Vk8MMlQfeMkBCSnk5q1SJ379KqhDg5EX9/2S+X9s2xYkWyfbvs6v6EqhkhIeS//0i9ekQgyD9T0ozw82diaEhiYylRq3yUvRnhn4SEhDx8+JAV1RkZpE0b8vffYoYVufViYkj16kSC+YBEKOmMMCsra8+ePVQ9EKTl7l1iaEgiI0sbI/OM0MeHzJ8vZkwZnBEC0NbG8uWYPZveEjBM1uD+/h3p6QoXKZNHq1bQ08P16+JHzpqF8eNhYUG/TSpYwt7evglLZQC1tHD2LHbvxsmTUly1cyfc3aGuTptZyoCGhsaQIUMqsBQs1K4d1q1D375ISqJeOIWBo0rmCAGMGgU+H2fP0qiib1/ExuLVKxpV5HHyJCpUgKEhE7pkQ5I8irt3cfcuFi5kxCAV5RJjY5w/D09PPHwo0fiMDBw7hnHjaDZLhTjGjMHAgRg0iLrAlv9DYeCo8jlCLhfr1mHOHOr/rHmoqWHECBw4QJf8ggQFoUEDJhTJzPDhuHsX8fElDsjJgacnNmxAxYoMmqWi/NG0Kfz8MGCARLW1jh1Du3aoXZt+s1SIY80a6OlR2QFRRJ06+PkTP39SIEr5HCGAzp1hbo7du2lU4e6OAweQk0OjChFPnqBDB9q1yEOFChg2DPv3lzhg2zYYGrLfUlhFeWDAAEyahH79xIezqcJkFAcuF/7+ePiQ4qIoHA4sLSFTWk1RlNIRAtiwAStW4McPuuRbWKBuXVy9Spd8EUIhvn6Fiwu9WuTHwwN79qDYTr1fv2LlymLq3KtQQROLFqFhQ7i6lhYo8PAhfvyAgwODZqkoFV1dXLiAlSsRFESlWKoCR5XVETZqhD59sHYtjSoYCJkJDQWHg7Zt6dUiP40awdy8+EK6c+di7Fg0bMi4TSrKKxwO9u7Fhw+lVZzZuRMTJ1LfNFiFPNSujVOn4OaGqCjKZFIVL6PE35Tly7F3b2l7V3IyZAju3MHnz3TJDw8P7969KdAuOZnmwnFUUGzIzP37uHULXl5sGKSiHCMKIt27FydOFPPpz584dw6urkxbpUIs7dph/Xoqg0hVjhA1amDKFCxeTJf8ihUxcCCNbev//ts3K8tNILA4Qp8O6hg4EFFRhbr2CASYMgXr10NHhz2zVJRXjI1x7hymTi0miPTgQTg5wciIDbNUiGP0aAwahIEDqYl2FAWOyp9Np8SOEMCcObhxA+HhdMkXrY7SlLNoYzMeeKClda937960KKAUDQ24uRWqSL5zJypVUtwMSBVlnrwg0oSEQuf37FGFySg0q1ejShVq/kcVK2bzeEdPnrwtpxzldoQ6OliyBLNn0yW/TRtoaeHuXVqEHzzYuV07vx8/ntZXqFaEJTNpEg4fzo3W+/YNy5Zh2zZ6u7OqUFE6AwZg8uRCQaS3bgFAu3YsGqVCDKIg0kePKGjWtmnTrpSUfydM8I6Rr9WFcjtCAO7uSEykOBKpIK6utITMBAYiIQFHjlTU1tamXjo91KqF1q0RGMgBMH8+xoyBtTXbNqko9yxciMaNMWZM7srNzp3w9GTbJhXi0NHBhQtYvVreR7exsaGmZhyf/0tXV1ceOUrvCHk8rF6N2bPpyvkbPRrnz1OTs1mQqVPRvTvq1qVYLN14eGD3bm54ODc4mMbdWRUqpGLPHnz8CB8ffP3KuXEDI0awbZAKCahdG4GBcHOTKxFw9Ohht29v0tK6pqlpIo8xSu8IAfTuDROT0jK+5cHAAA4O0lU4FMv+/UhKwt69VMpkho4d0x8/bt+p09ApU2IqV2bbGhUqAABaWrhwAZs3b2nUyFFff22lSmwbpEIy7Oywbt3/2LvLsCiXNg7g/10akbCxAwyUUAkFAUVRMbALFRBU7O46nqMH6xjYYoKCYqAYWJioKNiBCYqKgYEgiMTuvh/WVxFh86ld5vfhXJe7zzNzK4e995m5Z0Zoa+thamp/6tQpxRqxs2vYrVvF0FClIlGHRAhg2TL89Re+fqWlcXf3xGnTBq9atZmS1oRCTJ2Knj1hakpJe4wKCdkuENQSCDqdOUP2FSU4pEIFZGWtKijYlZy8nu1YCDlUr34uPz/n69dFo0eXvCxUGvER4spUNapJIrSxQbt2WL6clsYPHFicmekzf/7mzMxM5VtbvRqZmVivmr+tnp6eurpXNTUDR4zwYjsWgvhNp05t+Hy3unVrsB0IIQcHB4eyZd/zeHMNDCYqnMlatkTZsj/qpBSjJokQwKJFWLsWKSnUnz85YEBnY+OZubk1AKXmYwEIBPjrL3h7o3x5SkJjWo0aNTIzH2dmJvYhyyYIjjl4cOuXL/eePYtlOxBCDgYGBl++PExLu6Kr213CVkFSDR2q1O7T6pMIq1dHjRozGzVyHDFiJrUt+/j0+/Tp2sCBByZNUnatQGAgvn/HqlWUxEUQBKEO9PRw5AhCQrBrl4IteHvj9GnFNwJTn0QI4MOH0zk5IZGRpylvmc/nBwUhNhb79ineSF4eFi9GQACUK/QlCIJQNxUqICoKkyfj/HlFbi9bFj17IiREwd7VKhFu3brE3n5hfv4SOs7ULVMGYWEYP17xLx2zZ0MoxOLFlIZFEAShFiwsEBEBLy88farI7QEBCA6GUKjIvWqVCDt0aHvt2q5//23r6Qkq6lqKsrXFsGFSzn8pSU4O1q7FpEnQ16c+MIIgCDXQujUWLICHhyK7ctvZoVw5nFZoQFCtEqHYqFFwclIwXUk1d66CNZ8TJ0JDA3/9RX1IBEEQasPfHz17omdP5ObKfe/w4QqWzKhhIgSwdi1SU2k5rVBTE7t2Yf58+XZDyMjAtm2YMwfa2tSHRBAEoU6WLEG1aoo8zHh54cIFpKbK3aN6JkIdHRw6hLVrcfw49Y3Xq4eFC+HtLccxIqNGQU+Pxs3BCYIg1AaPh23b8OKFpLOXi2VggL59FdllTD0TIQBTU0REYMgQ0FE4ExCA6tUxf75MF3/8iL17ERgITU3qIyEIglA/eno4fFiRBRUjR2LzZgjkXE+utokQgKMj5s5Fz57Izqa+8a1bERIiU6Wvvz+MjMgBaQRBEHKoWBFRURg79lzz5l4HD8p6SoWVFapWxYkT8vWlzokQwOjRcHCAtzf1hTMVK2LrVvj4ID1d0mWpqTh2DCtWgK/m/9IEQRAUs7CAru6Mmzf/GjZsjux3BQTIXTKj/h/Pa9fi9WssW0Z9yx07wsMDEydKusbXF5UqYfBg6nsnCIJQe4MGeZqYDM7JaTd+PPLzZbqlf3/ExSElRY5e1D8Rigtn1qyhpXBmxQpcvVriIU1JSTh3jhzjThAEoaBly2Z//Hj1xYtlT57A3V2milA9PQwYIF/JjPonQhQqnElKorhlfX2EhWHsWLx6Vcy7gwejZk307ElxpwRBEKUHn8+vWBHR0ejRA3Z2iI6Wfou4ZEbGJ0iUkkQIwNERc+bQUjjTvDnGjIGfX9FpyAcPEB+v1IboBEEQhBiPh/HjcfAgxoyB1GHSRo1Qrx6OyVphU2oSIYAxY2BnR0vhzJw5yMtDUNBvLw4ahPr14e5OcV8EQRClloMD4uORlAQnJzx/LulKuUpmSlEixP8LZ/77j+Jm+Xzs2IHAQNy9++OVa9dw754i6zoJgiAICSpUwJEjGDgQDg7Yvx+CEtYM9umDmzdlXUdeuhKhri4OHcLq1dQXztSpg0WL4OWF798BwMcH1tZo0YLiXgiCIAjxMOnhwxg79pKRUcsGDVwyMjKKXKOjg0GDsH27TJWKpSsRAjA1xZ49tBTO+PvDwgLz5uHsWc1nz8jjIEEQBI1atICX18Vv34Y8fWrq6vp07lycOYOcnF8XDB+OHTt4suyFWeoSIQAnJ0yfnm9vP6V9e+93795R2PL69Vixwq57d5fKlcOsrSlsmCAIgihq+vShffrcmzTJes2a5gYGWLwYFSvC1hYzZiAmBqam+d++9d+y5bbUdkrp9pe2tlezsr6cPu04c2bY9u2TqWo2J+eVQJABRL596wkMpKpZgiAI4k+VKlWKiPhxKp6zM6ZPR3Y2rlzB+fP46y/cunUpJyfzw4f6UtspjU+EAKysrBo0eFGhQsjly+2aNUNEhNybtBaRn4/du9GtW2UgG+hZqZIuRZESBEEQsipTBu7u+PdfXL6Me/caaWk9NjG5J/WuUpoIjYyM7t6N+fAh7skT69WrsWULzM0RFPTb+LKMvn5FUBDMzLBuHWbO1E5KurJq1ejXr+9Kv5MgCIKgTb16VV68iPX1bST1ylKaCAtr1QqnTyM0FDExqFMH8+fjyxeZbnzxAjNmwMwMly9j715cuoQ+fVC3bi0/Pz9NcuQSQRAE26pUqWJoaCj1MpIIf2jVCkeO4PRpJCejXj2MH483b0q8+OZNeHujeXPk5CAhAXv3wsGBwVgJgiAI6pBE+BtLS4SG4sYNAGjcGN7eePr017tCIY4cgbs7evdG48ZISkJQEGrWZCtYgiAIggJkBK8YtWsjKAizZiEoCI6OsLKKefx4ftWqzT5/Xl2xIiZPRo8e0NBgO0qCIAiCCuSJsESVKyMwEElJePVqR2rqurt3b2zYkBkXh969SRYkCIJQHyQRSmFoiODggIYNxw8Z4uruLn3SlSAIglAtZGhUutatnR8+PM92FARBEAQtyBMhQRAEUaqRREgQBEGUaiQREgRBEKUaSYQEQRBEqUYSIUEQBFGqkURIEARBlGokERIEQRClGkmEBEEQRKlGEiFBEARRqsm0s8yTJ08uXbpUtWrV9u3b8/kkdxIEQRDqQ3pWO3ToUMuWLRMSEmbPnt27d28GYiIIgiAIxkh/Ipw9e/aaNWu8vLyysrLMzMyuXLni6OjIQGQEQRAEwQApT4TPnz9//Phx9+7dARgYGHTo0OHYsWOMBEYQBEEQTJDyRJiammpiYqKvry/+Y7Vq1d68eVPSxXl5eaGhoVWrVi38oqurq7W1tfKBqpD8/Pz8/HwWAxAIBPTFUFBQkJ+fr6lJzi35TX5+vkAg0FDrkyrz8/OFQiHbUUhC668e67/XCuN45CKRiL4IhUIhj8eTepmUjzOBQFC4FT6fX1BQIOHi1NTUvLy8wi9aW1sLBAKpcagTgUDA7l9Z8H8q17jqKg3/LOK/HZfL5Wj9Eajuz5fjkYtEIvoiFAqFsnw9lZIIq1atmp6enp+fr6WlBeD9+/dFHvgK09PTmzlzZml7/vtTfn6+rq4uiwGIv77QFENBQYGGhga7f0EO0tDQEAgEav/PoqGhIf4o4CZaf/VY/71WGMcjF4lEQqGQpgiFQqEsKVbKl7u6detWrVo1JiYGQEFBwZkzZ9q0aUNNgARBEATBAVKeCDU0NGbNmjVs2LDx48fHxsZWqlSpffv2zERGEARBEAyQPtwfEBAQEhKSkZHh4eFx5swZLs8QEARBEIS8ZKr9a9u2bdu2bekOhSAIgiCYRx7vCIICQqGwSL00QRCqgspEmJxc186uG4UNEoRKWLduna5uXSOjJj4+PmzHQhCE3KhcFl1QEJif3+vt27empqYUNksQXPPpE+7dw8OHP/4bG7sP6ADUCQ9/vno1jIzYjo8gCHlQmQh5vCBg6uXLpmRrbkI93Llzd8iQ6VZWVmPGLE5M5CUm4sEDJCbi82fUqwcLCzRuDA8PODquDwzsCJQzMjpavjy6dcPOnfj/dkwEQXAdTyQSUdWWpWWz5OSE/HyNly9RpQpVraqer1+/li1blsUABAJBbm6uPj2fxOIt1vT09OhonGtcXcdevNiDx/uvadMNzZvXEme+xo1ReFeJV69gbY0yZUQVK4pWrOC/eoXJk/HlC4YMwbp1UKet6L5//87xBfW0/uqx/nutMI5HLhKJvn37VqZMGToaFy+ol/o/LZVzhDxewejR73V14eREYasEwYLr19G6NV696lu58pwOHcrHx1cPDsaECXB3/y0LFhSgf380aQI/P2GPHsIDBzB4MNLS8N9/CA+HoSH++ou9vwNBELKhuGq0Z89POjpITYW/P7UNEwRDXr9GQAC6dEGPHnj61PnduyvHj+8sabvCv/+GgQGePkW/fqIePQQHDkC8K/W4ccjMxJQpWLIEhoZYvZrRvwJBEHKhOBHq6QlHjkTbttixA4cOUds2QdArKwvz58PGBiYmePoU48dD8m69sbHYuhVDhqBWLZibi8zMROXL49q1H+/yePjnH2RmwtsbU6agfHmEhOQsXbr07t27DPxdCIKQHfXrCMePx7Vr6NcP/fsjLY3y5gmCevn5CA6GuTmSk3HvHhYvhtQplfR0DB6M4GBER8PL68eLvXrhwIHfLtPWxtq1+PgRrVvD13fG9OnPmzbtlJWVRctfgyAIhVCfCMuXR//+qFMH1avD1ZXy5gmCYjExaNoU+/bhxAmEhkLGtT8jRqB3b7Rti2PH0LfvjxfFifDP+jNDQxw4gNq1bwJlhELDtDROn+pHEKUNLTvLTJ6MTZtw6hSeP8eoUXT0QBBKSU1N3bhxc3R0qqsrxo/HkiU4fRqyHyC2cSOePUNgIKKiYGf3q0ba0hI6Orh1q/i7bt06Mn26Vq1a++3sDDMzKfhbEARBCVoSYZ06aN8eBw9i584fGZEgOMXNbcCoUYLu3b2GDMG9e+jcWY57ExMxbx7CwqCtjfBwDBz427s9ehQdHf3J2Nh48eJFT55Y6OujUSOQHdkIgiPo2mt0xgysWoVu3dC3L7p3R0YGTf0QhNzy8pCWZqitfb9ZM0NfX8h1nkpuLry8sGQJGjbE58+4eBHdft9VsGdP7N8vqQVtbTx8iNxcNG2qSPAEQVCOrkRoZYXGjREWht27UaECXFxo6ocg5CMQwNsbLi6RZ84Mio2NlPf2yZNhbo4hQwBg7154eMDQ8LcLbG2Rm4vEREmNGBjg3j2kpMDZWd7+CYKgHo2nT0yfjiVLIBQiLg4PH2LSJPq6IgiZiEQYMQKfPmHvXm0npxbybpJy/DiOHMGmTT/+GBb2q170Jx5P0ujoT6amuHkTCQnw9JQrBIIgqEdjImzTBiYmOHIE1aohJARBQTh3jr7eCEK6qVNx/z4OHoSOjtz3vn+PYcOwezfKlQOAly/x8CE6dCjmyj8XURSrfn1cvozjx+HnJ3cwBEFQiN7zCCdPRmAgAAwYAE9PdO2Kr19p7ZAgSjRvHk6fRnQ0DAzkvlcoxKBBGDECjo4/XgkPR58+0NYu5mJHR6Sl4ckT6c02b47oaISGYto0uUMiCIIq9CbCnj2Rno7YWACIjIShIdzcaO2QIIq3ejX27MGpUzAxUeT2pUuRm4uZM3+98me96E98Prp3l3VnJXd37N6N5cvJNmwEwRp6EyGfjylTsHQpAPB4uHIFt29j1ixa+ySIokJCsHw5Tp9G5cqK3H79OlauxM6dv3Zcu3sXGRmSNpeXcXRUrE8f/PcfJk7Erl2KhEcQhJLoTYQAfHxw8ybu3weA2rURHIwlS3D5Mt3dEsQPkZGYOROnTqFWLbnv/fDhw9ix8z09I1ev/u12cZkMj1fija1b48ULvHwpa0cTJ2LmTPj4IDpa7iAJglAS7YlQRwdjx2LZsh9/HDIEHTqgQ4fnf/+9KDs7m+7eiVLu9GmMHInoaDRooMjtM2YsXreuYnr6CmfnNz9fFImwd28x9aKFaWigSxccPChHXwsXYtgwdOuGq1cVCZUgCIXRnggBjBqF6GikpPz4Y1hYZna2x/z5r+3sOjHQO1FqXbmCQYMQGQkbGwVbsLVtCkRUrpxvUmhqMTYWhoawtJRyr1yjo2IbN6J7d7RqFWFq6nzsGHk2JAiGMJEIDQ0xZAhWrfrxRx5PyOdnA9rfvzPRO1E63b6N3r0REaHUMdFVqgxydj6UlHRFT0/v54thYSWWyRTm7o4HD/D2rXw9hoXlCQTz3r2b17//FDmDJQhCQQylokmTsHMnPn4EAGNj49jYCBMT6wYNjjHTO1HaPHgADw8EB6N1a6XaOXAA/fuXK3wqb14eIiN/HTchgZYWOnbE4cPy9aitrW1oCGBmXl4/shkpQTCDoURYpQq6d8eGDT/+6OjouGqVb0yMfn4+M/0TpUVKSsrOnWc6dhSsWIEuXZRqKi8PJ04U3Ur0+HFYWKB2bZlaUGB0FEBGxuPXr08ZGv5lZSX3vQRBKIC5wcmpU7FmDX7Wx3h7Q0fnx3J7gqDE58+f7ey6+/qGW1svGzBA2dZOn4aFBapW/e1FCcsH/9SxI+Ljfwx66rgRAAAgAElEQVSEyKVatXIPHyI1Fe3ayX0vQRDyYi4RNmgAJyfs2PHrlQEDsG4dY/0T6k8oFGZm8jQ0dCwtBcq3duAAevX67ZWvX3HyZNEXJdDXR7t2OHpUkd4rVEBCAi5ehL+/IrcTBCE7RstVZs7E8uUoKPjxx+XL8ekTOa2QoEx8fAVT0z0REZ0XLJiuZFMFBTh6FD16/PZiZCRat0b58nK0o9joqFjDhjh+HCEhZOCEIOjFaCK0t0eNGti378cfDQ1hb082miGokZ2NMWOwdWv9Hj06a2pqKtnauXOoWxc1a/72YrHHTUjWpQtiY6HwefRt22LzZsydi/BwBVsgCEIqphcwTJ+OxYshEv3448qVuHkTaWkMR0Goodmz4eZG2Wa2f46LpqXh+nW5C3DKloWzM44pUR89ZAhmzcLgwThzRvFGCIKQgOlE2KkT+Pxfw6EtWqByZUydynAUhLq5fh179mDxYmpaEwpx+DB69vztxfBweHpCX1/u1pQZHRVbsAADBsDDA48eKdUOQRDFYmFJ++TJP7bhFpswAXv3Mh8FoT4KChAQgJUrUaECNQ1evAhTU9Sr99uLctWLFta9O2JioOR+grt2oVkz2NkpUoNKEIRkLCTC/v2RmHht/Pi16enpAKZNg0iEjRuZD4RQEytXomJFKL9e4qc/x0WfPUNKCtq0UaQ1Y2PY21NQFHblCipXRpMmIAvtCYJaLCTC3Nzsr18D1qzJGzFiFgAeD507UzaoRZQ2KSlYtgzr11PWoEiEqKii46K7dmHAAChcgtOzp7KjowD4fNy9i4ICNG2qbFMEQRTGQiLU0dExMdEAYuvU+VGTt3YtXr7EzZvMx0KovDFjMHUq6talrMG4OBgZoWHD317cvVvBcVGxHj1w7Bi+f1cyNOjr484dPH+u7KY5BEEUxkIi1NTUfPr0cr16i/v2/XHgt6kpGjXCFLLJMCGn8HC8fIkJE6hs889x0YQECIWwtVW8zcqVYWmJs2eVDA0AqlXD5cs4eRIuLhGbN2+moEWCKPXYOf9BV1fX3b3BxYu/Xlm0CBcv4ts3VsIhVNLnz5gyBRs3QkuLymYPHSqaCMPCMGiQpGN4ZaF87ehPTZvCy2tNbOyOgIBdy5cvp6ZRgijFWDsIycUFsbG//ujpibJlMWcOW+EQqmfqVPTti5YtqWzz+nVoaPx21qBAgIgICipxevXC4cO/tlVSko1NAY/3VSTS2r+/EjUtEkQpxloidHVFbOyvlfUA/PywfTtb4RAq5uJFxMRgwQKKmz1wAL17//rj169fx4zZYGwcW7++si1Xr466dXHhgrLtiE2cODEkJGDGjOE3bgy2sYFQSE2zBFE6sZYITU1hZISHD3+98u+/+PqVsuEjQo3l5mLECKxejbJlKW754MHfxkWnTw8MDv706tWU9+/fK984haOjAAYPHrxoUd8HD/D8OWrVwpcvlLVMEKUNm2fEu7ig8DShri5cXTFvHnsBESoiMBCNGxc9KVB5d+8iNxfNmv16pW7dGjzeVQOD3DJlyijffrlyJ4KDHQcMGKN8Uz+Zm+PVK/D5qFmT7DtDEApiMxE6O/82TQhg9Wo8fIiUFJYCIlTB48fYsAGrVlHfsnhctHBRjI/PqDJlVj57dtnAwED59iMidgoEG06cuJGVlaV8az8ZGuL5czRrBisrcpYLQSiC5SfCIlMmjRujZk1MmsRSQATniUQYORJ//41q1ahvfP/+ovWi8fFo0cLcwICCx0EA8+aNMjae1Lx5B0rSamF8Ps6f/7EfKTnjkyDkpexpNcqoWxd8PpKTf1sNPX06xo+HQAANDfYiI7hq82Z8+4aAAOpbfvIEGRlwcPjtxYQE2NtT1oWzs9P8+WeePqWswSJCQtCwIcaNw6NHWLOGrl4IQv2w+UQIwNn5t2lCACNHQksLy5axFBDBYe/eYc4cbNoEPg3/2+7di169ii4WjI+HnR2VvVhb484dKhssYuZM7NuHjRvRujWNvRCEmmE/ERaZJgTQuzeCgtiIhuCw+/fvDx2aPHw4rK1paf/PDWVEIiQkUJwIbWxw9+5vq4Yo17Mn4uORkIAGDcgOFQQhE5YTYZHCUbEVK5CWVkyCJEqtY8eiHR0nnTjh1a3bXTraf/4cb9/Cyem3F5OToa8PU1MqOzI2hpERXrygss0/NW2K5GRkZKB69Y99+gScIiU0BCERy4mwUSNkZuL1699eLF8eTZti2jSWYiK4Jz39a3Z2RX19o4IC5Y71K8G+fejRo+i0dHw8lROEP9E9OipWuTJevEBm5uD9+806dhxWQNWWNgShjlhOhDweWrXCpUtFX1+xAvHx+PyZjZgI7vn0qa+NTf+oqJktqd1R7f/+HBcFqB8XFWMmEQLQ1UXNmunAPZGonLm55p9DLwRBiLGcCFFcvQwAFxeUL4/p09kIiOCYL1+waBEvJKRrmzat6Wj/9WskJxdTXUJtyehPjCVCAE+eXNq82eXBg5ONGqFNG1hY4P59hromCBXCfiIsdpoQwOjRCA9nPBqCewID0a0bmjShq/39++HpWfTQXYEAd+78tssMVZhMhJqamkOHDrWwqBQdjYcPYWQEKyu0bIm3bxkKgCBUAvuJ0Noab94gLa3o63PmID+fbMNd2r14gW3bMH8+jV0UOy567x5q1oShIfXdmZnhwwdkZlLfsmT16yMuDpcuiYto0LMnBQcFE4R6YD8RamigZUtcvlzM6x06YOFCNmIiOGPWLIwbR3HpZmHv3iExEW3bFn2dpkoZAHw+LCxw7x4tjUvl6IjEROzYgfPnYWSESZPIyRUEwYFEiBJWEwJYswbPn+PBA8YDIrghPh4XL2LyZBq7iIxEp07Q0Sn6Ok2VMmJMjo4Wa/BgfP6MBQuwcSOMjODuvsLKqm1CQgKbMREEeziRCEuaJqxdG+bmZOvR0mvKFCxYACoOfihRseOioPOJEBxIhGLTpiEzE927f4yJ2Xnvnq+T01+TJuHUKeTlsR0ZQTCLE4nQzu7HTo9/mjr18+nT/wwcOJLxoAiWRUYiIwPe3jR28fEjbt5Ehw5FX//2DcnJv51TTy2OJEIAmprYvt1YW/szjxdkaGgdEYEuXaCjAwMDWFigXz+sX4/U1B8XL168eNu2bazGS8jty5cvI0aMSCFn+kjE5qbbP2lpwc4OV67Aw6PoW1FRI0WiBuHhL9u23ebn58dGdAQL8vMxcybWrKF37/VDh9C+PfT0ir5+/TosLaGtTVe/VlZ48ABCIS2bpspLU1MzJ+d5SkpKnTp1xK+kpuLIEVy4gDt3EB2NMWOgoQE+f0Be3hMg68WLF4sWLWI3ZkJ25co1EYncw8MdRaJU6VeXVhz4RQRQ8jRhq1Z2PF4k8Kpp06aMB0WwZsMG1KuH9u3p7YWVcVEAZcuiUiU8e0ZjF3Lh8/k/syCAatUwYgR270ZiIr5+RV4eoqNhaJgCmAF1Dx6sVezgDcFBZ86IN7a1A6ru28d2NBzGlURY0jTh9OlTHj06oqWVEB9PEmFp8eULAgOxZAntvRQ7CAGaK2XEuDM6KpWmJtzd8erVWWvrt5qarq1bj6hfH0FBpNyU096+hbc3fH1hablJT29f48b/LVyINm1I7WHxuJIIW7bEnTvFb5Zfv36dXr10/vmH8ZgIlohX0NM3RScWFYW2bVG2bDFvkUT4J11d3atXT3XvPs3GBsePIyICDg6Ij2c7LOIPBQUICoKlJUxMcP063rzpfO7c0dRU1zNnMGAA2rXD+PEsLGPlOK4kQj09WFnh2rXi312/Hu/e4fhxZmMi2MDACnoAGzeGjh7tYmKy68+3PnxAejrMzekNQOUSodiQIcKtW9GsGS5fxtix6NYN3t74+JHtsIj/i41Fs2Y4fBixsQgKwqFDaN0aFhbCdu0QGYnhw3/ssde4MUJD6T0OTLVwJRGi5NFRACYmcHWldz0ZwRF0r6AXW7AgKDv72IkTK/98SzxBWOSEXsqpaCJs21b04QNu3waPB29vPHwIExNYWCAoCAIB28GVbu/ewdsbAwZgyhScOYNGjQAgOBgjRwKAvz+2bgWA8uURFISDB7FuHdq0YW1jB67hUCIsqV5GLDgYjx6p5GcHITsGVtCLubv76eq2GzfO/8+3aNpru4jatZGZiU+faO+IWnw+fH1/bXxobIygIJw6hX374OCAq1dZDa60EgoRGgobG5iY4OHDXyuOrlxBZibc3ACgfXu8e/fr89PWFnFx8PND+/ZkpBTgVCJs1Qrx8SUu5jUzQ5MmGD2a2ZgIZjGwgl6sTJnRCxZcmz591J9vxcfTPkEIgMeDpaVKfh/398fu3b/tU2pjg9hYjBuHnj3h7Y2UlOzr168LSS0NzQoKCnx8JlhZdW/WLHnbNpw5g6Cg3+a8N2zA6NE/xjb4fPj4YMeOX+/y+fD2/lE7Y2GB0FA8e5b08uVLRv8OnMGhRGhoCHNz3LhR4gVr1yIuDm/eMBgTwaCDB2lfQf/TiRPo1KmY10UiJiplxFR0dLRGDdjYICrqtxfFI6WJiTAxQb16rdu0We/jM4GlAEuLu3fvRkZ+vHdvQPXqoefPo3Hj3979+BHHjv322+Tnh7Aw5Ob+dlm5cggKwoEDCAw826iRn41Nn9u3bzMRPcdwKBFC4jSh+N1q1TBmDIMBEUzJz8eMGVi2jN4V9GKPHiE/HxYWxbyVnAx9fdpnKMVUNBGi0IRTEcbGWLlSaGycm53d5MEDVRv2VTU1ajTKy8uqVm3N/Pld/3x361b06IFy5X69Urs2rKxw+HAxTTk44O+/P2toVM/IqPzmzRfaQuYubiVCydOEABYvxpEjxa+yIFQaMyvoxY4fL/5xEPQvpS9MdRNhjx64cwfFbtrF5/MvXdo7bZrh69drnz9nPLLSZNUqPS+vQ69fX7K1bV7kLZEIW7b8KJMprKRvMAD69u0VFta7Q4eRISGtqY+V87iVCF1ccPmypPIzLy+ULYtp0xiMiaAfMyvofzp+vPh19GBkBeFPlpY/nk1VjrY2+vcv8azQhg0bLl48dNYsk379ig7EEVR59QpbtuDff4t/9/hxmJjA1rbo6z174tat4r/B8Hi8Xr16HDzokZSEtWspjpb7uJUIK1RA1aq4e1fSNZMmYetWsquFWlm0CJ6etK+gF8vOxrVrP0rp/sRMyaiYnh5q1MCTJwx1R62hQ7Ftm6TvrBMmoFYtsuSJLlOmYOxYVK1a/LsbNhTzOAhARwd9+yIkpMRmdXRw4AAWLsSVK9TEqSq4lQghbZoQwOzZALBsGTPhEPTKzs7esGH/li0pf//NUI9nz8LOrvgNZQQC3LmDZs0YigSqPDpqaYnKlXHmjKRrxKWMO3cyFVOpEReHuLgSz6d7+RJxcejXr/h3hw3D9u2SHiRq1cKWLRgwoHTtk8C5RCh1mlBcn0YSoXrw9p4wduydvDzPSpUYWo8tYVz03j3UrAlDQ2YCAVQ5EULihJNY2bLYuxeTJyMxkamYSgGRCOPHY+lS6OsXf8HGjfDxKfFdKyuUK4ezZyV10aUL+vSBj08pGnjjXCJ0dcXFi1L2/lm1ChkZ2FXM9liEisnK0hGJvhga8nh0b+XyfydPlpgImayUEVPpRDhgAE6dwocPkq6xtERgIPr2RXY2U2Gpu9BQaGqW+MCXl4ft2xEQIKkFqd9gACxejK9fmZu2Zx3nEmG1aihbFo8eSbpGTw8eHpg7l6mYCHqIREhPXzltWvc7d2L4jBzN9/BhiQsnwGyljJhKJ0IjI3TtivBwKZcNHQo7OwwbxkhM6u7bN/z1F1atKnELwAMHYGmJ+vUlNTJwIE6ckDLyqamJvXuxbh1iYhSPVoVwLhFChmlCAMHBePlSyiAqwXG7doHP1woMbFuhQgVmejx+HJ07l/gu80+E1aujoADv3zPaKYX8/bFli/TLNmzAw4fYvJn+gNTd4sVwcZH0f2lJZTKFGRmhc2fs3i3lsipVsHMnfHxKxR4mXEyEUqcJAVSpAjs7jB/PSEAEDXJyMGcO/vuP9u2tC5MwQfjtG5KTGapcLczSUkqZNJe5uCA3FwkJUi7T1cXevZgzBzdvMhKWmnr9Ghs2YOHCEi9ITERSErp0kd6Uv79M30vatMHw4fDyQkGBHHGqIi4mQhcXXLgg/bJNm3D7NpKS6A+IoMGiRXB2RqtWzPWYnY34eLRpU/y716/D0hLa2szFI6bSo6M8HoYMkT7hBMDcHGvWoF8/kNPtFTZtGsaMQc2aJV6wbh0CAqClJb2p1q2RkyNpP8uf5s6Fri7mzZMjTlXExURoZgahEFK3pbC2hrk5RoxgJCaCUq9fY/16BAYy2umZM7C3L37hBNgYFxVT6UQIYMgQ7N8v02ZPffuiY0cMHkyOwVPE1au4dAlTppR4QVYW9uyBfzHnqRSDx4Ovr0zfYPh8hIdj9+6iu8uqGS4mQsg2TQhg1SqcPYv0dPoDIig1dSrGjpX03ZYOEsZFwUaljJiqJ8IqVdCyJfbvl+niFSvw6RNWFnMKJCGJSIQpU7BokaSDWXbuhJsbqlWTtU1fX0REyPQNplw57NmDgADpDyeqi6OJUJZpQgAeHqhYEePG0R8QQZ24OCnfbWly6pSURMjKE6GFBZ49U+2tyPz9sW2bTFdqaWHPHixbhkuXaI5JvYSFoaAAXl6Srtm8WXqZTGHVqqFFCxw8KNPFDg6YNg1qvGceRxOhjE+EAObPR0SESm7YWDqJRJgwAUuWMHHoYGHihRPiY7v/9OED0tNhZsZoSGI6OqhXDw8fstA1Vbp0wZMnUpY8/VSjBrZuhZeXlAWIxE+ylJVduoScnBLnv0siy4LCnyZORM2aLHx/ZYZMiTA1NfXChQv3798XMTW637gxvnyRqWx3xAjo6pI1hSojNBR8PgYMYLrf6GhJ1XTiCUIm61cLU/XRUU1NDBqE0FBZr+/UCQMHon9/SVuVEj8tWYJWraSUlW3YgBEj5P4f2NMTjx7h2TOZLubxsH07Tp9Wz51MpCdCPz+/pk2bzp8/39PTs2XLll++MHFaFY8HJydZlwmOGIH162kOiKBCVhZmz0ZQEAspR+oEISvjomKqnggBDB2KHTvkKLJfuBB5eS+trEauWLGRzrhUXmoq1q0r8ZQJsY8fcfy4Iidaa2rCy+u3Y+slE++ZN3r0ITc3/1u31Or8XumJcNiwYampqefOnXvy5AmPxwsKCmIgLMg8TQggMBC5udiwgeaACKUtXgw3NxZSTnY2EhIkDRzFx7NTKSOmBomwfn3UrYsTJ2S9XkMDVaosT0x0/fvvva9evaIzNNU2YwZGj0atWpKu2bIFvXrBxESR9ocNQ0iIHI/mlpYioXDuuXND/fxmKtIfV0lPhC1bttTS0gKgqalpYWHx+fNn+qMC5Jkm1NRE79745x+aAyKU8+oVNm6UtByYPuKFEwYGxb8rErFWMipmY6PyiRByTjgB6NvXzcRk+bdvfIGgMm1BqbYbN3DhAqZOlXSNUIjgYMVXkTVogBo1cPKkrNfzeDwnJ2td3QlCYUcFu+QkTdkvffXqVVRU1KFDh0q6oKCg4MyZM8nJyYVftLS0rFu3rgKRWVvj1Sv+hw/C8uWlX7x6NSpV4kdFCbt2VaArigmFQiGr27YL/49TjU+ezBs3DtWri5j/t4mO5nXsCKGw+BnupCTo6/MrV1bqH0yZf/MKFaCpyX/1Sih77TsrhEIhj8cr6e/YuzemTOGnpgpNTWVqrVevrp6eHf/7T9vLi3funFCWZeCyREjfrx7Dv9d5efkTJugsXCjS05P0KxMdzatUide0qaTQJEfu58fbsgUdO8pa/xEdHZqZmd+6tc62bSJfXwqqRkQiEa2fV7JcxhPXv3z58mXw4MF/vj1jxgwnJycAGRkZ7dq169Chw8KSv9KbmZlVrVrV8PdjbHx9fTt2VPC7Q8+eev7++Z07yzTz0LWrXmoq7+ZNGdbF0CwrK8ugpKcPRggEgtzcXP2SDmJRTkFBQX5+vp6enlx3Xbum4eure+NGNj1BSdG4cZmDB3Pq1y/+V2LfPs3DhzV37vyuTBf5+fkCgUBXV1ex27t31xs5Mq9DB05Xj3z//l1DQ0Or5JQ1ZoyOublo/Pg82dsUiTBggG79+sJ//pHjrpLQ+qvH5O/1vn1REyasys83evlyr65uidsdZWdn9+uXM2BAlYEDJX1ISo48O5vXsKH+9evfKleWI6s9e8Zv317v4MEca2tlE5hIJMrJyaHp80ooFOro6Ojo6Ei+7McTob6+/qhRo/58u169egCysrI6d+7s4OCwYMECCW3p6+uvWbPG2tpa0ZiLat0aCQkaJR04UsT27ahbF0+eGDB5sGqxRCIR64lQS0uLO4lQKMTMmVi2DJUqsfDPkpgIPh/NmpX4r3HvHlq2hJI/MiUTYbNmePJEr1cvZUKgnaampuREOGIEfHwwa5a2XMVQO3eieXO4uGh3765shLT+6jH5ex0TczUzc4qR0fq8vO8VKpQr9pq0tDQrqw4fP5r4+o4yMOgjoTXJkRsYoGdPHDpUZvJkOSK0scGaNRgyRP/6dWWP8BSJRHw+vww9C6qEQqFAhinQH4lQW1vbo4Sium/fvnl6ejZo0GD16tWMHRon5uJS4inMf6pdG5UqbbW3D5k2rX1g4Bw64yLkExIi6QQ1ukmuFwWQkACJX/CYYG2NI0dYjkF5LVpASwtXrsDJSY67TEywZw88PWFtjTp1aAtOpXTvPvXo0cX//DO4atWqJV3z+fPnrKxyQMsXL5Td8cXfH35+mDRJvnLufv1w/jyGDUNEhJL9s096sYyvr+/9+/fNzMyWLVu2ZMmSgzJuRUAFe3s8foyvX2W9/sOHfwSCPUuXkuNeOCQrC3PnSjpBjW6SE6FAgDt3wPooghoUjor5+clXMiNmb4+ZM9V54xJ5HT5cc8GC9ePGDZVwTcOGDS0tJw8cWHX69LFKdufoCA0NXL0q942rViEpCWvXKtk/+6QnQmdnZz8/v4yMjPT09PT09KysLAbCEtPWRrNmuHJF1uutrS0Ad01NpUdYCOoEBqJdO9ZW6UldOHHvHmrWVHZsR3kNG+LVK5k2fuQ4Hx9ERSEzU+4bx49HrVpqu3GJXL58QXQ0Bg2ScplIhGfPOi9ZMlreCftiyXiKSBE6OjhwAAsXIi5O+RDYJL1qdOxYZb9uKMPFBbGx6NBBpotv3jz+4oXAzExj0yYEBNAcGSGD589/nJbFlpgYODiUuHAC7B06UYSmJurXx4MHbK7ioET58mjc+Fz37jFbtw6tI+dA57ZtsLfHrl3Sc4B6CwlBp04oV/zM4C/37qFiRchYoyvVgAF55uZratfWmTVrFJ8vx9abtWphyxYMGIAbNyBLhT83cXSv0Z+cnWVdTShWu7bG0KGYOBF5FNSgEcqaNg2TJqFGDdYCkGWCkCO5Rz1GR4VC4Z07486da+7tLU/pBYD/b1wyebJq77yqvO3bMWyY9MsuXICLC2WdnjkTkZv76t9/r5+UfVHh/3Xpgt694eOjwgdscT0ROjri1q0379/LMdSyYQN0dODjQ19QhHQCgSAi4um1a5g4kc0wTpyQkgg58kQIwMpKHRIhn8+vW7cKj7fe1tZWgdstLfHvv+jbVx1GiRUTF4ecHDg7S7/y4kW4ulLWr4VFIwODyzzeXTOF9p5fvBiZmVi6lLJ4GMb1RHj2bPT37/0tLFxfvnwp4y08HnbuxN69ePCA1tAISWxs2g0aNKtOnYmsLBwUe/AAfD4aNizxgm/fkJwMS0sGYyqZejwRArh163SLFns8PGYpdvvQoWjeXKZHIrW0eTOGDZNeWSYSITZWpnwpI1tb2/Pnz1WseNnc3FyB2zU1ER6OVatw5gxlITGJ64nw2bMXIpFtTo5pWlqa7Hd16YKmTdGzJ31xEVKkpHwWCLrm5SWxGIPUcdHr12FpCe0S1yszytoad++q8OBSYW3aVLh8WfHbN2zAgwfYvp26gFRERgYOHZJp++xHj2BgQPGkg42Nwffvegrv/Fq9OsLC4O2Nt2+pjIoZXE+EI0f6+/g0sLAYJu9Iy7FjSEoiO3Gz48EDaGuHTJyYduAAm2cLSE2E3BkXBVC+PAwMkJLCdhxUcHJS6uhdPT3s3YsZM3DrFnUxqYLwcLi7o1Il6VdeuEDluKiY+MwfZX5wbm4YOhQDB6reAVtcT4Q6OjoLFgS8eNFD3m/KlStj5EhMnozvSu2cRcgtNxdeXli50mb58ikSlgPTLTsb16+jdWtJ13CnUkZMbUZHnZwQH69UwVr9+li9Gv36ISODurA4b8sWWceEqa2U+alVKyjzKA/gr7+gpYW//6YoIKZwPRECqFoVZcogSf4xttWroaeH4rZQJWg0Zw7q1WP/n/30abRoIWnhBNg+hvBPapMIjYxQr56yy2b69UO7dqVosvD6daSnw81NpouprZT5qVUrpZ4IAfD52LkTO3bIcSYXF6hAIgTQooUiux7weNi1C5GRuHuXhpiI4ly6hD17sJkDe/tIHRf98AHp6VCoRI4uapMIQcVHKoCVK5GcjLVrUSD7mb8qS1wmI8sSvmfPoKGB2rWpj6FZMzx/DiUPX69UCRER8PVFUpLK/NTUOREC8PCArS2pmmFIRgYGD8amTZxYV3vypEwThGxt/FYsdUqETk7KDrIB0NFBaOi3iRPbVKnS8tgxude3qZCsLOzfL+uirwsXpIz5K0xTE7a2FGwT07IlOnQ41qhRSyurtt9VYXZKzRMhgCNHkJKiDrvhcd+4cejcGZ06sR3H/xdONGgg6RqujYsCMDfHu3eK7E/GQeI9oZQvgtXQeKWnp//p0+idO89TEBZX7dkDV1fIOKVO0wShGCWP8gC0tc/n549//FgzJeUNBc3RTDUSYbNmePQI2dmK3FupEsaMwdSppGqGXocO4fJlLF7MdhwAgOhoKfn43bt3kVPIx4oAACAASURBVJEhZmbvmIpIJhoasLDA/ftsx0GFatWgr48nT5Rtp0GDBtOnt3dyij93buy1a1RExknicVEZ0TRBKKZk4ehPf/01vkePS/Xrd504sW5ODgUN0ko1EqGODpo0wY0bCt6+ciXKlMGAAZTGRBSSlobRoxEWJqU4hTFSJwjd3Prdv58ZGNifqYhkpU6jo1Q9W8yePf7SpfXbt1ft1k2OLfhVyN27ePsW7dvLdPGLF8jLg0Kr3mXi6IgbNyg4BqR69eqRkRtv3x5TsSI8PMDgYQ2KUI1ECKBlS8VHRwGEhSEqCgkJ1AVE/J9IBD8/DB8OBwe2QwEAZGXhxg0pX5m1tMrweMkmJrScBaoMdUqElEwT/tSpE/bsQc+eiImhrE2OCA7G0KHQ0JDpYjpWEBZmYIAGDRR/6ihCQwPbtqFePXh4cHrMX2USoYODUomwQwe0aIG+fakLiPi/DRvw7h1mKbijFvViYqQvnJg48WDr1j0vXIhkKihZqVMipOqJ8KfWrbFvH7y81OEQ459ycn7UWMqI1nFRMWp/cBoa2LIFzZrBzQ2fPlHWLLVUJhEq+UQI4NAhvH6NZcsoCogAACQl4e+/ERYGLS22Q/k/qeOiAG7c0OnSxVlHR4eRiORgbY379yEUsh0HFZo0wadPFG+45eyM6GgMHw4GDwin1969aNECNWvKej3dT4Sg+lEeAI+HVavg6gp3d3z4QGXLVFGZRFirFng8pTagqlQJU6ZgzhwFi26IPxUUYOBAzJ8vpT6TSffv3z94cH/btlKmOK5eRYsWzEQkH0NDVKigyPYRHMTjoWVL6mf1bG0RHY1Ro7BrF8Uts0KuMpm3b5GZKWkfeUq4uODSJYq/jfF4WL4c/frB1RVvuFdGqjKJEEqPjgJYtAjGxujPuQoJVbVwIYyMMGIE23H834cPH1xcfD5+vBwaukTCZTk5ePgQNjaMxSUfNRsdpfbZQqxpU5w9i5kzsWMH9Y0z6dEjJCdLH8D46dw5uLjQvvi1ShWYmODRI+pbnj4dPj5wc8Pr19Q3rgxVSoTKrCb8ae9eREdDjeuwGXPjBoKDsWMHh9aka2lp5ecLtLU/GBhIOvzpxg00aQI9Pcbiko+aJUJqpwl/atQIMTGYN0+1lwgHB8PfX45pBQYmCMXo+8FNn46xY+HszK1hj1KXCF1d4ehIqmaU9e0bBg7EqlUwNWU7lEKMjY3NzU/899/oefMkHY/O2XFRsejoIcuW2a5bt4PtQChgZ4eHD+kqnW/QABcvYtUqrFpFS/t0y8tDWJgcZTJgZIJQjPJpwsJGj8bMmXB1RWIiXV3IS5USoZ0d7t+nYF384cN4+xZLJA2eEVJMmcLFKtz375GSUiUgoCVP4lPqtWtcWenxJ4FA8OzZzdzcjbt2qUM1iI4ObGxoHICpXRtnz2L9evzzD11d0OfAAdjYoF49Wa9PS0NaGpo0oTOm/2vVCrGxNLY/fDiWLEH79rh3j8ZeZKdKiVBPD/XrU3BEmYkJpk/H7Nnh/foFZHF8nSf3fP/+/cCBL0ePYuVKtkP5w6FD8PCQPtDE5SdCDQ2NBQsmaGgs/fffuWzHQg26P1Jr1kRsLPbvx4wZSE1NpbEnqslVJgPgwgU4O8u0K7fyGjZEVha9M3kDB2L5crRrh7i4gnfvWN7jSZUSIahYRCHWsGGYQBC6d28ZHx/OVHqogvv37xsbWwwe7Nm16zYTE7aj+UNUFLp1k3LNmzfIy0OdOowEpJAxY4bY2+/V0pLvJGrOonWQTaxyZcTEICioo4VFf2PjRvR2RpGkJDx4gK5d5biF1i1Gi+Dx4OhI+w+uXz8sXfrRyamuuXmXUaOm0NuZRCqWCJUvHBWrXbs2n/8IuPrgAVP/Z6mFq1evCgRNgK7JycfZjqWorCxcvowOHaRcduUKWrZkJCAlWFmpz9lhTk64dg10H6NUqRKApyLRssxMpTcHY0RwMHx9IddCVsYqZcQY+AYDoF69RzxeJZFo3IkTD2nvrGQqlgipeiJ0cnJ6/Dhm27YVycnDySFNsjM19dPXL2dtfWnbttVsx1LUiRNwcoKhoZTLuDxB+JOlpfokQhMT1Kql7CG9sjhwYE316ssMDILT0mjvS0n5+di1C/7+ctzy+TNevWJ0zQ99haO/99JqyBCnmjVPZWZuo2pfNwWoWCKsVw85OaBkIsDMzGzIkBZnzuDIEbKyUCbZ2Rg7ln/gQEhcXJQpp6pFAcg2LgpuTxD+ZG2tPokQTH2kdurUKTHxwPDh7SZLKhnmhKgoNGiA+vXluOXiRTg6yrofKSVsbfH0KTIyaO9oy5agxMRdmzdX7tqVtaNXVCwR8nhwcKCyCM3ZGTExOHBAviLm0mnGDLi5oV07pY+Yo4FAgBMn0LmzlMvy83H7Npo3ZyQmJVhaqs9Ga6BtWX2x/vkHcXGc3ox0/fodw4f3cXWVb2iLyQlCMS0tNG9OzQicLHr0wJo1aN+enTUVKpYIQdFqwsJcXXH6NMLCMGQIlc2qmWvXEBnJ3Z1aL15E3bqoXl3KZXfuoF496cOnrDMyQoUKSE5mOw6K0F04Wpi+PjZvxtix+PqVoR7lIhKJ5sxZkZ4+JypKvgUfjK0gLMzZmblvMAB69cKqVXB3p2VTG8lIIgSA1q1x8iR27sTQoRS3rB7y8uDvj9WrwcFKUTF1GhcVU6d6mVq1oK2NZ88Y6q5NG7i5YfZshrqTC4/HMzVtoa/v7+vbS/a7MjLw7BmaNaMvruJRdUiv7Pr2xb//okMHpr8FajLaGxXs7XHrFvLyoK1NZbNubjh0CN27Q1sb69dT2bIaWLgQdeqglxy/uUw7fFim0bBr1+DmRn80VBAnQrWp5BJPE5qZMdTdypVo0gR9+6JVK4Z6lJ1AEHz+POzs5LglNhYtWlD8iScLR0ckJFD/YSuZry8EAri54fx51K7NUKeq90RYtizq1KHly3KXLjh4EMHBmDiR+sZV1/372LABGzeyHUfJ7twBn4/GjaVfqUJPhJaWXNl0gxIMP1sYGWHlSgwdSsFGVNS6cgV8vnxZEIwvnPjJ0BBmZhTsYSIvf39Mngx3d2rqImWheokQ1C2i+FPXrjhwAGvWgPuFZ8wQChEQgEWLUK0a26GUTMZx0U+f8OkTh06MkszKSn223gZThaOF9e4NCwssXsxop1KFhsLHR+67WJkgFGNyfrewsWMxejTatKH4PMuSqGQipGpZfbG6dcO+fVi1CtOm0dWFClm9Glpa8i14Yp6MiTAuDvb2DO1QpTxzc7x7h8xMtuOgiKUl3r8Hwyv81q3Dhg2sVeT/6ft37N8PLy/57srKQmIibFnaaIiZZfXFmjABAQFo1er5tGnLHtFcP6Minwq/o++JUKxHD4SGYvlyzJhBYy/cl5KChQuxcSOHDlr6U2oqXr6Eo6P0K69dU5lxUQAaGmjUCA8esB0HRfh8tGzJ9EeqqSkWLIC/PwQCRvstyeHDaNYMNWrId9fly7C1ha4uPTFJ4+yMS5cgYmnN1OTJ+Pp1yLJlNdq186a1I5VMhA0b4vNnvH9PYxcDB2L7dixbhgULaOyF4wICMG0a7cdhK+nQIXTuDE0Zqr6uXlWBPWUKU7Nl9aw8WwwbhrJluXJmYWgovOX/PGdrglCsWjUYGODxY9YCsLOrWbbskbQ003//RX4+Xb2oZCLk8WBvj/h4envx9sbWrZg/H4GB9HbETaGhePdOBeqGZBwXFQqRkAB7e/oDoo6a1cswP00IgMfDhg1YsID9RZlpabh8Gd27y30jixOEYqz84H46ejQkIWFecvLBGzfQvDldH/sqmQhBz2rCP/n6YvNmzJ37wNjY3c9vPO39ccbHj5gxA1u3ynF2NisyMhAfj/btpV+ZmIgqVVC+PP0xUUedlhICsLfH/fv49o3pfs3NMXUqhg1jbXxPLCwMPXrAwEC+u3JycOcOyyMZTG4M9Ccej9egQYPq1fmRkZg7F56eCAig/qhnkgil8PODvv7wjIwp27dfWLcujd3fJWbk5uaOGQMfHxXYiuz4cbi4oEwZ6Veq0MKJn8RDo2rzv5yeHqytaTykV4LJk5GRgZCQfCF729YpNi4aFwcbG+jr0xCQzNh9IiysTx88fAgAVlY4dYrKllU4ESYk0H62i1jnzk14vLF8vu7EiRV1ddG6NaKjmeiXFRERB6tUcT540HXaNI6twCqOjOOiUJFDJ4ooVw4GBnj5ku04qMPWR6qmJiZPvuHv71izpsObN2+YD+D+faSnK7JZKPNbjP7JwgLp6QwtY5DKxASbNmH9egwfjr598fEjNc2qaiI0Nka1agzV1O3Zs+nTp2sCwdXv33lbt+L7d3h6QkcH7dpx5YsShU6ejP3yZbSursbXr1w/zCY/H6dOSd9oW0wVnwihdqOjLNbif/qUAHT5+NEykY1NnXfsgLe3Ikt32K2UEWPmkF65dOyIxETUrYvGjROrV2/n7j7gu3JbJ6hqIgT9iygKMzExAcDnY9AgXL2K3Fxs2YJPn+DqijJl4OGB69cZioRugwdPLlPm5sKFA2vWrMl2LFKcP4+GDVGlivQrv35FSgosLemPiWrqt6z+6lV2FjP4+g7q0+e7jk59V9c2DHddUIDwcAwcKPeNubm4cYMT50gzv+moVPr6WLwYHTpEpqYOO3tW39Pz7sqVuHJFwb2EVDgRUnsek1w0NDB4MG7dwtevmDcPL17AwQHGxhg0CI6OnpUr22zfvp2dyJR28mS1sWODxo7l9hJ6AHKOizZrJtMSC66xslKrwtFy5VCtGjup3cDAYM+eRbVqzbh6lcEz/QAAp06hTh1FtjSKj4eFBcqWpSEmOXFnmrCIOXP6WlltbdOmwNfX+sULTJuGChXQuDECAhAaigcPMGrUlLNnpR/4q8KJsEULxMWxHIO+PqZPx8OHePcOQ4bg1KnXcXF3c3KWDh2qkssPRSJERMi98wUrRCIcPixrIlTRcVGo3dAo2P5I9fJCWBjTnSpWJgMOLJz4ydYWjx9z8WSr+vXr37lzKiYmxMtLJygIly7hwwds2oT69XHkCNq0ebRp071Ll6Tv9a7CibBJE7x9i8+f2Y4DAFCxIlauxOvXlYBcYKlI5PvpE9sxyS82Fvr6qjGEePMmdHVl/ZatipUyYg0b4uVLFpYc0IfFaUIAXl7Yvx95ecz1mJGBU6fQt68i93KhUkZMRwdNmzI3FaUMPT20aoXJk7FvHxISyvJ4T3Jynkq9S4UToYYGjesrFaOtrZ2bmxIePsXcfJ6ZGb1739AhPFzBr67Mi4pCjx4yXSkSqXAi1NRE/frsnNlNk1atcPEia73XrAkLC5w8yVyPERFo106RgzwLCpCQACcnGmJSCLurCRWzb1+1Vq0SZ8yoL/VKFU6EYHA1oey0tbW7dOn44AFMTWFujlev2A5IZvn5OHAA/fqxHYdsZJ8gTEqCvj6qVqU5INqo2eho3brQ1GRzn5eBAxkdHVXsuAkACQmoVw/GxlQHpCgO1stI9vIlFi/Gpk16RkaGUi8miZAWmpq4e/fH18+UFLajkc3x47CwYO4kTGWkpODdO1kf8lR3glDM0lKtEiEAR0c2P1L79sWJE8jIYKKv58/x7JlMOx8VkZSU5Ofnp6/PjT1SAQBOToiPp3G3T8qNHYuJE2WdPVHtRNiyJa5dA3ubRUiiqYl799CwIRo3Zn+fQ1mEh6tGmQyAQ4fg6QkN2ar/VHdcVEzNngjB9iCbiQlcXBAVxURfO3bAy0uRfQrnzl356FHXu3cPpDJ2NK00xsaoUwe3b7Mdh2z278fjx5gyRdbrVTsRVqiA8uXZ3BldMh4P8fFo3BiNG//YGYizvn7FiRPo3ZvtOGQj+7goVP+JUM2WEoLtwlEAAwciPJz2XkQi7Nql4KR7p07tebxAMzO9SpUqUR2X4lj/wckoMxMTJ2LDBujoyHqLaidCcGMRhQQ8Hq5eha0tmjbl0AGhfzp4EK6uqrEn9ZcvuHkTbm4yXZyTg4cPYWNDc0x0qlIFWlrgzIMBBays8OYNPnxgLYBu3ZCQgHfv6O1FXIOt2P97urqeHTpcvXEjWotL296zW/Eru9mz4eGBNvJsnKAOiZCtZfUy4vEQG4tWrWBry92v9mFhKjMuevQo3Nxk3Yb4xg00aQI9PZpjopmajY5qaMDBgc3vr7q66NIFe/fS20tICHx9Fbz38GF068b0wn+pXFzYPKRXRgkJ2LdP7rPz1CERcvmJ8KeYGLRrBzs7Llb3pKUhIQFdu7Idh2xK1biomJolQnCgBJHulfU5OYiKUvDLpUCAEyfQpQvVMSmtenXo6OCp9FV5rCkoQEAAli9HhQry3ajyidDaGs+fIzOT7ThkcPQounSBiwuuXGE7lN/t3o2uXVk+6kVGubk4c0bWjbah+pUyYmp2Qi84MNvUrh1evsSTJ3S1f/AgHBxgaqrIvbGxqF0b1atTHRMVWP/BSbZqFYyMFPn+ofKJUEsLTZsiIYHtOGQTGQlPT7i64tw5tkMpRLEdgVlx9iwsLeX4ukeeCLmpRQvcvcvmjjkaGujbF3v20NW+wtuqATh8GJ6elEZDHS5PE4oXDm7cCB5P7ntVPhGCw6sJi7V/P/r3h7s7o9tbSJCUhBcvZK09YZ1c46Jv3iAvD3Xq0BkQIyws8OwZcnPZjoM6enpo0oTl769eXti5k5aW37xBQoLiyezoUTn+J2cYl58Ix4yRY+FgESQRsmDnTvj6onNnTJ8eeY7tZ8Ndu9C/v2qczCAS4dgxOeYyr1zhxBE2ytPRQd26XF+BIy/WP1IdHMDn03KA2q5d6N1bwRKt+/eRn8/d/X6bNMHHj7QX3Cpg3z48eSLHwsEi1CERtmyJuDiu1zIVsWULLCyWLV26o23bcUeOHGExkt27VWZcND4eRkYwN5f1evWYIBRTv9FRLgyy9e9Py4JChZcPAoiKQvfulEZDKT4fLVpwrsohMxOTJmHjRjkWDhahDonQ1BRlyiApie045NSiRQqgKRIZpadnsRXD9esoKICdHVv9y0eucVGoywShmJodTAjA3v77hQu7b95kc6uSwYOxZw/FBwXfuIGsLDg6Kng7lycIxSwskoODQ798+cJ2ID8IhcJZs9CpE1q3VrwRdUiEUJ1FFIVt2BA0b551mTL/7N8/gK0YwsIwaJAic8uskCsR5ufj9m00b05nQAxSvyfCZcsCc3JutWs39AN7S+vNzFCtGsWVa6Gh8PVV8HfqzRskJcHZmcp4KLdjR89Tp9J69QpgOxAAiI29XKVKi02bXKdPV2r3WPVJhBxfVv8nDQ2Nv//+68QJt2PHcOoUCwEIhdi3DwNYy8JyEAqFGzcefv/+nOwPr3fuoF49GErfd141qF8iNDEx1NZ+lZ+fr62tzWIY1C4ofP/+0549XwYNUvD2I0fg4cH1CXtDQx0gWU+PE79aJ09e+PjRT0fH9NMnpZbCqE8iVLknQrFWrdClC/r1Y2Hr8DNnYGqqYJEVwyIi9k2YsDcra8nVq7L+mNVpXBRA9erIz0daGttxUGfevMnLl09s0CDGyMiIxTD690dUFHJyKGgqPj6+QQOP9PS2374puJtiVBTXx0UB3L4d4+bm1a3bRrYDAQAvr6FaWndGjLBqrtzgj5okwmbN8PgxsrPZjkMh+/ejoABDhjDdrwotHzQ2NiooeK+rm25gYCDjLepUKSOmZucx8Xi84cPtnz+vyO42qqamaN4cR49S0NTz5y+ysqw1NBq9fPlSgduzsnD5Mjp0oCASWpUtW9bfv9WRI5zYAS4mplL//hv++28Wn69ULlOTRKijA0tL3LjBdhwK0dJCeDh27WJ0WdX37zh8GH36MNejMsqV61i1atCtWxGWMteVq9kTIdRxdFRLCx4eOHaM5TC8vKipHbW07KWv77B8ubuHh4cCt588CUdH1RjM79QJFy9y4sFj925qJnfUJBFCBVcTFta1K1xcGN1d8MgRNGuGatWY61EZmzdjzBiLOnVqy3j9p0/49Ek1Rn1lp2ZPhGLdujF0NKAEvXrh3Dl8+qRsO6tXa0yZMnTUKB+eQqUy8hZFs8jICPb2OH2a5TCSk5GcjHbtKGiKJEKuOHYMmZmYMIGh7sLCVGZcNCsLkZEYPFiOW+LiYG8P5QZLOMfaWg0ToYcHLl3C169sxmBoCHd3REYq1cjnz9i3DwGKllJydqPtknDhG8zu3ejTh5raIvX5qODgMk+56OsjOBhr19K4EfBP6ek4fx49etDeESXCw+HmJt/+xdeuqdu4KIAmTfD4MfLz2Y6DUgYGcHJif7vBgQOVrR3duBE9eqByZQVvj41FrVoc3Wi7WN274+hRFBSwGcOePZQVvatPIqxVCxoaSElhOw4lDB4MGxsoNL8gn3370L49WC3Wk8PmzRg2TL5brl5Vt0oZAHp6qFGDie9JDOPCs0WnTrh/X/FPj/x8bNyIsWMVD+DwYZUZFxWrVg21arH57HHnDjIzFd+4oAj1SYQAy6d9UuL0abx6hfnzaexi0qR/xo2z19VdR2Mf1Ll7F+/fyzcNIBQiIQH29rTFxB71q5cB0L07oqNZftLV1BTq6Q22sWkRHa3Ikt69e9GgAaytFQ/g6FEVWDhRBLvfYMR7Q1K1GYhaJUJVXFZfhIkJlizBwoWgr6Z8587I3NxjsbGhdHVAqU2bMGwYNOQp1U5MRJUqKF+etpjYo34HEwKoXBnm5oiNZTOGDx8+ZGWlfvmydPPm/Qrcvnq1UrP79+8jLw9WVoq3wIpu3ZSdWFWYSISICCo3A1GrRPjhw65Nm1yDg3exHYhSJk6EmRmNA6T16k2pUmXgsmVz6OqAOjk5iIiAj498d6nfwomf1PKJEGw/WwCoXLmyv38bA4OlDg6j5b03Nhbp6Ur9wh4+zOmNtktiZQU+H/cV3DxAKZcvQ0+PyjM61CoRhoevzMk5+s8/K9kORFkxMUhMRFAQ9S0nJeF/7Z17VFNX1sB3SEJ4JAgooFngKsUHAjJ8CggWklBXQbGigxSxi1GL6FhRXFpG64NaHWuRZUdZfNblIPqJWpGB8YEoSkcTdVqePrAw1EhAQNASK0ggIc/vj9SUQV5J7oOE8/srN9yz905uuPuec/bjyZP4//znRkzMiLsZkUdeHgQHw+TJ+o0yv1R6Hb6+8PAh2UbgwOLFcPEiyQ1kDh5M/cc/ruTk/EHfGk+HD8PmzUaFKJtQ4kQ/oqLIeYI5dw4MrmM3IGblCDds+MTScl5UFOE1WrDG1RV27oStW+HVK4wl//WvkJwM9vYYi8WJ48chMVHvUWY8I3znHZBIMMh4G214eYGVFfk+fv58cHSEggI9hjQ2wu3bhjddAoDWVhAKISTEcAkkQspUXqmEggKMiySblSPcvn1DQkK5p+cGsg3BgD17wMUF47wikQiuXDEqto1Ifv4Z6ushMlK/USkpaXV18TY2T/AximQoFPD2NsNtQgCIioKLF8k2AmDHDtizR4/CvxkZkJgItraGaywshMhIoNMNl0AioaHQ2AgtLYQqLSkBDw9wd8dSplk5QgCYMwfKy8k2AiOKi6G0FMvS+F99BRs2gIMDZgJx5e9/h08+0e8G0djYePw4X6mM2b//f3Gzi2TMMq0eRsE2oZbISLCxGalL7uqC06dh/XqjNJpc4kRfqFRYsAAI7iyOVVm1vpibIwwMNPnAUR1eXpCYCGvWQE8PBtKamuDyZUhOxkAUAcjlcOaM3oXI2Ww2k0lnMr9atkzPiaTpYJaBowAwdy48fw4NDWTbAZCaCnv3jmjD8vhxCA8HNzfDdZlKoe0hIPgJRiaDoiLsiySbmyOcMQNevgSxmGw7MOLYMWCxIDoaA1H79sG6deDoiIEoArhwAXx9YcoU/UZZWlpOmVJYUFARGRmOj13kY66BoxYWEBmJTRcII1m0COh0uHx5mNNUKjhyBDZtMkrX9esQFGQahbYHIyICfvwRCOtXf+kSBAQYXsFnMMzNEVIoMGsWoW0c8KawEEpK4No1o4Q0N0NBgclMB8GgajIA0NsLVVUQHIyDQaOGmTOhpgZUKrLtwIFRsjoKADt3wp49w0wKL10CFxdj45NNel1Ui60tcDjEFcnDY10UzM8RgnltEwJAYCAsXQpxcUaV9du/H/78Z3Byws4sPGlogOpqQ24Q5eXg5QUsFg42jRpYLHBxgfp6su3AgQ8+gMpK7COlDWDxYlCrh+kPpc2aMAaVCq5dg0UmkMc0DIQ9wXR0AJ+PS86lGTpCc9om1HLuHPT2/p+V1azExM8MGN7SAnl5xPW1MJ6sLFixAhgMvQcKBMDl4mDQKMNcswmtrSEsDK5eJdsOAAoFUlNh9+5BJ4X37kFTk7F35Lt3TazQ9mBERUFxMcjluCvKz4cPPsClSLIZOsLgYCgvJzk5F1uoVFAo9qhU3584YchzV1oaJCaCszPmduGCUgmnTkFCgiFj79wBDgdrg0Yfvr7mGS8Do2l1NDoalEooLh74r3/7G2zcaGwDoMuXTa++6IA4O4OnJwgEuCvCaV0UzNIROjuDra25rR1xuYEWFgEUysfOzlBZqcfAtjbIzYXPDJlJksOVK+DhAV5eeg9UKqGszFQTk/XCXONlAGDRIigpAZmMbDsAKJTfdgrfprUVrl2D1auNVVFYaPIbhDoIeIJpa4P79/GqPYmlI1Sr1Wp9yxPhw5w5pK2OqtXq7u5uzMXevHlepaqXSPb6+8OcOcN0qVUoFFKpVPs6LQ1WrcJyOthXOB4YFiYDAFVV4O5OWtEcuVwuI+r+TZYjlMlkcpzXv8aPB19fuHXLwOESiUSD3VpQTAx0d8ONN+0odMKPHIH4eGN/aTU1xBXa7sK/8XF0tOFF8pRKZc8IUsRyc+GPfwRra1yEj9QRfvfddxQK5eDBg0OcU19f/+TJqKjoQWK8THV1dXg4XrH71tZw9Srk5cE//wnOzlBV2ng+gwAAC0xJREFUNfBp586dS05OBoDnz+H0adiyBUsbTp8+/RluE8yWFigrg6VLDRl7+zaZG4TffvvtF198QYwuDw8Qi+H1a2K0/U5qaurRo0fx1mLM3GLWrFltbW1YWWJhATt2/N4TbebMmWKxWCqF48cxCMC+dIm4QttTp06VSCQ4qwAmE+7dM2Ts9evX/zT0oz0AGLouWlRUNBLhI3KEYrE4LS0tyHQKOJpfvExfli4FsRj8/SEwcJipYXo6rFwJbDZRlhlNdjYsXw42NoaMvX17TGwQAoCFBXh5mfM24cWLehQ5w5XYWHj1Cv71r9/fycmB4GDw8DBWstlsEOpYsgTH1dH6enj6FHg8vOSPyBEmJSXt2rXL0VSSsQFmz4affoLeXrLtwI2RTA3FYsjJgZQUwo0zFLUaTp40cOtFrYZ//3tMbBBqMeNtQg8PcHIaLanAVCrs3Am6qb5GAxkZxmZNAMCLFyAUQmiosXJGFbhuE545A8uXGxudNATDCy4sLJRIJLGxsadOnRr25B/f6hDv7OxsT8a+javrtPz8Fh8fLKqT6cPjx4+lUulDQsLbp0yB77+3SEp6NzDQdvny9r/8pVX7fnNzc0dHR0pKe3g4RSx+hm2dnZaWll9//RWPD3j3LovFmkihCA2QXVdn7eAwua3tZ+wWxvSjra1NLBYTc90BwNFxwq1bVnPnElrtWCwWW1paEvAZg4ImZmVRrKz0vpYKhaK2tra9vR1DY7y9Ka2t07OzW5RK5cmTbUoly97+ZyO/g4KC8UFBtrW1TRjZOAxqtfrRo0c2hq20jBgrK2ht9SoqeuLqqt9GckNDQ1dX19C/q5yc6Xv2ND98qPf9vKGhYSSRKxTt9m9ubu4PP/zQ72+urq5r164NDg6+ceOGm5vbwoULw8LCUgafYixatKikpIRCofR908nJiRRH2NPja2nZQqP9SrBelUollUqZTCaRSl+/nqdQTBg//rz2UKFQKBQKgCBLy2YaDeP8ZLlcrlQq8finksvZarWtlZXQgLFKpYNc7mZjQ9osSS6Xq1Qqa3238g1FoXBSKl2srQltiiqVSqlUqqWlJd6KFAonpdLJ2rpW34GvX79msVj9bkHGI5NNo1IlUmkdgzFdo7G1sjI2EqK39x0ADYPxFAvrhqezs9POzg7zr+Vtenr+h8EQUamdeo1SKpW9vb22Q7XwsOjqeo/FugugdzSOQqHYvXt3XFzc0Kf95ggFAkFtbf+fnZOTU0VFxYMHD5YuXQoAmZmZnp6ea9aswS8YBIFAIBAIgvltaZTL5XIHCrmTSqUdHR1VVVUA0NnZ2dLSUm9mCXoIBAKBGNtQRp52M+zSKAKBQCAQJoceCfUhISHTp0/HzxQEAoFAIIhHjxkhAoFAIBDmBzYl1ur+m85O/aKGzIPGxsa+xdVevXrV2tpKjOru7u7m5mbd4cuXL589e6Y7bG9vFxudQtFPhUQiEYlECoXCSLFmQFtbW98ff99vyWzQJuSQbcWgPH/+vK95Uqn06VMMAjKbm5t1/9EqlUokEumKCyoUCpFINMpnEQqFoq6ujlwj1Wq1SCRSvukhp/0aX79VFenzzz/vW6Hw7a+33yUeod7ePrnk7e3tv/zyy6ADNEajzdKYMmWK9xsuXrxovFiTY8uWLQsXLtS+lslk3t7eJ06cIEZ1dXU1g8Ho6enRHkZFRdnb2yuVSu0hl8vNyMgwUkVBQcH06dO1r58+fTpt2rRt27ap1WojxZoB8fHxEyZM0P34ExISyLYIe3g8XmZmJtlWDMqHH3749ddf6w5LSkpcXFyMFxsREZGenq59ffv2bQA4duyY9rCwsNDNzc14FbgiFAoBQCaTkWiDdlLU1NSk0Wh6e3tjYmK4XG5nZ2e/0xwdHTs6OnSHjx8/BgCpVKp7p98lHgkxMTFJSUna1+3t7S4uLsXFxYOdjFnR7UuXLv30hsVmU1NdH/bv39/Q0KAtO7Br1y42m71q1SpiVPv4+LBYrLKyMgBQq9Xl5eWenp7V1dUAIJPJysrKeNjVJhKJRGFhYR9//HFaWhoBmUkmQUJCgu7Hn52dTbY5CGzgcrmCN72F+Hw+j8fTHQoEgrCwMPJMMz16enqioqKkUum1a9fs7OwI0Hj06NH8/PySkhIAWL9+fXR0dERExGAnm2EbJrJgMBjZ2dlbtmzJy8s7ceJEdnY2YX6CQqFwuVw+nw8ADx8+nDp1akREhPawtLTU1tbWx8cHE0U1NTUcDic5OXn37t2YCEQgRi08Hu/OnTsqlQoABALBjh07tPNCAODz+QPmmyEGpLOzMzw83MHB4cKFC4QVnZgwYUJGRsbatWuzsrIqKioOHDgwxMnIEWJJUFDQihUr4uLivvnmGzc3NyJV655etf+iHA5Hd8jj8SwsMLjQL1684HA4Bw4c2LRpk/HSEIhRTkBAgEajuX//vlwuf/ToEY/Hc3FxEQqFXV1dDx48wHCVxeyJiory9PQ8e/YsnU4nUu+yZcv8/PzWrVuXnZ3NYrGGOBM5QizRaDR1dXV0Ol1D+AY1j8crLS2VyWQCgYDL5QYHB5eVlanVau0hJioYDIaDg0NFRQXxnw6BIB4ajTZ37lw+n19RUeHn50en0zkcDp/Pv3PnzqRJk959912yDTQZ2Gx2bW0t3q2g3kYul9fX14/khowcIZZkZ2c3NjZev349JSWlpYXQgsg+Pj52dnalpaWlpaVBQUHW1tYeHh4VFRVlZWVYbWbY29vz+fyioqLNmzcjX4gYPdjY2PRtvtrd3Y1VLVztviCfz+dwOACgXWgRCATvv/8+JvLHCDk5OS4uLvPnz387XnRAtJevbxC+RCIZshjpwOzdu9fJyenkyZOrV68eujsxcoSY8ezZs+3bt586dYrH461cuXLdunVEaqdQKKGhoZmZme7u7tqS3xwOJz09nclkent7Y6XF1dX11q1bV65c2YJtt18EwgimTZt2r09P2KqqKk9PT0wka7cJb968qV1W0c4I0QahvtDp9Ly8vJH7QjabzWKxdNe0t7e3pqZG32t6//79I0eOZGVlLV++3N/ff/v27UOcjFt/pzGGRqNJTExMTEwMDAwEgH379vn5+Z05cyY+Pp4wG3g8XnJy8tatW7WHHA4nLS0tOjoa25gdrS/UzjIPHTqEoWTEaObq1au6bFQ/P78lhLVXHwFJSUm+vr4bN24MCQl5/PhxZmbm5cuXMZHs7++vVqtLS0sDAgIAwNHR0dHRsaKiIjc3FxP5BLBv3z4qlap9nZCQMHnyZFLMoNPp58+fj42NXbBgQXFx8dA7dhQKJTU1dc2aNdu2bbOzszt79uzUqVPnzZs3cnVyuXzlypVpaWnaFeyjR4/OnDnzo48+GuwJhvrll1/q83EGNppGo3G5XLz7XY1m2traJBLJ1q1baTQaANDp9JCQkKamptmzZxNmA5vNdnZ2jo2NdXJyAoCJEyfa2dnFxcVhErZDpVLZbLb244wbN27JkiWNjY1sNtvBwcF44SYNjUbz8vJyd3cn2xAcodFoLBbL4g2TJk0aVdUWmUzmihUrhEJhZWWltbX1oUOH/P39MZFsYWHh6uo6f/58nUBXV9eAgIAFCxZgIh9XKBQKk8mk0Wi6C+fn5zdu3DiCzWAymaGhoQwGg0qlRkdHS6VSpVLZ7/8lPT1906ZNVlZWunfee++9GTNmVFZWCoVCLpd7+PBhvQJthEKhvb39p59+qp0G2NraBgQEvHjxYrDlMVRiDYFAIBBkMn78eJFIRLyT1oH2CBEIBAIxpkGOEIFAIBBkkp+fb0BQKIagpVEEAoFAjGnQjBCBQCAQYxrkCBEIBAIxpkGOEIFAIBBjGuQIEQgEAjGm+X8ApFc3Q9rM3AAAAABJRU5ErkJggg==", "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.0004494364987582204, 0.00044944118256774946, 0.00044943725507162996]\n [-0.000449441182568017, -0.0004494364987580068, -0.00044943725507160426]" }, "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.3" }, "kernelspec": { "name": "julia-1.5", "display_name": "Julia 1.5.3", "language": "julia" } }, "nbformat": 4 }