{ "cells": [ { "cell_type": "markdown", "source": [ "# Tutorial" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "DFTK is a Julia package for playing with plane-wave\n", "density-functional theory algorithms. In its basic formulation it\n", "solves periodic Kohn-Sham equations.\n", "\n", "This document provides an overview of the structure of the code\n", "and how to access basic information about calculations.\n", "Basic familiarity with the concepts of plane-wave density functional theory\n", "is assumed throughout. Feel free to take a look at the\n", "[Periodic problems](https://docs.dftk.org/stable/guide/periodic_problems/)\n", "or the\n", "[Introductory resources](https://docs.dftk.org/stable/guide/introductory_resources/)\n", "chapters for some introductory material on the topic.\n", "\n", "!!! note \"Convergence parameters in the documentation\"\n", " We use rough parameters in order to be able\n", " to automatically generate this documentation very quickly.\n", " Therefore results are far from converged.\n", " Tighter thresholds and larger grids should be used for more realistic results.\n", "\n", "For our discussion we will use the classic example of\n", "computing the LDA ground state of the\n", "[silicon crystal](https://www.materialsproject.org/materials/mp-149).\n", "Performing such a calculation roughly proceeds in three steps." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "using DFTK\n", "using Plots\n", "using Unitful\n", "using UnitfulAtomic\n", "\n", "# 1. Define lattice and atomic positions\n", "a = 5.431u\"angstrom\" # Silicon lattice constant\n", "lattice = a / 2 * [[0 1 1.]; # Silicon lattice vectors\n", " [1 0 1.]; # specified column by column\n", " [1 1 0.]];" ], "metadata": {}, "execution_count": 1 }, { "cell_type": "markdown", "source": [ "By default, all numbers passed as arguments are assumed to be in atomic\n", "units. Quantities such as temperature, energy cutoffs, lattice vectors, and\n", "the k-point grid spacing can optionally be annotated with Unitful units,\n", "which are automatically converted to the atomic units used internally. For\n", "more details, see the [Unitful package\n", "documentation](https://painterqubits.github.io/Unitful.jl/stable/) and the\n", "[UnitfulAtomic.jl package](https://github.com/sostock/UnitfulAtomic.jl)." ], "metadata": {} }, { "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "n Energy log10(ΔE) log10(Δρ) Diag Δtime\n", "--- --------------- --------- --------- ---- ------\n", " 1 -7.899836980555 -0.70 4.8 \n", " 2 -7.904927493139 -2.29 -1.52 1.0 1.92s\n", " 3 -7.905182173903 -3.59 -2.55 1.4 790ms\n", " 4 -7.905210975759 -4.54 -2.92 2.5 232ms\n", " 5 -7.905211374737 -6.40 -3.22 1.1 153ms\n", " 6 -7.905211526104 -6.82 -4.57 1.0 85.2ms\n", " 7 -7.905211531382 -8.28 -5.27 2.4 122ms\n" ] } ], "cell_type": "code", "source": [ "# Load HGH pseudopotential for Silicon\n", "Si = ElementPsp(:Si, psp=load_psp(\"hgh/lda/Si-q4\"))\n", "\n", "# Specify type and positions of atoms\n", "atoms = [Si, Si]\n", "positions = [ones(3)/8, -ones(3)/8]\n", "\n", "# 2. Select model and basis\n", "model = model_LDA(lattice, atoms, positions)\n", "kgrid = [4, 4, 4] # k-point grid (Regular Monkhorst-Pack grid)\n", "Ecut = 7 # kinetic energy cutoff\n", "# Ecut = 190.5u\"eV\" # Could also use eV or other energy-compatible units\n", "basis = PlaneWaveBasis(model; Ecut, kgrid)\n", "# Note the implicit passing of keyword arguments here:\n", "# this is equivalent to PlaneWaveBasis(model; Ecut=Ecut, kgrid=kgrid)\n", "\n", "# 3. Run the SCF procedure to obtain the ground state\n", "scfres = self_consistent_field(basis, tol=1e-5);" ], "metadata": {}, "execution_count": 2 }, { "cell_type": "markdown", "source": [ "That's it! Now you can get various quantities from the result of the SCF.\n", "For instance, the different components of the energy:" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "Energy breakdown (in Ha):\n Kinetic 3.1020977 \n AtomicLocal -2.1987839\n AtomicNonlocal 1.7296072 \n Ewald -8.3979253\n PspCorrection -0.2946254\n Hartree 0.5530398 \n Xc -2.3986215\n\n total -7.905211531382" }, "metadata": {}, "execution_count": 3 } ], "cell_type": "code", "source": [ "scfres.energies" ], "metadata": {}, "execution_count": 3 }, { "cell_type": "markdown", "source": [ "Eigenvalues:" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "7×8 Matrix{Float64}:\n -0.176942 -0.14744 -0.0911692 … -0.101219 -0.023977 -0.018408\n 0.261073 0.116915 0.00482507 0.0611643 -0.023977 -0.018408\n 0.261073 0.23299 0.216733 0.121636 0.155532 0.117747\n 0.261073 0.23299 0.216733 0.212134 0.155532 0.117747\n 0.354532 0.335109 0.317102 0.350436 0.285692 0.417258\n 0.354532 0.389829 0.384601 … 0.436925 0.285692 0.417516\n 0.354533 0.389829 0.384601 0.44929 0.627662 0.443806" }, "metadata": {}, "execution_count": 4 } ], "cell_type": "code", "source": [ "hcat(scfres.eigenvalues...)" ], "metadata": {}, "execution_count": 4 }, { "cell_type": "markdown", "source": [ "`eigenvalues` is an array (indexed by k-points) of arrays (indexed by\n", "eigenvalue number). The \"splatting\" operation `...` calls `hcat`\n", "with all the inner arrays as arguments, which collects them into a\n", "matrix.\n", "\n", "The resulting matrix is 7 (number of computed eigenvalues) by 8\n", "(number of irreducible k-points). There are 7 eigenvalues per\n", "k-point because there are 4 occupied states in the system (4 valence\n", "electrons per silicon atom, two atoms per unit cell, and paired\n", "spins), and the eigensolver gives itself some breathing room by\n", "computing some extra states (see the `bands` argument to\n", "`self_consistent_field` as well as the `AdaptiveBands` documentation).\n", "There are only 8 k-points (instead of 4x4x4) because symmetry has been used to reduce the\n", "amount of computations to just the irreducible k-points (see\n", "[Crystal symmetries](https://docs.dftk.org/stable/developer/symmetries/)\n", "for details).\n", "\n", "We can check the occupations ..." ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "7×8 Matrix{Float64}:\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 2.0 2.0 2.0 2.0 2.0 2.0 2.0 2.0\n 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0\n 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0\n 0.0 0.0 0.0 0.0 0.0 0.0 0.0 0.0" }, "metadata": {}, "execution_count": 5 } ], "cell_type": "code", "source": [ "hcat(scfres.occupation...)" ], "metadata": {}, "execution_count": 5 }, { "cell_type": "markdown", "source": [ "... and density, where we use that the density objects in DFTK are\n", "indexed as ρ[iσ, ix, iy, iz], i.e. first in the spin component and then\n", "in the 3-dimensional real-space grid." ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "Plot{Plots.GRBackend() n=1}", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzdd2BTVdsA8OfcJN17pS3QXXbpoAVaZgFlKaC4EHDxvjgQBEVBZYuKLx/T160o+KoMEcSFAi2zpUBLoXSwpdCRdDfdSe75/ojWUlpa2tzc9fz+asLpzekh5zz3zEsopYAQQgjJFcN3BhBCCCE+YSBECCEkaxgIEUIIyRoGQoQQQrKGgRAhhJCsYSBECCEkaxgIEUIIyRoGQoQQQrKGgRAhhJCsYSBECCEka1ILhN9//31SUlI7ExuNRk4zI0lYaB2AhdYBWGgdgIXWMVILhEeOHElNTW1n4pqaGk4zI0lYaB2AhdYBWGgdgIXWMVILhAghhNBdwUCIEEJI1jAQIoQQkjUMhAghhGQNAyFCCCFZw0CIEEJI1uQbCP/33fZH//Xifza8bzAY+M4LQjzIzMx87F+zZ817LT8/n++8IMQnmQbCH/b8+OL7O3/rOXdFws03Vr7Dd3YQsrSKiop7Hpqx3W3K58bY+Psf5js7CPFJpoFw7/7DFUNnQ7d+Nfe8tuv3Q3xnByFLO30uq9I/DroPpeH3aVj7wuJSvnOEEG9kGgjjY6MdT30N5QXKY5+V+MRE7zF8foGtxiFSJAPnSukLx40Pnw80Xj4B+dlwJbleVxqxz+H1U8ZrOsp37hDigUwD4ZPTH186uX//X5+fE1xbsG3F6hjFHzep33f6Z48Zz5ViW4AkqN4IO6+x9/xmGLfP6GQFZ59QH/z6/eGnV93355fnfvvu2EQVBRi013DPb4ad11gDy3d2EbIgQqmk2v25c+eGhobOmTOnPYl1Op2jo2Pjy4Ia2HqJ/TCb9bSBWT2ZGSGMrZKzjIpWs0JD7cFvoV2ppJ9dYL+6yIa5kVk9mQf8GWUrN8D1Rtiby36aw2aVwYxQMrs3082eWDaz/8BvWgdgoXWMTHuELfKxg4XhzLVHlatjFD/lskHb9YtwsAiJFkvhQB595KBx4I+GWgMkTVTuH6d8OLDVKAgA1gp4OJDZP055YLyi1gBhuwz3/2E4kCetm2WEboNdnuYYAqO7kNFdlJcr6ecX2EF7Df3auo9GSFAKa2HLRfajbNbDBmb1ZLYMV93t2EYvF7IxVrEqWvHdFXZBirHWCM90Z/7Vk3G35ibHCPEKm/ZWhTiR1TGK3MdUs3oyn+awAdsNi04Zb1TjzTESrmOF9JGDxj7f66/q6I/3Kk5PVs7q2fERfkcVzOrJpD+o/DZecVVHQ3foHzloPJCHVQBJDfYI22AaLHo4kMkupx9nsxE/GEb5MrN6MqO6EN4mTxC6VUUDbL/KbspkWQrP9mQ2D1M5qMx5/f4e5JMhivcGKLZcZJ87bnSxglk9mWkhjD22H0gSsEfYXqbBoj8fU43uQhakGHvuNLx3li2p/+tfjUZjaSnuxEKca2hoqKioaHyZWkyfPWYM3K4/kEc3DFJkPaR8qS9j3ijYyMUKXurLXHxYuTpGcSDvr1XWGU1WWZeWluIT0pEFFBUVmXeZJ97R3R3TYNGsnkxqMf00hw3doR/ty8RVJv/n9TlGe/euDszR3/bY2dnxnU0kTd9s/37+kreInUtUqN9jq77cmAU6PfyrB3PxYZWHjYXy8PckuqKgRrH1EnvfH0ZPG3jSr+azFydr9UpGV/TDlo8HDRxoodwgmbl58+bwCVOqrV1ta0sS9u4IDAw0y2Vx+0SnlhoX1cGXF9klU0c1PL0FXHysEjauH+b0wnOzOnNNgcP12R1grkLr0isqf04iWNmR7+YNGTtp6eP38D5Er2fhx+vsG++9f1lHaPzzUHoj/Kdn04/u7/yV8ZvWAZIvtKdemL/VaiTtOxYuHnmoaOfOLz82y2U57BHu37//p59+cnNzmzVrlq+v7+0Jrly58sUXX9TW1j788MNxcXGmN2tqarZs2ZKVleXl5fXEE0/4+/sDwOXLl3ft2tX4iw8//HBQUBB3OW8/Txt4rR/zmY3+so0DABhsnKtqavjOFJIsI6WgtAYAO0fnBT3qRnfhf55axcBDgUxWQMPyK14AADaOdfUNfGcKSVZ1bR11dgYAsHWurqk112W5miPcvn379OnTe/bsqdFoYmNjdTpdswT5+fkDBw5kWdbf33/8+PGHDh0yvT9hwoQ9e/YMGjSovLw8KioqLy8PALKzs99///2yv+n1eo6y3TErXpvn8fEE2x3zbI98/NT0qXxnB0nWw49NYzbe57bjhSBN8r333st3dv7x7yend03a5LxrPrPp/kULXuY7O0iyFs9/XrXtJYddr6i3z1rxartG/tqFciMiImLr1q2mn4cMGfLxxx83S7B06dKHHnrI9POaNWvGjRtHKTWtN7lx40bjRb7++mtK6d69ewcPHtyez50zZ86mTZvamcnKysp2pmxTXl7ewcPHum7RnS5izXVNYTJjocmHuQrtwf2GRb9cOnHihF6vN8sFzaimpubo0aNTd99845TBLBfEb1oHSL7Qvr1sjPqm6MjRo6WlpWa8LCc9Qp1Ol56ePnr0aNPL0aNHHzlypFmao0ePNk1w9OhRAHBycgoICDh58iQAFBQU5OXlhYWFmdJotdpVq1Z9+OGHubm5XOS5k3x9fUcOG/xmjN3rp3DVHOLEqSJ6soguuTdk4MCBSqXglrnZ2toOGTJkzT2+n2SzN3G7LeKAnoWlqeyaEe5DhwxxdXU145U5qU4FBQUA4OHhYXrp5eV1+PDhZmkKCwubJqiqqqqsrHRyctqzZ8/YsWNffvnl4uLitWvXhoeHA4Cjo2NcXBzLsseOHVu0aNHPP/88bNiwFj/6zz//PHXqVHp6uumlQqFYtWqVg4NDi4lra2sVCkWn/9x/TPODDRnKX67WxXtL9tBisxeaHJil0F49oXgzjEJDQ42A5+BcCcwIYlacZjfGdPaOEL9pHSDtQvvkIhPkQAa5NNzVSgwrK6s2bxw5CYRWVlYAYDQaVSoVADQ0NFhbNz+aSaVSNT4avqGhgRBiZWVVW1v72GOPzZ49e/r06ZmZmTNnzoyKioqJiRkxYsSIESNMiZctW7ZkyZLbI6uJi4uLra1tTEyM6aWtra2zs3NrpdBixjppRX+67CwZ46/ifxkDN7goNMnrfKHtuwmFdfTpHkT45/wtjoJeu9iX+ip6uXSqEuA3rQMkXGhVeliTxf50D2NtfXdhi2HarjOcBEIfHx+GYW7evBkSEgIAeXl5t68a7dKli2khjCmBu7u7jY3NgQMHSktLFy9eDAABAQETJkz49ttvG6OaycCBA7/66qvWPtrFxSU0NPS5555rTz4VCoXZ754eCYYNmYYfrpNHggTfYnUIF4UmeZ0sNJbCkjTD6gGMtUoEXyp3W3g5jKxMpztGdep7gt+0DpBwoW08x470Zfp7cvLXcVKvrK2tx44du337dgCoq6v78ccfJ02aBADV1dX79u0zdQQnTpz4/fffm86h2L59uymBp6dnRUWFRqMxXefixYtqtdr0i6Z3KKV79uzp168fF9k2CwKwOkaxOJXVS3ZwFFnad1dYFQOT/EUQBU3m9WVOaOkJLc4UIvMoroP3M40roriqAlxNua9YsWLcuHFnz569dOlSYGDghAkTACA3N3fcuHFlZWUuLi5PPPHE5s2bhw0b5unpefLkSdNimfDw8Mcff3zAgAGjR4/OysqqqKj497//DQAzZ868fv26n59fTk5OVVXVr7/+ylG2zWK4Dwl0gC8usM/1Ek3LhQSrgYVlaexnQxUiGmy3UcCSSGbRKeOhCYJb1IPEaFW68fFgJtiJq0rA4ckyJSUlR48edXNzGzJkiGmUtr6+Pjs7OywszNR51+v1hw8frq2tHT58uJOTU+Mvnj9//tKlS56enoMGDTJN71VXV6elpRUWFvr6+sbExJjmIFtk4ZNlWnO2lI79zXDxEZUjN6c+8kjyR1dwoTOFtimT3Z/H/nSvyCKKkUK/XYb1sYp7O7rrH79pHSDJQvtTR6P3GDIfUqltufoIPGKNqy/N44nGPq7kzQipdQolWdO41uFCq9JD9536X8coI9xF1CH8yw9/sivT2LQHlEyH8o7ftA6QZKHNOGQMcSLLOBsXBXz6BHfejmbWZxi1ZjsDCMnR/2UY7+nCiDEKAsCDAYydEnZcxdly1HEZpXR/HvtyGLehCgMhVwIdyWPBzHvncH896qCiOvggi+X0Rphrqwco3jjNNmAoRB216JRxcaSC6zkmEdcx4VsWpdhykb2qk9TgM7KYlWnGGSFMkKMou4Mmw7xJTxf4LAcjIeqIo4U0qxz+3YPzOIWBkEOeNjC7N7MyDVsBdNf+1NHvrrCLwkW/J+y9GMWqM0adsM7JR+Kw6JTx7WjGmvtKgIGQW6/2U/yRx54txU4hujtvnmbn9VV4cbZMzmLC3MioLsz683g7iO7O7j/ZKj08ZpGTSTAQcstBBQv7KZacxlYA3YVzpTQhn32pr0Sq56r+zKbzRg0uHEPtZqSw+DT7nwGKji05vlsSqWlC9nxvJrOMHirATiFqr9dOGpdGcb5AwGICHMm0EObds7hwDLXXlxdZDxsY09VCE+QYCDlnxcDy/syiU0aMhKg9DhfQixUwk/sFApa0JFLxzWX2SiVWAtS2OiOsTGNXD7DcBLmkKptgTQtmGozw03UcIEVtoACLThnfjWGspFU1PWzgxd6K5bhwDLXDpkx2oBeJ9bLcemlp1TahYgisilYsPMUasB1Ad7TrGltrgIcDJVgxXwljDuazZ0qwU4jupLwB1mYYV/a3aBWQYH0TpvHdiK8dfH0ZIyFqlZHC0lT2/wZaaIGAhTmo4I0IxeLTOFOI7uTddOMD/kwnH2Z5tzAQWs67MYqlqWytge98IKH6/AKrtoXRHT2lWvie7clcqICEfOwUopbl19AvLrCLIy0dmDAQWs4ATxLjST7Mxk4hakGtAd4+w66OEf0O+jtQMbCyP/M6LhxDrViWys7qyXS1t/S9IAZCi3o3hnnvrLG0nu98IOHZkMnGqclACy4Q4MXUYMZAYc+feDuImrtQQX+8zi7ox8O9IAZCi+rhTCb5M/+HJ3GjW5XVw4bzxrcsu0CAFwTgnWjF67hwDN3m9VPsa+EKN2sePlr6FU9oVvRnPs1hb1bj4BD6x9vpxocCmVBniXcHTcZ0JV3t4atLGAnRP04W0VNFdHYvfkISBkJL87UjM3swq85gK4D+kldNt15il0RKeXawmdUDFMvT2BpcOIb+9vop48r+jK2Sn0/HQMiD1yMUu6+z2eXYKUQAAEtS2Wd7Md7iP1+7/aI9SKwX+W8W3g4iAIBfb9D8GpgRwls8wkDIAxcreCVMsTQVWwEEOeX0txvsgjAZdQdN3o5m/u8cLhxDwFJYfNr4Xgyj5C8cYSDkx9w+TIqWntBip1DuFp5iXwtXOFvxnQ+L6+5MHgxg3sOTuGXv2yuslQLu9+czGGEg5IeNApZGMYtOYSsgaylaeq6UvsDTAgHeLYtSfHGBvYELx2SsgYXlaezqGAW/68RkWgOF4OnuTFEt/H4TWwH5WnTKuCLKEg/gFiYfO5jVk1mJJ3HL2EdZbG8XMsKH5/XSGAh5oyCwKpp57aSRxVAoSz/lsiX1MI2/BQJCsDBcsTeXzcKFY7JUpYf3zhnfiua/CvCfAzl7IIBxUMH2q3hHLDsshSWn2XdjFDwPCfHN2Qpe66dYfBqrgBytOWe8twsT7sZ/HcBAyLPVMYo3T7P1OFcoM1svsY5WMKEb/00A72b3ZlKLaZIGO4XyUlQHH2SxS6MEEYMEkQk5G+pNernAZxfwjlhGGlh4S+rna7efjQKW48Ix+VmRZnyyOxPkKIh7QQyE/Fsdo1h1xqjT850PZCn/zWT7uZHBakE0AULwRChTWg+/3sBOoVxc09FtV9iFfJyv3SIMhPwLcyP3dGHWZWCnUBZ0evi/DOMqASwQEA4FgXeimUWncOGYXLx5mp0fpvASzGlKWBsF4a3+zPuZRk0t3/lA3HvvrHFcV6aPK3YHbzHRn3FSwTdX8HZQ+s6W0sR8dm4fAUUfAWVFzgIcyfRQ5p10nCaRuIIa+DibXSaMBQJCszpGsTQVF45J32spxmVRCkcV3/loAiukUCyOUHx7hb1SiWNDUrbyjPGZHoyfA3YHWzDEm4S5ko+zsVMoZYcL6OVKeKaHsEKPsHIjZx42MKePYjmesiFdlyro99fYheFCWSAgQP8ZwKw+a6zEhWMSRQEWnTK+G8NYCSzyCCw78vZyX+ZgPnumBDuF0vTmafaVMIU7Hw/gFoueLmRsN+b/zuHwqDR9f401sPBwkODijuAyJGcOKngzQvHmaWwFJOh0MU3SUkEtEBCmVf2Zj7PZQlw4JjlGCstS2Xf5Pl+7RVgthWVWT+ZSBfx8sbKsrIzvvCDzMBqNhYWFr580LItk7Hh6ALeIdLEnM0KZlafrCwsL+c4LMhutVvvRudqu9jC6iwDjIGC9FBYVA+GpGx94fZurndWkkXGfbVzDd45Qp+Tk5NwzZVqVnXd1WdEXh/YA+PKdIxG4j00fPX3mDi9fT1p5/I+9bm5ufOcIdZzBYLhn8qOZhboSrWbV8qUAU/jOUQuwRygsNTU1h/d+Z1h0rOilQz+eyL5y5QrfOUKd8sqyd29Ofr/82T2G8W+8s3Yj39kRh9eXrWCf/a7k2R8vRcxc/+GnfGcHdcq+ffvSSEDRc7+wCw9/uG4139lpGQZCYTEYDERlA4QBALBxqKur4ztHqFPq6hvA2gEAqJV9dS3+b7ZLfX09WNkDAGvlUIuFJnL19fVGK3sAAKWV0SjQVfE4NCosTk5O9w0dsPfzKaXULsIXevfuzXeOUKesWvTSmOlPNQQOcruRvPjHbXxnRxzefXPBk/OmlPpEu+alvLT/R76zgzpl/PjxXd/bdOnrfLfySy+/8C++s9MyQilXi/VLSkpycnKCg4O9vb1bTNDQ0JCWlubk5NSsub969apGowkKClKr1Y1vGo3GtLQ0a2vrsLAwQlqdbp07d25oaOicOXPak0OdTufo6Ni+v8aiLly4MP7nmp1PRkR5CG5iWbCFJlj37S6OqD2/aGK0g4MD33kRjeLi4gU/Znn2iFozBAvtLgizen59seHL/We2TurStWtXvvPSMq56hLt27Xr22WcjIiLS09PffvvtZ599tlmC69evjxw50tvbW6PR9OvXb8eOHUql0mg0Tp069eTJk2FhYSkpKXPmzFmyZAkAFBcXjxw50srKqqamxtvb+5dffrG1FcxxrRzo0aPHuBJjQgEVYCBEd4WlkFLtvGEMRsG74+Hh8eiI/isycdOlFBzRKB6Mj+naVbgzcZzkTK/Xz5kzZ8uWLQcOHPjjjz8WLFhQUVHRLM2qVatGjRp1/Pjxc+fOZWVl7d27FwASExMTExPPnz//008/JSQkLFu2zLSLYO3ataGhoadPnz579mxNTc3WrVu5yLagxPuQxHyBjqej9ksvoV42RG2DhyTctQEebFYZLW/gOx+o0xLyabyPoO/pOQmEx48fZ1l2/PjxABAVFRUcHPzrr782S7Njx45nnnkGAOzs7B599NGdO3cCACHE1tbW1Ntzc3NTKBQsy5oSP/XUUwCgUqmmT59uSixt8b7McQ3VYygUuYQCOtJX0E2AYFkxMNCLHCvEewhxy62iVQbaW9iPW+FkaPTGjRsBAQGNM3mBgYG5ublNE1RUVFRWVgYEBJheBgQEHDx4EABGjhw5efLkiRMnRkZGHjp0aNOmTe7u7pTSvLy8polv3LjR2kdXV1dfuHDhwIEDppdKpXL48OF3mFMULDdrCHQkp4tprJf4Mo8aJeazMwV2vrCIxPswiQXsfX54OquIHcynI30ZgbdinATCmpoaKyurxpfW1tY1NTXNEgBAYxobG5vq6moA0Ol0OTk5zs7OdnZ2SqXy3LlzlFKDwdDQ0HB74hbl5eVdvXo1KyvL9NLW1rZHjx6tzR5XV1cLOUYO9VT+/qchzM7Ad0ZuIfBCExQDC8cKrT6MrsdC64Dq6upBror5p5VVffC8tfYS4Ddtf64qzpOtquLt5EgbGxulso1Ix0kg9Pb2LikpaXxZUlLi4+PTNIGXl5dCoSgtLTWdGVFcXGxK8Pnnnzc0NJhGPl9++eWuXbs+/vjjQ4cOdXd3Ly0tNf1uY+IWde/efcKECe1cNUopFfIShnv96MZM43IHG74zcguBF5qgJGtpkJPR391ep2Ox0O4WpXSYl931w/o6pYOHsCqBcAmweh4vNqwcYOUg7EePcTJoExUVdfny5aKiIgCor68/efLkgAEDmiZQKBRRUVFHjx41vTx27JgpQVVVlaurq+lNa2tre3t7nU4HAAMGDDh27FizxJI33IekaGkdHsEtWon5OEHYKUoGBqvJkUKcKherSxWUUghxEnot4KRH2K1btylTpkyfPv2ll17aunVrVFRUVFQUAHz55ZdfffXV4cOHAeDll19++eWXnZyc/vzzzz/++GPt2rUAMGnSpNWrV69duzYyMvKHH34wGo2DBw8GgPnz5z/66KNdunSprKz8+uuvk5OTuci20DiqoLcrSdHS4cJecIVak1jAvtQH57c6Jd6XScynDwbwnQ/UIWJZLMbVNP4XX3wxdOjQzZs3BwUF7dmzx/Rmr169HnjgAdPPjz322MaNG3fs2JGTk3Po0CHTRsvw8PCjR49eu3bto48+cnR0TE5OdnZ2BoDRo0f/73//+/nnn1NSUvbt2yef81bifUhiAd4Oi1IDCylaOsRbBK2AkMX7kIR8XDgqVon5NF4MgZDDk2V4IY2TZRr9kUdXnTEeuU9AJ+EJv9AE4nABXXjKeGKiErDQOsRUaCwF9Tf6jCkqbykfoWE2gvqmUQCfb/Qpk5T+wp4gBDx0W+CGqMmZElotrHWjqF0SC1iBbyIWBYbAEDWDh0uIUWYZdVAR4UdBwEAocHZKiHQnSRpJ9dplIiGfxvti/TKDeF+SWIBVQHxEtFgMK6rQjfTFs9bEp8YAZ0roYLU4WgGBG+lLEnGaUIQSC4R+slojDIRCF+/DJODtsNgc19BId2IvoLldEevjSnR6er0Ka4GYsBSOFLBiWfGOgVDoYtUEjx4WncR8ViyDQsJHAEb4MIfwdlBU0kuo2pb42omjFmAgFDorBgZ44tHDIpNQQON9sHKZTTyOjopNQoE4Nk6YYF0VgXhfBncTiohOD1lldCCelm4+uJtQdBLzxbRqGgOhCOBiAXE5XEAHeBIbPFLGfLo7EwpwuRJrgTgYWEjS0OHiGRQRTUblLMaDXNXR4jq+84HaJ7GAxY0TZhfvg7eDonGqmAY4EhEdlY7VVQTw6GFxEdH2KRHB3YQiIroqgIFQHExHD/OdC9S20nq4UkmjPcTUCojCKF9yMJ/FOiAKiQWsuBaLiSmvcoaLBcTiUAE7xJuosGKZm58DcVCS7HKsBUInxuPmsb6KQ6Q70dbRQnxSt+Al4slqnIn3xdtBEUjW0N6uxMWK73zcDayx4oBHD4tFQr5ojpUSHVwvIwpiPG4eA6Fo4GIB4dPWQkEtjXAXWSsgFiN9mcMFLM4TCpwYj5sXWXblDHcTCt/BfHa4N6PAOMgNHzvwtCVnS7EWCFeNAdJFeNw8BkLRwKOHhS9RVMdKiRHeDgrcMQ2NEOFx8xgIRQOPHhY+0W2fEp14H4LHDQqZSI+bx0AoJnj0sJDdqKaVetrHVXytgIjE+zJHC6kBQ6FQJYrzuHnx5VjOcDehkCXk03gfBsMgp9ytIcCBpBZjLRAi8R43j4FQTPDoYSFLzMcJQkuI9yX4qGphOlxAB3qJ8rh5DIQig51CwTpUgDsILSHeh+CGWmES3clqjUSZaTnD3YTCdLmS6lno7oyBkHPDfZgTWlpv5Dsf6DYJoh0UwUAoMqN8SQIePSw8ifl0lDibANFxtoKeLuRkEVYCYSmth6uiPW4eA6HI4NHDwoQ7CC1pJB46KjyiPm5enLmWN2wFhIYCJOaL73xF8Yr3YXA3odCI+rh5seZbznA3odBkl1M7JQlwxEBoIUO9SWoxrTHwnQ/UREI+HSnae0EMhOIT78McwqOHhUS8awREyk4J4W4kSYN1QChMx82Hi/a4eQyE4uNjB1549LCQJOKjlyxupC+etSYgYj9uHgOhKOHRw8LBUjhSyI7AQGhZ8b4MVgHhEPtiMQyEooRHDwvH2VLqYUO62Iu4FRCjOC9yvoxW6vnOBwIA8R83j4FQlExHD+sxFAqA2JsAkbJWQIwnOVaInUL+SeC4eQyEomQ6ejgNjx4WgMQC3DjBj3gfBs9aEwIJHDePgVCs8OhhITBSOK6hw8V5vqLYYRUQCAkcN48VWKzw6GEhOF1E/eyJpw3f+ZClgZ7kSiUtqec7H7IngePmMRCKFR49LAQJIl8sJ2pKBmK9yBFcNcYraRw3j4FQrExHD6fg0cO8wpPV+BXvy+DDWPiVIInj5jEQihjuJuRXAwspWjoMJwj5g1WAd2LfQWiCdVjE8OhhfqVoaQ8X4mLFdz5kLNKd5NdQTS3f+ZArCnBIEoMiGAhFDI8e5lcC7iDkm4LAUG/mEN4O8kQyx81zFQjr6+tfeuml4ODg6Ojon376qcU0W7duDQ8PDw0NXbJkCcuyAHDz5s17bvXrr78CwIkTJ5q+efLkSY6yLS52Sohwx6OHeZNYwMbjuCjf4n1wdJQ3krkXVHJ03VWrVqWlpR06dCgjI+Oxxx47e/ZsYGBg0wQnT56cN2/e3r17fXx8Jk+e7OPj88ILL7i5uS1cuNCU4ObNmzNnzvzss88AoKioqKCgYMOGDaZ/CgoK4ijbomM6a210FwXfGZGdOiOkFtM4tRRaAVGL9yUfZGOPkB+J+XRKoBSqACf3s5TSzz77bMWKFd26dRs/fvz48eO//PLLZmk+/fTTJ554YsiQIcHBwW+88cYnn3wCAHZ2dqP/dvPmzZEjRwYEBJjSu7i4NP6Th4yjW+oAACAASURBVIcHF9kWo3hfBh/Sy4tjhTTcjTiq+M6H7IW5kcoGerMaa4GlmY6bH+6NgbAVxcXFGo2mf//+ppdRUVFZWVnN0mRmZkZFRTUmyM7ONo2OmlBKv/rqq2eeeabxnZycnEGDBo0fP37Lli2U4pf+L3FeJBOPHuZDYgErjUEhsSMAw7xxEwUPpHTcPCdDo8XFxYQQR0dH00sXFxetVtssTUlJiZOTk+lnZ2dnvV5fUVHh6upqeicxMbGkpGTy5MmmlyEhIR999FFwcHBOTs6CBQvKysrmzZvX4kenp6d//vnnS5cuNb20srI6c+ZM4wc1U11dTYjo/xej3FT7r9WM8bXQ6JA0Cq3zDtxQLQs3VlW1q9ix0Dqg/YUW56744zp5wBuXjVn0m7bvmmKoJ6mqqrPMx3WYjY2NUtlGpOMkELq6ulJKq6qqTBGosrLS3d399jQ6nc70s06nUygUTcPV5s2bp0+fbmtra3rZq1evXr16AUBUVFR9ff2GDRtaC4T9+vWbMGHCrFmzTC+bXbYZSqmDg0MH/0jBGN2VPVGmnNLdQtOE0ii0TqrSQ1alPt7P2rZ9FQgLrQPaX2jjAun6HKODA550Z9FvWlKp4clQRhrFzsnQqKenp5OTU05Ojulldnb27ctbTN27xgSBgYEKxV9NeUVFxe7du5uOizbl4uJSW9vqviGGYezs7Fz/docoKBkj8ehhiztSSAd4knZGQcS1Hs7ESOGqDmuB5RhYOFZIR0hl1TQnf4ZCoZg2bdqaNWsMBsOlS5d27979xBNPAIBGo3nxxRdNYezJJ5/cunVrQUGBqYdnSmDyzTffhISEREZGNr5z9OjR4uJiALh+/fq77747duxYLrItUgPw6GGLS8zHjRPCMsKH4KoxS0otpv4OxEMKvUEA7vYRrlq1qqqqysPDY9CgQcuWLQsPDwcAnU63Z88evV4PAGPGjHnqqad69uypVqu9vLwWLFjQ+LuJiYlz585terX9+/cHBQXZ29tHRkZGRUW98847HGVbjJQMxKnx6GGLwrO2hQZ3E1pYQoFEdhCaEE5XYNbX11tbW98hAcuyBoPByqpdp1TV1dXZ2LRxBzJ37tzQ0NA5c+a054I6na5xRY+orTnH3qimm2ItMU0omULrsPIG8PtOXzxDZdXu20gstA64q0K7pqODfzLkPa6STtvcIRb7pt37m2FOH+Z+P4mMi3D7Z9w5CgIAwzDtjIIA0GYUlK14XxwXspxDBexgNWl/FEQWEOhIrBXkQjnWAkswHTc/1Fs6dUA6f4mcRbqTghpaiEcPW0RiPo33xYojOCNxmtBSTkjuuHmsz1JgOnr4ME4TWoRkzleUmHhfgtvqLSNRclUAA6FE4GIBy9DWws1qGukuqVZAGkb6ksR8lsVKwD3pHTcvqT9GznA3oWUkFrDDfBgFxkHh8bUj7jYkowxrAbdqDJBaTIdI4ojRRhgIJaKvG6lsoLlV2ApwKzGfSuAxpFKFD6y3gCQNjXAn9tI6TQIDoUSYjh4+hJ1CjiVKa/uUxMT74DQh5xILpPBI+mYwEEoHLhbgWn4NLamjfV2l1gpIRrwvc6SANWIl4FKCFFdNS+3vkbORuJuQYwfz6UhfhsE4KFSeNtDNgaQVYy3gik4PmWV0kKfU6gAGQukwHT18pRJbAa4k5uPJakIXj7sJuXSkkMZI8bh5DISSgnMknEoswJUyQhfvSxJxQy1npHrcvAT/JDmLx1VznLmqo/VG2tMFA6GgjfBhkjW0AUMhN6S6WAwDoaSM8iUH83FLMScS8ulIya0RkB4XKwh1Jie1WAnMr7wBrlTSGMlNEAIGQonxdyB2SpKDRw9zAHcQisVIXD7NjcR8Nk5NVFIMGlL8m+QNFwtw5HChNAeFpCfeh0nMx7FR80ssoJKcIAQMhNKDuwm5kFNOFQQCHTEQisBQb3KqmNYa+M6H5CRId9U0BkKpGeXLHMKjh80tIZ+OlmgTID0OKujnRpJxmtCstLWQJ93j5jEQSo2PHXjg0cPmllgg2XthSYr3wU0UZpZYwA6X7nHzGAglCI+YMS8KcLiAHYErZcQj3pfBKmBe0l4shoFQgnA3oXmdK6Wu1qSrvWRbAekZrCYZpVSn5zsfEpIg0R2EJhgIJWikL3O0EI8eNhvpPY9b8mwU0N+DHNdgHTCP/BpaVk/7SPe4eQyEEuRuDd0cSCoePWwmeLKaGMX74iYKszmQR+N9pHzcPAZCacInlJqLkcLRQna4RLdPSdhI3FBrPpJfLIbVW5pw1Zy5pBXTrvZEbct3PtBdGuBFLlbQsnq+8yEJhyQ9QQgYCKVqhA+TpKH1Rr7zIX4Jkl4sJ2FWDAzyIkcK8Xaws65U0noj7eEs5VqAgVCanK2ghzM5VYRDQ52VWMBKe1BIwuJ9GTxlqfMSC+goqR83L/E/T87ifUkCtgKdo2chWUOHe2M1ESWcKTcLOTyPGmu4ZOHRw52XoqWhzsTVmu98oA7p70FuVFNNLd/5ELlDMlg1jYFQsoZ6k1Q8erhzpPoYUplQEBiiZnCasDOyy6mKkf5x8xgIJctBBWFuJAmPHu6ExHxWqs+dkQk8ZamTEvLpKBncC2Ill7KRvgRHRzuszgini+lgb+m3AhKGj+fsJMnvIDTBQChl8T64aq7jkjQ0zI04qfjOB+qEcHdSWk/zqrEWdAQFOFLASn6CEDAQSlscHj3cCYnyaAKkjQAMw9vBjjpXSt2sSRcZHDePgVDKbBQQ7UGOFWIr0BEJ+TRe6tun5CDeB6cJOyhBNsfNYz2XuHhfBs9a64BqA2SU0jgvWbQC0jbSlxzEQNghcthBaIKBUOKsM/d9PG/qk8/P02q1fOdFNLbt3DVo/EOOP75eX1XOd15QZ9lV5Go/nz34/seOHT/Od15E4+LFi5Omzdy34gnXwnS+82IJbQTC6urqixcvVlRUWCY3yLyys7P/85/VuvtW/I8ZMnHaTL6zIw4pKSmzV396Pv6dQre+j/1rNt/ZQZ1175RpteEPJMW8PmXW/MLCQr6zIwIsy45+YOrebk/oR8+f8cy/dDod3zniXKuB8KuvvurevbuDg0OPHj1cXFz8/Pw2bNhAKY4wiElaWpqu7yTwCmEjJ13PK+A7O+KQnHKqLPxR8AxkB007n32R7+ygTqmvry+vp9BrFPj2bgiNz8jI4DtHIqDRaBocfSAkDrqFG/xjLly4wHeOOKds8d01a9a89tprAwcOXLNmjVqtLikp+eWXX+bPn5+bm7tu3ToLZxF1WExMjNO7TxWHDCP5mSH+3fjOjjgMHRzruvnVUv9I1eVjkf368J0d1CnW1tYedsqisz9RJ7X1xYPh4a/wnSMR8Pb2tqnWQOZ+sLZTXT/ds+davnPEPXqb+vp6JyenZ599ttn7b731lkKh0Gq1t/+KcMyZM2fTpk3tTFxZWclpZoQgIfHQwMlP+j78emlpqVkuKIdCW/2/X12Gz1iweIVOpzPLBeVQaGZnrkLLz89/+sWXVcOf+SMp1SwXFDJzFdq1a9ds7pk9YfqszMxMs1xQ4FroERYVFVVWVs6e3Xx25IUXXliyZMn169c9PT3bE1937tyZlpbWq1evadOmKZUtfFBqauqPP/7o4ODw5JNPqtVqAMjJyTly5EjTNI888oiLiwsAZGVl7dy5U6VSTZ8+3c/Pr/2RXubiRwzfP3i477d6e2fcGd5etO+Yp1ffu2aQgu+MIDPw8fHZ/P7awt8N1bgZpt0Mbv6eT2/4eWrLQ4bS08I3w8HBgRBSUlLS7P3i4mKFQuHv79+e6y5YsGDVqlVqtfqTTz558sknb0+QkJAwatQoa2vry5cvR0dHl5WVmT4i9W87d+585ZVXTBE0LS0tNjaWZdmioqL+/fvfvHnzrv9QGXNUQbAjSS/B+d32StbSWLUsVo3LR6wXk6zBKtBeSRoaJ6sq0GI/cfTo0X379r18+XLjO3l5ecOHD587d257upnFxcW2tramXy8tLbW1tb106dLtH7Fu3TrTz/fee+/atWubJZg1a9bTTz9t+nnq1Kmvv/666edp06Y1/nw7HBpt0fPHDBsyjGa5lOQLjaXU8+uGG1WsGa8p+ULjgnkL7UAeO+QnvRkvKEzmKrTnzNdiiELLYwVffvllQ0ND9+7do6Kixo8fP2DAgMDAwNTU1Pz8/Ef+lpCQ0FpwPXHiRNeuXYODgwHA1dU1JiYmMTGxaQKWZQ8fPjxmzBjTyzFjxhw6dKhpgtra2h07dsyc+deK/8TExLFjxzYmbnY11KZYNUnGx1C0z+UKaqskXWVwrJSsDPIi6SW0Ac+WaJ9kmfUIWx4C7tq1a2pq6tdff33o0KGbN2+yLBsWFgYA165da0xTWVnZ2kULCgq8vLwaX6rV6oKCW9buFxcX6/X6xjRqtTo/P79pgh07dnh5ecXFxQGAwWDQarWNE5O3J24qNzf31KlTZ8+eNb20tbVdvny5vb19i4nr6upUKlnMnEU5w6JCRV2dGU4dlXyhHc4jAz1IXV2dGa8p+ULjgnkLTQEQ5MCk5NfHeEj5jtAshabTwxWdood9vVkrAW9UKpVC0cZ8f6tzoQ4ODs8///zzzz/fgQ9WKpUs+8+tl9FobLZYxvSyMY3BYGj2n7d58+aZM2cSQgCAYRiFQnGHxE3Z29v7+fn179+/8YNsbW1bKwWFQtFmAUlDdxdgKRTWK7rYdfZSki+0k8UQ6wXm/RslX2hcMHuhxanhZAkZpDbjJQXHLIWWqoFIN7BVSeQba4ojd8bJoiAfH5+8vLzGl/n5+ffff3/TBK6urjY2Nvn5+aZOYX5+vo+PT+O/Xrp0KTk5edu2baaXDMOYeoF9+vS5PXEz7u7uoaGh7YzfKpVKPvfpA72MJ0vII86dXTgn+UJLKTbM6q1Qqcw5LiT5QuOC2QttsDf7Uy5VSaV9b5FZCu1kKRvnLfGCaoaT9cRDhgypqKg4ffo0AOTm5qanp5umAwsLC02DloSQCRMm7Nq1CwBYlt2zZ899993X+OubN28eN25c02h33333/fDDD6aff/jhh6aJUTvFeuE0Ydsq9XBVR/u5yWh2RD5ivfBJLO2SrGFjZXbcPCc9Qnt7++XLl0+aNGnSpEm///77vHnzTFFt9+7dn3zySXp6OgAsXrx41KhRN2/ezM3NNRgMU6dONf2uwWDYunXrRx991PSCr776alxcXEVFhU6nu3Llyv/+9z8usi1tcWqyIAWXCrQhRUv7exAr3G8mRcFOhKX0RjXthiuhWkcBUoroV8PlVQe42i85b9684cOHnzlzZsaMGbGxsaY3J0+ePHDgQNPPERERmZmZ+/fvd3R0HDdunLW1ten9+vr6LVu2jBgxounVgoKCsrKy9u3bZ2VlNXbsWAcHB46yLWHRHiSzjNYYwE4ue2Q7IklD5XYvLCuDvJgkDX00CP+LW5VVRt2tidqW73xYFoeNYmRkZGRkZNN3fHx8mg54ent7z5gxo9lv2dvbjx49+varubm5Pf7441zkUyZsldDHlaQW06He2Aq0KlnLvtBLXvfCsmLaR/RoEN/5EDDZbaUHAHweoazE4W7CO6IAJ4torBorhWTFeZEkPF/mjpK1chwUwTovI7FeBE+ZuoPMMuphQzxt+M4H4ky0J8kqozUGvvMhYMla7BEiSRusJklaXC/TqmQNjZPfvbCs2CigrxtJLcbbwZaV1UNBDe3jKrtagIFQRrrYE2uGXKnEVqBleNa2HODo6B0kaekAT6KQXyXAQCgvsWqShNOErUjCHqEM4Lm7dyDDHYQmGAjlBacJW1NSD4W1tLf8BoXkZrCaHNewWAdalKSR6WIxOf7NcoYLR1uTrJHpoJDc+NoRWwVOELTASCG1mA70lGMdwEAoL5Hu5HIlrTTDUyikJlnLynCxnDzFqXGasAVnS2g3B+JqzXc++ICBUF5UDES6k1NF2Ao0l6ShsV5YHWQBz91tkTw3TphgzZedWFw1dxsDC6nFdIAsB4VkKE6NM+UtkOdWehMMhLITqybJuJvwVmdLqb9cB4VkKMKdXNHRiga+8yEwcj5oFwOh7MR5MSe0FJfNNSXnQSEZUjEQhRMEt9LUQmUD7eEi01qAgVB2vGzBzZrkVGAr8I9kGd8LyxNuqG3muIYd5NWOR7lLFAZCOcLDNZpJwjNlZCbWiyRrcILgH8ly3UFoIt+/XM5icbFAEwU1oGug3Z0xEMrIYDWTUoQTBP9IkvfsAAZCOcJt9U0ladhYtXwHheTJwwbcrUl2OdYCAIB6I5wtoTEe8q0EGAjlqK8rya+hxXV850MYkrW4g1COcFt9o7QS2sOFOKj4zgd/sP7LkYJAjCc5iavmAECuj+RGuK2+ET6ADAOhTMV64W5CAIB6I5wrpdEyHhSSLZwgaIQPIMNAKFOxXgyOCwFAajHtKe9BIdnq60oKcIIAAOR9powJBkKZilWTU0XUIPs+IW6lly2GwABPkiL7CYLrVdTA0kBHWdcCDIQy5WIF3exJRpncWwG8F5Yz3E0IAEkaOljGOwhN5P73yxmumgOAExgIZSxWjRMEOEEIgIFQzmJlv1jgTx2lFALkPSgkZ4O8yOliuU8QJMl+yShgIJSzOC+5ny+DJ6vJnIsV+DmQc6XyrQW1BrhQQaNkv2oaA6F89XAhlXpaUMN3PviDZ22jOHmfvp1SRMNciY2C73zwDQOhfBGAgZ7khIx3E8r8fEUEf62XkW8gxFXTJhgIZS1Wzch2mrDaABfKaZQ7tgKyJvNt9TgoYoKBUNZiZfw8plNFNNydWMt+UEjmujsTnZ7m18ixFlCAE1p2EAZCDIQyN9CLpJfQeiPf+eADLpZDAEAABnmRE7LsFF6qoPYq0sUeawEGQnmzV0J3Z3KmRI6tQLKWxSWjCABivRh5ThMm4bjo3zAQyp0850gowAktHeiJrQCCWLkuHMVjlRphIJQ7ea6au1hBHXBQCAEAwABPclaWEwS4ZLQRBkK5i1WT4/ILhDhBiBrZK6GHC0mT2QRBpR7+1NF+blgLADAQoiBHQoHmVsmrFcDzFVFTMjxl6YSW9vcgKowAAICBEAHAIPk9mxAfyY2akuG5u0kaFicIG2EgRKan1cuoFajUw/UqHBRC/5DhhtpkDQ6K/AMDIZLd85iSNTTagyjxu4/+FuhICIHrspkgYCmcLKKDvLAO/AULAkG0B8mpoDUGvvNhKbiDEN1uoKeMbgczy6iXLfG04TsfgsFhIDQajZmZmQUFBXdIk5ubm5OTQ2nz719BQUFGRkZVVZXppV6vL2tCr9dzlWlZslZAX1dyqkgurUCShsbivTC6laymCfG4+Wa4ag6uXbvWq1evadOmRUREvPjii7cnMBgMjz76aFxc3IMPPti/f//i4mLT+1VVVZMnT+7Tp8+MGTP8/PxMcXTfvn1eXl7Bfzt06BBH2ZYt+Wyr/3tQCFsBdIs4OU0T4lnbzXAVCBcvXjxmzJj09PSsrKwffvjhyJEjzRLs2rUrIyPjwoULWVlZISEhq1evNr0/Z84cQkhBQUF6evrNmzfd3d1N7w8cOLD0b/fccw9H2ZYt+ayXOV9GvW2JBw4KoVv19yA55bRKHoNNuJW+GU4CoV6v37Vr16xZswDA3d19ypQp27dvb5Zm27Zt06ZNs7e3B4BZs2Zt27YNACorK7/55ps1a9bo9fr6+no7OzsrK6vGXykqKmpoaOAiwyhOTZI0rBwiYZIGmwDUAmsF9HMjp4ulXwmK60BTS3u5YC34ByeBsLCwsL6+Pjg42PQyKCjo+vXrzdJcv349KCioMUF+fn5DQ8OVK1esra3feeedqKgoX1/fWbNmGY1/HXx0+vTpyMhIZ2fnhx56qKysrLWPrquru3HjRurf0tPTb5+ARLfztSP2SnK5QvplhecrotbIZIIgWcsO8iIKrARNKLm4aHV1NQBYW1ubXtrZ2el0utvT2Nj8NT5la2tLKa2pqSktLa2qqvLz87t48WJZWdmAAQO2bNnyzDPPDB48WKvVOjk5lZSUPPjggwsXLvz0009b/Ohr167l5OTs37+/8co7d+50cnJqMXHjYhwEADHuVgm5Dd4BbRy5KPZCO15oPTu4QaezaHsn9kLjheULLcKR+fZP5YtBIh52ak+hHb6hjHIGna7WAvkRAhsbG5VKdec0nARCtVoNAGVlZR4eHgBQUlLi7e19e5rGjl1paamNjY2zs7Mp2RNPPAEArq6ukyZNOnbs2DPPPOPm5mZK6e7u/sorr8yfP7+1j+7Vq9fEiRPnzJnTzqw6Ojre3d8mXUN92TNlymcd235SrXgLrbgOShv00V0cGIvfDou30Hhk4UIbGUDnpRocHB1F3Vlqs9BOlxveiFCI/K80M06GRl1dXYOCgpKTk00vk5OT+/fv3yxNVFRUY4KkpKSoqChCSHBwsKura9MA6ezs3OwXCwoKXFxcuMi2zMXJ4GE0SRp2oCexfBREouBrRxytyCVJTxAYWDhTjA8ga46THiEAzJ0799VXX7W1tT1//vzx48e/+OILALh8+fKoUaMyMjKcnJyef/75mJiYuLg4tVq9YsWKdevWAYCNjc3s2bPnzZv37rvvXrly5fvvvzctN924caOzs7O/v392dvbSpUvfffddjrItZxFu5LqOVjSAs1XbiUUqWUtj1biDELXKdNZad2fJxon0UurvSCRcxzuGw0CoUqnWrVvn6uqakJDg5eUFAA4ODuPHjzcN1/bo0ePnn3/+73//W1tbu3r16kceecT0i8uXL9+0adM777zj4eGxf//+iIgIAAgICNi5c6dGo/Hx8dmyZcuECRM4yracKRmI9CApRfTeLpJtBZI0dHEkBkLUKtM+oqe6850PzuADyFpEJLaocu7cuaGhoe2cI9TpdDhz09Qbp4zWCrIs6k6hQryFpmfB7Wv9zakqy98Oi7fQeMRLoaUW06cOGzOmcNVD4FqbhTY10Ti2K3kyFG8Hb4HFgf4RqybJWpbvXHAlvYQG4aAQuqNwN3K9ipbV850PzuADyFqEgRD9I9aLOaGlUt1Xj6dpoDYpGejvQU5JdFt9fg2tNtAQ6c6AdhgGQvQPDxtQ25Kscmm2AriVHrWH6ZQlvnPBiSQNjVPjoukWYCBEt5DwE0rxcDXUHrFeTLJEqwDeC7YGAyG6hVQfRpNfQ2sMNNgJWwHUhlg1SSmiRglWAtMDyLAKtAADIbpFnBeR5O3wcRwUQu3jbg1qW5JVJrVaUG+EjFIajVvpW4KBEN2ijyvR1NKiOr7zYW74ADbUfpI8Zel0Me3lQuzFujGEWxgI0S0YAgM8SYrkWgF8JDdqv1gpjovgquk7wECImpPebsJ6I2SW0WgPbAVQu0jyeUw4KHIHGAhRc7FejMQWjp4qor1ciB0OCqH26e1CiuqkNkGQrGVjsUfYCgyEqLlYNUktpgYJ9QlxUAjdFdMEwQkJjYtc01ECxN8Ba0HLMBCi5pxU4O9AzpZKp1OI26fQ3ZLYbkLcRHtnGAhRC+LUktpWn6xhMRCiuxIrrYWjeC94ZxgIUQtMD6PhOxfmcVVHGUL8cFAI3Y1BXiS1mOqlMjiKPcI7w0CIWiClVXNJGjoYmwB0l5xUEOgokQmCagNcrKCR7lgLWoWBELUg1JlU62letRRagWQtxcVyqAPipHLuboqWRrgTawXf+RAwDISoBQRgkBdzQhKdQnwkN+oYyZy7i6um24SBELVMGq1AtQEuV9JI3EqP7p5kzt3FxWJtwkCIWiaN5zGd0NIId2KFX3N090KcSY2B3hT5BAEFSCmig7ywDtwJlg5q2UBPklFG64x856NzknFcFHUUAYhVi36C4EI5dVIRHzu+8yFsGAhRy2yV0MOZpBWLuxXAY6VQZ0hgH1ESLhZrBwyEqFVifxgNBTihxUEh1HESmCDAs7bbA9sI1CqxP4wmp5y6WhNvW77zgURrgCc5L/IJAtxK3x4YCFGr4tTkuEbER2sk4b0w6hxbJfR0JqminSAob4Ab1TTMFWtBGzAQolb5OxAlQ67pxNoK4PmKqPNEfe5usobGeBIlNvNtwRJCdyLqxQI4KIQ6T9RVIFnLYhVoDwyE6E7EO01Y3gA3q2lfHBRCnROnJkminSBI1tJYXCzWDlhG6E7E+zCaZA0dgINCqNP8HIhKnBMELIVTRXSAJ94Ltg3bCXQnUe7kYgWt0vOdj7uHg0LIXES6iSKjjPraEQ8bvvMhBhgI0Z1YK6CfGzklwlVzSRocFELmIdJzd3HVdPthS4HaIMajh40UThfTgdgKIHMQaY8wWYNnyrQXBkLUhlg1SdaKbLFARin1tSNu1nznA0lCfw9yqVJ8EwT49KX2w0CI2jBYzSRrqLjuh3HjBDIjFQPhbuRkkZgqgbYWSuppT2esBe2CgRC1QW0LTlbkYoWYWgHcSo/MK05s04TJWnaQF2GwErQPBkLUNtHNkeCgEDKvWC+SLKrdhLiD8K5gSaG2iWvVnLYWSutpDxwUQuYTp2aStJQVTSXAJaN3BwMhalucqHqESVo2FgeFkFmpbcFVPBMEehbOlNABGAjbDQMhalu4O8mtomX1fOejfZJxByHigIhOWTpTQoMdiZOK73yIB7YXqG0KAtEeJEUkq+bwkdyICyI6dxdXTd8tDISoXeLU4lgsoGchvYTG4PmKyNxE9DymZLwXvEtcBUKNRjNjxoywsLApU6Zcu3bt9gQsy65atSoiImLYsGE///xz4/sNDQ1vv/12TExMZGTka6+91vj+xo0b+/fvHxsb+91333GUZ3QHsWpGFOtl0oppiBMOCiHz6+dGblbTUjFMECRraBxOEN4NJUfXnT59ur+//+7du7/44ouJEyeeO3eOkFv+Yz744IPvvvtu27Zt165dmzZtWkpKSs+ePQHgX//6V25u7rp161xcXM6fP29KvG3btvXr13///feVlZUPP/xwQEBAbGwsE6E6XwAAHHlJREFURzlHLRrkRR7XUqPgQ2ESbpxA3FAQiPYkKVo6rpugv2B51bSepcFOgs6k0HDSI7xw4cLRo0fXr18fEhKyatUqrVZ75MiRZmk+/PDD5cuXh4WFTZw48aGHHvrss88A4Ny5c7t37961a9fQoUPDwsKmTp1qSvzBBx8sXLgwOjp65MiR//73vz/66CMuso3uwN0afO1JZpnQI2EyrhpHnIkTw3GDxzU0DheL3SVOyuv8+fM9evRwdHQEAIVCERUVlZGR0TRBQ0PDhQsXYmJiTC9jYmJMCU6ePBkTE/Ptt99OmDDh3//+99WrVxsvGB0dbfo5Ojq62dWQZYhiWz1upUfcifUSwQQBThB2ACdDo0VFRS4uLo0vXV1dNRpN0wTFxcWUUmdnZ9NLFxcXrVYLADdu3Dhx4kSfPn1WrFixZ8+eIUOG5OTk2NjYlJeXNyZ2dXU1JW7RuXPnvv3223Xr1pleKpXKQ4cONf5uM9XV1c0GbNEdRDkrjuYzk12FW2g3qqHeaOVF6qqq+M7KrfCb1gECLLR+9uSkVlWhq1IIK1//qK6uPppv9U6ksapK6D1Xi7GxsVEq24h0nARCZ2fn6urqxpc6nc7V1bVpAlOYrK6uNr1fVVVl+sHJycnGxmb9+vVKpTI6Onrbtm0HDx584IEH7OzsGi/YmLhFvXv3HjFixBNPPGF6aWdn5+3t3VpiSqmDg0PH/06ZifejGy8Y7e3tBVto57TsEG8h/p/iN60DBFhoDgBd7A1/6u3D3QQaCWv0NKeSGdLN2o6r5R/SxMnQaGBg4NWrVw0Gg+nlxYsXAwMDmyaws7NTq9UXL15sTBAQEAAAQUFBjo6OjdHbxcWlqqrKdMFLly41S9wipVLp7u4e9Lc7REF0t3q7kpJ6qq0TaBMAeNY24p7AJwjOlDF9XAlGwbvFSSAcOHCgh4fHli1bAOC3334rKSkZN24cAKSkpLz//vumNNOmTdu0aRPLslqt9ttvv50+fToAjB8/vq6u7sCBAwCQmpqanZ0dFxcHANOnT//ggw/0en1lZeXmzZtNiZGFEYCBnuRUiXAjDZ6viLgWqxb0tvqUYgbnyDuAk0BICNmyZcvKlSv9/Pyefvrpr7/+2sbGBgAyMzN37NhhSrNkyZLq6mpvb+8ePXpMnz591KhRAGBtbf3NN98888wzISEh999//6effhocHAwAL730kpubm6+vr7+///Dhwx955BEuso3aFKTL/OyLzenp6XxnpLna2trNW785t29bHwcx7PNCotXXuvz37V/u2bOHZQU3CZeYmPj9N1+F6HP5zoj4EMrZI1dZli0tLXVxcbnDRGVZWZmNjY2trW3TN41GY3l5ubu7e7PEFRUVSqXS3t7+Dh86d+7c0NDQOXPmtCeHOp3OtLQVtcfv+w88Mm95ZfQM97PbPlv58gMT7+c7R39hWTZ88MiLXUYZjGy/oiNpRw8IbZEFftM6QICFVllZGRY38kbY4/YVufHO5Xu/+4rvHP1j2TtrNv2WWh4wzD3lsxO/7gwJCeE7R2LC4VgywzAeHh53TtPisheFQnF7FASA1hZ/Isv4cMu2yoc2QLfwkh7DN32xXDiB8MqVK1or74axiwCgYEvG9evX7zCLjFCHJSUllYfeS0fNqQJIWT/EYDC0uRzRYrZu31U+OxEYZamd87ZdexYvXMB3jsQE912i9grx76LKTQUA5vrpYP+ufGfnH56enqz2Cuhrob6aLbra4l0UQp3XpUsXVf45YA1QqbGiAoqCAKD28oSb5wHA7ubpID8BVU9RENB/JBK45YsWpE+feW7NJxWOfst//ZLv7PzDxcXlnqdf2bM23sOGWb3yTaGNpyHJCAsLm/vgyI/XDy4xWs95ZwPf2bnFe2vXjZ76nJu+ZOyo4Y89iqso7g6Hc4S8wDlCrul0ulmn7aLcyav9hDKcUGeEkB2GX8YoBLu7C79pHSDkQtt1jX33LHtqslI4X7iZR4x+DuTl0GrBFpqQCaUtQyKyNJL5vwyjTs93Pv72YRYb50UEGwWR9DwYyBhY+O2GUHoR16voT7ns3D7YnncQFhy6a71cyChf5qNsQSwfrzPCuvPs4kj8JiPLIQBLo5ilqUJ5HMuKNHZ2b4WrNd/5EC1sPlBHLI9i1gmjU/hBFjtYTfphdxBZ1gMBjIGFXwXQKbxSid3BzsKyQx3R3ZmM8mU+zOK5U1hrgPXn2Tcj8GuMLM3UKVwmgE7h2+nsnD7YHewUbEFQBy2PYtaf57lT+EE2dgcRbx4IYIwUfsnlMxRid9AssPhQB4U6k9FdmA/46xTWGGDtOeNi7A4inhCApZHMsjQ+O4Wr0tm5fRQuVvzlQBKwEUEdtyyS2cBfp/DDbHa4DxOG3UHEn8kBDAO8dQqvVNKfc9k52B3sNCxB1HGmTuF/+egUVhtg7Tkjzg4ifhGAN/nrFL51hn0Ju4PmgO0I6pRlkcxGPjqFH2axI3yxO4j4N8mfYQB+zrX07eCVSvrLDfZF7A6aAxYi6pRQZ3KPxTuF1QZYl2F8Ixy/vYh/BGBxJLMslbVwp3DlGXZeX+wOmgc2JaizlkVZeqbwgyw2HruDSDAm+TNKBn66brnbwcuV9Ncb7Oze2ICbB5Yj6qwQJzKmC/N+poVagWoDrM8wvoGzg0hIFkcwy9Ms1ylcmYbdQXPC1gSZwdIoZmOmsdIincL/ZrIjfZm+rtgdRAIy0Z9RMrDXIp3Cy5X09zxcLGpOWJTIDEKcyNiulugUVhtgw3kjniyKBGhxBLPCIp3ClWnsS30UTiruP0k2sEFB5rEkktlw3ljewO2nvJ/JjvRlerlgdxAJjqlT+CPHncJLFfQP7A6aG5YmMo8QJzKhG/NfLjuF1QbYeN64BLuDSKiWRHLeKVx5hn2pr8IRu4NmhW0KMpslkczGTA47he9nsqO6MD2xO4iE6n4/xorLTuGlCro/j30RF4uaGxYoMptgLjuFpu4gniyKBG5JpGLJaZajXuGKM+w87A5yAJsVZE7cdQo3ZbKjsTuIBO8+P2Kn5KRTeKmCHsjDvYOcwDJF5hTsRO7zM//y0So9bMLFokgklkYplqaav1O4PA27g1zBlgWZ2eIIZpO5O4WbMtl7ujA9nLE7iERgQjdip4Q9Zu0UXqqgB/OxO8gVLFZkZsFO5H4/ZpP5OoVVeng/0/gmdgeReCwzd6dweRo7H7uDnMHGBZnfmxHM+5nGsnrzXG1jJntvV+wOIjEZ3404qGD3n+a5HbxYQQ/msy9gd5AzWLLI/IKdyEQzdQqr9PDfTDxZFInPskjFsjTzdAqXp7Evh2F3kEPYviBOLI1iPsgyQ6dww3l2DHYHkQiN60YcVPBDpzuF2eX0YD77fC9sqzmEhYs44e9A7vdjNmYaO3ORSj1syjS+jt1BJE7LIhXLO90pXHmGXYDdQY5hE4O4sjSK+TCL7UyncON5dnw37A4isRrXjbhaw65OdAqzy2kCdge5h+WLuOLvQCb6d7xTWKmH93F2EIncmxGKFZ3oFK5IY1/tp3DA7iDHsJVBHFoS2fFO4cbz7AQ/pjt2B5GYje3a8U5hVjlNLGCf64mtNOewiBGH/B3IJH9mw/m77hSauoOvh+P3E4ne4gjF8g7tKcTuoMVgQ4O4tTiS+SibLb3LTuGG8+x92B1EkjCmK3Gzge+v3V2nMKucHsbuoKVgKSNu+TuQyQHMxrvpFFY0wH8zjYuwO4ikYknkXc8UYnfQkrCtQZxbHHF3ncIN59n7sTuIJOTeLsTdBna2u1OYWUaPFuJiUcvBgkac83MgDwS0d6awogE+yMLuIJKaJZGKle3uFK5IYxeEKeyUHOcJ/Q2bG2QJb7a7U7j+vHGiPxOK3UEkLfd0Ie42sONq253CzDJ6TMM+h91BC8KyRpbg50AeDGDWt9UprGiAj7JZ3DuIJGlJpGLlmbY7hcvT2Ff7YXfQojgs7LNnz/7+++8uLi5Tp051dHS8PUFxcfH27dtra2snTZoUGhpqevP333+vrKw0/ezp6TlixAgAyM/PP378eOMvDhkyxMfHh7ucIy68GcFE7ja81EfhYdNqmvXnjff7MUGO2B1EEnRPF+JpA9uvslODW73VyyyjSRq6ZbjCkhlDXN1679u3Lz4+vrKy8rfffouLi6urq2uWoKSkJCoqKiUlRavVxsTEpKammt6fP3/+Rx99tHPnzp07dx46dMj0Zmpq6gsvvLDzb4WFhRxlG3HHz4E8HHinmULsDiLJWxKpWJ7GGlvvFC5LYxf0Y7A7aGFclfeqVaveeeed5557jmXZmJiYnTt3zpgxo2mCzz//vG/fvlu3bgUAOzu71atX79y50/RPb7311uDBg5tdsEePHjt27OAot8gy3jB1CvsqPFvqFK7LME7E7iCStNFdiLct7GilU5hZRpM1dCt2By2Ok7vvmpqa48ePT5gwAQAYhhk/fvz+/fubpTlw4MD48eNNP0+YMOGPP/5o/KeEhIQtW7acPXu2afry8vLNmzfv2bOnoqKCizwjC/BzII8EtbynELuDSCaWRimWtdIpXJbGvordQT5wUuT5+fkAoFarTS+9vb2Tk5NvT+Pt7W362cfHp7KysqqqysHBISAg4PLly5cuXZo/f/5TTz21bt06AFCpVN7e3ikpKTk5Oc8999y+ffsiIiJa/Oi8vLzMzMzc3FzTS5VK9dprr9na2raYuL6+3srKqtN/rrx0stAW9IYBP5Hnuxs9rG9pCf5zltzfDXytGurN9Fx7QcFvWgdItdCGuIO3DfnfhfrHAm95P7MckgqZzwaxnakCUi20zlAqlQpFG51sTgIhIQQAKP2rpWNZ9vZ8EEKaJgAAhmEA4NdffzW9eeXKlb59+z755JPh4eFjx44dO3as6f158+YtXLjw999/b/GjlUqlnZ2dq6tr40ulUmm68u0Yhmntn1BrOllofg7wcCD8NwdWRv4zBFreAJ9dhOMTgGGkOS6K37QOkHChLQ6H2SfII4GgbPL3vXUWFvQFB6tO/ckSLrQOM8WjO+MkEPr6+hJCCgsL/f39AaCwsLCx89c0TeOal8LCQhcXFzs7u6YJgoODg4ODc3JywsPDm74/ZsyY3bt3t/bRarU6NDR0zpw57cmnSqVSqfAIo7vT+UJ7M4pG/GB4uZ+qcaZw0znjAwEQ6irZqRH8pnWAhAvtXj/okmH44QYzLeSvoHW+jJ4sNn47UqnqXJMs4ULjFCf3Dra2tsOGDdu7dy8AGI3GX375xdSfq6+vz8jIMPX/xo4da0oAAHv37jUlMP2TSW5u7tWrV0NCQpq9n5iY2LjXAolRN3vyaBCzPuOvmcLyBvgkm8WjZJCsLItSrDjDGv5u2JamsgvDGVucHeQJVwW/dOnSKVOmXL16NTs7W6lU/n97dx8UVb3GAfyHiyDKrLgLsS+8LOMtbpopSphatHZZA0QZ3sKGK/EioKIUadN47d6h0plg0KIpNaY79qLBZKgECSwoAWogapFJDijsgLuBxIvAAgu759w/zsxeBxBfRvcc9nw/fx2cR+c7zyzncc/vd84JDw8nhLS0tDz77LO9vb1OTk4JCQkHDx6MiIiQSqXffvttZWUlIeTy5ctJSUkrVqygabqgoCA2NnbZsmWEkISEhMHBQQ8Pjz/++OPChQulpaWPKTZYxq4lM5YcN6YvErjMIvuumMIVM7ywWRT4ZLXURjab5LdQ//zbjN976bpb9FGl1V4R4b7/L9Q9cs3NzWVlZc7OzqGhocx2lcHBwaqqqldeecXW1pYQ0t/ff+LEiaGhoZCQEHd3d0LI6Ojo2bNnr127Zmtru3TpUl9fX+afam9vr6mp6erqkkqlKpXKvAQ4UVpa2v1fGh0YGJj0Tn+YwqNqWup50xzasMmj/4UqcV2orXUPQnzSHoLVN+1sB51QbaxY2ZX2u/gf7jO3L3wEF0WsvmmPyWMchKzAIHzcHlXT/nusOHnHbnuRVGgvaKs5Yd1b3fBJewhW37Suri7PF9bZODoburXnT3zp57P43n/nXqy+aY8JFmaAHXv27KHePjO8rahfvuz43Xc/AVirzJwDI8rtQ1tOmBK+fPvfe9iOw2sYhMAOk8lEbO0IIZTd7JERa7x5EGBKI4ZRYjeHEELs54xY5f2z0wd2KQE7/pWe+p+cIFryd+eea5ER2P0EvLMzNelkULjhetkMTX3mof1sx+E1DEJgx+bEuPVBa3Q63eLFi3HnE/CQQqFounT2ypUrTz65VyQSsR2H1zAIgTUymUwmk7GdAoA1s2fPXr58OdspgN9rhGq1+s5b9eF+4CbOh4CmPSiTyTTxSf0wtZGREeaGbHhQvB6Eqampt27dYjvFNLNx40YDFvYfhF6vj4uLYzvFNPPnn3+mpaWxnWKaaWpq2r17N9sppiVeD0IAAAAMQgAA4DUMQgAA4DVr2zV6/fp1tVptfq/F1Pr6+qKjo6374V6P3NjYWHBwMN55dv8oijIYDCqViu0g04nBYOju7kbTHsjg4KBGo0HTxgkLC9u6devUNdb2rNH6+vqbN2/e59P2Wltbvby87l0Hd0DTHgKa9qBomtZoNGjaA6Eoqr29nXkLLJh5eXnNnz9/6hprG4QAAAAPBBe4AACA1zAIAQCA1zAIAQCA1zAIAQCA1wQZGRlsZ7CEnp6ekydPNjU1eXp6Tnq/hNFoVKvV586dE4lEc+fOtXxCDtLr9UVFRQ0NDXK53MHBYWLB9evXy8vLr127Nm/ePLwXm6HRaAoLCzs7O728vKa4yaSxsbGxsVGhUFgwGne1tLQUFhZ2dXV5eXnZ2NhMWnP58uXS0lKtVuvq6jpr1iwLJ+Sg5ubmoqKinp4ehUIxadO0Wm1JSUljY6NIJMKv5z3QPNDS0uLq6hodHR0cHOzt7d3T0zOuwGQyrVmz5rnnnktMTBSJRGfOnGElJ6f09vZ6e3sHBgZu2LBBIpHcuHFjXMFnn30mlUqjoqLCw8OFQmFRURErOTlFrVaLRKLExERfX9+goCCKoiYt0+l0Li4uYrHYwvG46ccffxSLxZs2bfLx8QkNDZ1YQFFUcnKyu7v7xo0b161bl52dbfmQXFNQUCAWi5OSkhYtWvTaa69NLCgtLXVyckpMTIyPj3dycsI5bWq8GIRbtmxJTk5mjgMDAzMzM8cVlJSUKBSKoaEhmqYPHDiwatUqS0fknqysLJVKxZzKt2zZkpKSMq6gra3NYDAwx9nZ2T4+PpaOyD3Lly8/dOgQTdNDQ0MeHh7l5eWTloWFhe3cuRODkLFkyZLDhw/TND0wMCCVSqurq8cVHDlyZP78+b29vSyE4ySKory9vfPz82ma7uvrE4vFFy9eHFcTEhKSkZHBHO/atSsiIsLSKacVXqwRFhUVRUZGMscRERHFxcXjCoqLi9euXctc/YuMjDx37lxPT4+lU3JMcXFxREQEc8klMjJyYtPc3d3NF5mlUileSXHr1q26urqIiAhCiIODQ3Bw8MSmEULy8vLs7e3Xr19v8YBc1N7e3tDQwDTN0dExMDBwYtPy8vI2b97c3d195syZ7u5uNmJyS3Nzc0tLS2hoKCFk7ty5KpVqYtPEYrFer2eOh4aGnJ2dLZ1yWrG2R6xNRFFUR0eHXC5nfpTL5VqtdlyNVqv18/Njjl1cXOzs7LRaLc/fGa3Vau9sWkdHh8lkEggEEytHRkaysrI2bdpk2YCco9Pp7O3tzWccuVze0NAwruavv/56//33f/rpp6amJosH5CKdTicUCs0rWHK5vLW1dVzNjRs39Hr9sWPH3Nzcqqqq8vPzAwICLJ6UQ3Q6nbOzs3mhdNJzWmZmZkxMTFBQkMlksrGxOXr0qMVjTifW/42QoiiKosyLyQKBwGg0jqsxmUx37muYtIZv7uyJQCCgadpkMk1aFhsbq1Ao8Pa4+/kUbdu27Z133nF1dbVsNO5iTtPmHydt2sjIiEAgqK2tLSgo+OCDD9544w3LZuSc+2naqVOntFrtq6++umHDhtbW1rKyMstmnGas/xuhra2ti4tLV1fX008/TQjp7OyUyWTjaqRSqfkNvf39/cPDwxNr+ObOnnR2djo7O0/cbUtRVHx8/O3bt3/44YdJvyzyikQiGR4eHhwcdHR0JIR0dnZKpdI7C9ra2goLC4VC4c8//9zR0aHX61NSUjIyMsaV8YpEIunv7zcYDPb29mSyphFCZDKZv78/c+pXKpXbtm0zGo22ttZ/7robiUTS09NjvkLT2dlpvnhj9u677+bm5q5du5YQMm/evJ07d8bExLCQdZqw/m+EhJDVq1er1WrmWK1WK5VK5ri7u5v5lqNUKpl9DYSQsrKyBQsW4P/sSqXS3LSysjJz0/r6+kZHRwkhNE1v3bpVo9EcP36cOYvxnEwme+qpp5imURRVUVGxevVqQojRaGRWtkQi0ddff61SqQICApYtWzZz5syAgIA5c+awnJtVCoXC09OzvLycEGIymU6fPm1umnmd/uWXX25ubmaOm5ubJRIJn6cgIcTb21skElVVVRFCxsbGKisrmaaNjY319vYyNQKBgPk9JYQYDAb8P/UeWN6sYxEXL14UCoUZGRlvvfWWWCxua2ujaZr5lNTX19M0bTAYFixYEB0dnZ2d7erqevToUbYjs6+trU0sFqenp2dkZAiFQqZRNE27ubl99913NE1//PHHNjY2MTExycnJycnJqamprOblhK+++koikWRnZ0dFRT3zzDOjo6M0TdfW1hJCxt1KUV1djV2jjNzcXLlcvm/fvrCwMB8fH6PRSNN0VVXVzJkzmQKdTieTyXbs2PHRRx+5ubl9/vnnrOblhJycHE9Pz/3794eEhKxYsYL5dJWUlAiFQqZg7969crk8Ozs7KytLIpHs27eP1bxcx5e3T1y9evXYsWN2dnYxMTHMa0pomv7iiy/CwsKY3Q19fX2HDx/u7u5es2aNv78/23k5oa2t7ciRI6Ojo1FRUQsXLmT+MC8v7/nnn/fy8qqvr//ll1/Mxba2tgkJCSwl5ZDKysqKioonnngiLi6OeTJDV1fXyZMnk5KS7izr6OhQq9WxsbEsxeSWioqKyspKV1fX+Ph4ZuNMR0fHqVOnzJ8onU73zTffGAyGgICAlStXshqWK0pLS6urq+VyeVxcHHNd4ebNm6dPn3799deZgsrKypqaGhsbG6VS+eKLL7Ialuv4MggBAAAmxYs1QgAAgLvBIAQAAF7DIAQAAF7DIAQAAF7DIAQAAF7DIAQAAF7DIAQAAF7DIAQAAF7DIAQAAF7DIAQAAF7DIASwHikpKR4eHhqNhvlxdHTU39/fz89veHiY1VwAnIZnjQJYj4GBAV9fXycnp5qaGjs7u/T09IMHD54/f37p0qVsRwPgLgxCAKty6dKlVatWpaWlvfTSS+vWrcvJydm+fTvboQA4DYMQwNp88sknb775pqOjY0BAQEFBAfNudwC4GwxCAGtz+/ZtDw+P/v7+3377bdGiRWzHAeA6bJYBsDabN2+eMWOGu7t7amqq0WhkOw4A12EQAliV3Nzc/Pz8AwcOfP/997W1te+99x7biQC4DpdGAazH1atX/fz84uPjP/30U0LIhx9+uHv37tLSUpVKxXY0AO7CIASwEnq93s/PTyAQ1NXVOTg4EEJoml6/fv2FCxd+/fVXqVTKdkAAjsIgBAAAXsMaIQAA8BoGIQAA8BoGIQAA8BoGIQAA8BoGIQAA8BoGIQAA8BoGIQAA8Nr/AOvB5RM9b5xlAAAAAElFTkSuQmCC", "text/html": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ], "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "metadata": {}, "execution_count": 6 } ], "cell_type": "code", "source": [ "rvecs = collect(r_vectors(basis))[:, 1, 1] # slice along the x axis\n", "x = [r[1] for r in rvecs] # only keep the x coordinate\n", "plot(x, scfres.ρ[1, :, 1, 1], label=\"\", xlabel=\"x\", ylabel=\"ρ\", marker=2)" ], "metadata": {}, "execution_count": 6 }, { "cell_type": "markdown", "source": [ "We can also perform various postprocessing steps:\n", "for instance compute a band structure" ], "metadata": {} }, { "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Computing bands along kpath:\n", " Γ -> X -> U and K -> Γ -> L -> W -> X\n", "\rDiagonalising Hamiltonian kblocks: 7%|█▏ | ETA: 0:00:01\u001b[K\rDiagonalising Hamiltonian kblocks: 100%|████████████████| Time: 0:00:01\u001b[K\n" ] }, { "output_type": "execute_result", "data": { "text/plain": "Plot{Plots.GRBackend() n=44}", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOydeVxM+xvHPzMt2mnPEtmSUFTIvmUX+lFZExGufS27bDfEFSVbdiG7XLtwrUnIlrJfxBXSvs3M9/dHttIyU+ec7zSd9+u+7mt8O+f7fM6ZmfPM812eR0AIAQ8PDw8PT3lFSFsADw8PDw8PTXhHyMPDw8NTruEdIQ8PDw9PuYZ3hDw8PDw85RreEfLw8PDwlGt4R8jDw8PDU67hHSEPDw8PT7mGd4Q8PDw8POUa3hHy8PDw8JRreEfIw8PDw1OuUUxHGBwc/PDhQ46N7t6N0aM5tsk1EomEVk6+e/fu5eTkANi+HffuUZGQh9OnceJEnpbMTMyYQUnNd5YvR3w8A/2IxeJS9iASYfp0BpTIP6W/V7k8e4aAAEZ6Yh2mLpkbpHlkKaYjPHv2LJeOUCQSzZ07d9y4xc+eveTMKBWysrJofQe6du365csXAKtWQSSiIiEP2tqYNAk5OT9bKlTAqVM4e5aeJuDNGwQHM9BPenp6KXuIiEB4OANK5J/S36tcDh/Gs2eM9MQ6TF0yB0gkEpEUzwvFdIQcU6dOy2XLPiYn61+61P7p06e05SgyCQl4+xZNmtDWAbRqBTMz7N79s0UggJcX/vyTnibA3R3btkEeEumfO4cuXWiLKFNcuoR27WiLKK/wjpAB3rxJIGQT8AchnkFBQbTlKDKXLqFNGygp0dYBAFi0CEuX5glPBw3C27e4do2aJDs7aGvjyhVqAn5w7hw6d6YtouwgEuH6dbRpQ1tHeYV3hIygCvwLiIErDRs2pC1Gkbl4ER060BbxnVatYGqKPXt+tigpYepU+PrS0wS4uWH7dpoCAKSk4MEDtGxJWUYZ4s4d1KgBAwPaOsorvCMsLSEhAAKBTkBNJaVnzZo1o61IkZErRwhgwQIsXpwnKBwxAnfu4M4dapKGDsWRI0hJoSYAQHg4WrSAujpNDWWLy5fRvj1tEdLx+PHjBQsWHDt2jLYQJuEdYanYtAkzZkBDIw54CvyrrDzm0aNHtEUpLP/9h48fYWVFW8cvtG+PatWwd+/PlgoVMHkyVq6kJsnICG3b4tAhagLAj4vKzuXLZWOC8PLly40adVuzRsfJadmgQYqzSp53hCUnMBB//onLl6Gi0l4gOA78ram5r3Xr1rR1KSzXr6NdOwjl7DO7YAGWLMGva2n/+APh4YiLoybJ3Z3y6CjvCGVCLMb16ygTT465c1dKJKsIWUjIxYMHz9CWwxhy9lApOyxfjnXrcOUKlJVBSH0dnTd2dtd79NhZtWpV2tIUlmvX5GtcNJcOHWBsnCco1NTE2LE0g8ImTd7euOFsamq/YMEK7q3/+y8SE9GoEfeWyyp376JaNRgZ0dYhBQYGxkDuJo9XqqqqlNUwB+8IS8KCBdixA+HhqFYN27ejRg3B5MnjDh9eevq0ZWYmbXGKy9Wr8ugIAcyfj0WL8gSFEyfi8GH8+y8dPYMHj8vOnvj27RV//0enT5/m0nRo6GFn5wlmZtsACZd2yzRlZeNEfDyePAkUCM4JBBbA4P37/6KtiDF4RygbhGDKFJw4gX/+QZUqIAS7duH1awwfDlNTNGmCI0doS1RQxGKkpKBBA9o6CsLBASYmCA392aKnh+HDsWYNHT1v3sQDbQCV5OQu0dFPOLMbGnpk9Og9t24NevAgaskSxXlKsk2ZmCCMiUHLllBSUlu4MPzOnSiB4G6bNj1pi2IM3hHKgEQCT09cu4Zz574tdL58GdnZaNkSNWoAwKhR2LyZrkbFZPfu3UlJKcbGYQIBbSmFMG8eFi+G5JcoaPp07NjxZd++U69eveJYTLdu7TU1pwHHNDXX9u7dlTO7R45c+Pp1JtAiM3Ph0aPnObNbphGLce0a2ralraNIIiLQsSM6d4aSEry9Ubu2REsLgYG0ZTEH7wilRSzGiBF4+hQXLkBP71vjjh1QVcWoUd/+2bs3Hj8Gn1uGWdatW+/mtkIkUo2JOdCvnzttOQXTuTP09XHgwM+W1NRnmZmd3N2vN2065OBBTteaBwX5BgbaDB16T119i5lZfc7stmnTRFU1BPhSocL2li1tOLNbpomORpUqcj1BGBYGR0f89Rf+/htbtyJ3ZrBFC+zfT1sZgxBFxMXFZe/evQx2mJVF+vUj3buT9PSfjampREeHmJiQnJyfjdOnE29vBi3LF+np6Tm/Xi0n1K3bBrgEGAOxamq1ObYuPadPE0tLIhZ/++eUKfOB4wABEho1cqAiydGRrFtXkhOTk5NLdJZYS8vX3LzzuHGz0n/9qig0JbtXP1i9mvzxB1NamGfbNmJsTK5fJ/36kTlzvjUmJycfOkSUlakqkw6xWJydnV3sYXxEWDxZWXB2RnY2jhzJs0c4NBSGhhgxAsrKPxtHjcL27XlyMfOUkvr1TYFzAIArBga6lNUUTteu0NX9uYFPT09HReUdAOCdrm5FKpIWLsSffyIjgyNzGzYIe/b0io09GxCwTJ3fTi8d8rxSZvlyLFqEf/7Bmzd49Ahz5/780//+BwBHj9KSxjQc+GTuYTAiTE0lDg5k4EDyeyDUti2pWJE8fZq/vU0bcuQII8blDioRYVpampFRc0CpUqV6sbGxHFuXiZMnSYMG34LClJQUO7uuBgYtlJWb7dnzmJakXr1IQIDMZ5UgyklNJSYm5P59mW2VdUoTEYrFRF+ffPjAoBxmEInI2LGkUSPy9i359IlUrkyuXfv519xLtrIiXbtSUyglfETIAElJ6NIFNWti9+48YR+AV69w9y7s7FCnTv6zRo7Eli2caVR8NDQ0Ro++qamp/+TJZXNzc9pyiqJ7d2hq4vBhANDS0oqMPP327cVduyJ8fevTGiTw8cGyZeBgV09gIDp04LcPykZ0NIyMYGxMW0desrIwaBBiY3H1KqpWxfjxGDKkgMyxrq64cYOGPhbgHWGhJCaiSxc0boyNGwvIZrJ9OypWLLgSr7MzIiKo7SFTSG7dyv9DRG6ZPx8LF/5cPlqhQoUBA1C9Ovz96eixsUGTJti6lV0raWn46y/Mns2uFUVCJBINGTKhQwe7xMQ+8YwUU2aIr1/RpQvEYvz9N3R0EBaGO3fg41PAkePHIyWFZgYlBuEdYcH89x/atUOHDggMxO9L9gnB1q1IT0fv3gWcq64OV1f6FQAUBkJw+zZUVGjrkI6ePaGhkX/uJCAAK1bgJaWyzT4+8PVFVhaLJtatQ8eO4CuvSE9w8M4jR3SSkm7/99+EkSO9aMv5xvv36NABjRsjNBRqakhKwrhx2LKl4PzpOjqoXBmrV3OukgV4R5iHnJycnJycN2/Qpg369Su0ns7Fi8jMxMiRqFCh4AM8PbFlCyjVclc0nj2DtrbcpRgtgnnzsGBBnj2FZmaYNAkTJ9LRY2uLRo1YDArT0uDvj3nz2OpfIYmJeZWe3goAIa2eP39FWw4AxMSgRQv06QN//29ft4kT8b//FVUlsWtXnDzJmUAWKTtPF/ZZutS/SpWmJibNrK1Xjh2LBQsKPXL7dmRnw8Oj0AOsrGBsjPP8lmImuHULZau2laMj1NVx/Hiexhkz8OIFtcRDixbhzz/ZCgrXroWDAywsWOlcUenb10lJyQfYVamSu4eHK205uHULHTti3jwsXPit5dQpXLmCJUuKOmv6dLx9i9RU9vWxDO8Iv/H169c1a/Z++nTny5c7OTknBg36r7AjU1Nx+DAaNEDR6zb4JTNMERmJpk1pi5CROXPg4wNCfraoqmLDBowfj6QkCnpsbdGwISvD9ampWLsWc+Yw37MCI5Fg+fImw4dvW77888GDnjNnjqer58QJ9OqF4OCfP+6TkzFmDDZuhJZWUSdaWkJTExs2cKCRXXhH+I309HShUB8QAgINDYPUwn/khIZCUxPjxhXT4eDBuHABHz4wrLMcUuYiQgC9e0NZGWFheRrbtEHXrj9/cXPM4sVYtgzZ2Qx3u3YtOnfmw0HZmDcP2dkICrKcOXNyp06d6IrZuRMjR+LoUfTo8bNx2jT07ClVLa1mzXKLk5dteEf4jSpVqlSsWElVdbi+/qiGDYW1atUq7MhNm5CZ+W0/aRFoacHJCbt2MayzvJGTgwcPYFPW0nUJBJg9GwsX5gkKAfj5Yf9+OvXrbW1Rvz7DQWFqKtat48NB2QgLw+7d2LeP8kJoQggAf38sXIjLl/PsjggPx/nzWL5cqn48PfHwITsSOYR3hN+IiUFS0u7Q0NFhYSPOnw8VFJLd+eVLPHyI4cOhplZ8n7mjo/kehTwycf8+atUqZnxGPunRI/vFCw99fdumTXu8e5ebYgZ6eli2DKNH01lItXAhli5lMihcswZduqBePcY6VHiePoWHB0JCYGhITcOdO3fNzJoZGzc1MxuweXPOlSt53sG0NHh6IigI2tpS9ebqCkLK/JIZ3hECgFiMYcOwbJmgTx/7Fi1aFOYFAQQHQyD4mWW7aFq0gLo6rlxhTGc5pCyOi+ayYUNwRkatxMSoqKjpnp6zfrQPGwYdHQQFUZBkbw8LC+zYwUxvyclYu5bfOygDGRlwdcWSJWjViqYMd/cZr1/vT0i4/fZt/fHj9+UrJT5zJtq3R7duMnRobk7n88wgvCMEAD8/aGtjxIhiDiMEW7agbl0ZquK5u/NLZkpFWVwpk8vTp2+zs+0AENL01as3P9oFAgQFYdEifI8SOWXhQixZwkxQ6O+Pnj35cFAGxoxBgwbw9KQsIzU1DagCQCIxS0r6/Oufrl/H8eNYuVK2Dp2dcfUqgwIpwDtCxMbCzw+bNxewcT4f4eHIyMDkyTJ0PnQowsKQmFgageWashsRjhjR38BgrkCwRSAY7Ow85Nc/mZtj9GhMmUJBVYsWqFcPO3eWtp/kZKxbh1mzij+SJ5c1axAdjY0baesA+vQZLhT21tT0qVx5zZAhzj/a09Ph7o5166ArY2b7iRORlFS2y8+Vd0cokWDkSCxahMIXx/wkKAgiEfr3l6F/fX107449e0ossFyTkoJXr+S0Kn2x2Ng0uX49JCBANHas94ULHr/urwcwdy7u38eJExSE+fgwEBSuWYNevYrZQcTzgxs34OuLw4ehoUFZSXo6Tp709PNbvWdP40ePLlb9ZWB07lzY26NvX5n71NODsTG1JILMwH76bwpIX31i1SrSvj2RSIo/MiWFqKsTDw+ZxVy4QKysZD5LPuG4+kR4OGnd+ttrIyOjD3KYpV8KxGLSti3x98/ffvYsqVGDpKZSkOTgQLZsKeqAoisqfP1KDA1JXBzDqsooxVaf+PCBmJqSkye5kVMM48YRN7cC2m/cICYm5ONHqTr5/ZKHDiU1apRWGxvw1SeK5+VLLF8u1aAogH37IBBgvOw7Xzt0QEYGIiNLILC8U3bHRX9FKMT27Vi8GI8e5Wnv3BktWxaTuYMlfHyweHHJg8K//oKjI+rWZVSTgiISwcUFo0ahe3faUoALF3DsGNasyd+elQUPDwQElHwt69Sp+PdfpKeXUiA1yq8jlEgwfDhmzy6gjlKBrFmDatXQuLHMhgQCDB/OL5kpCWV3pUw+ataEjw+GDctfsfmvvxAcjOhorvW0bIk6dUo4Yp+UhKAgfu+gtEyfDh0dubhdSUnw8MDmzQVMAS5YgIYN0a9fyTtv3Bjq6mX4KVd+HeH69cjOljbCe/YMz55h2rQS2ho+HAcOICWlhKeXWxQjIsxl7FgYGeVfj2dsDB8fjB6NfDOIHLBoEZYuhUgk84mrV6N3b6nm1Hn27cOxY9i+XS5Sxk+YAEfHAvZF3L2LHTuwdm1p+2/aFLt3l7YTWsjB+0OD16/h44MtW6CkJNXxucWYBg4soTkTE7Rrh9DQEp5ePvnwAenpqFmTtg6GEAiweTP8/REVlad99GgIhazXC/ydli1RvbrMT67ccJBfLCoNDx9iwgQcOgR9fdpSgGPHcP06/vwzf3t2NoYNw19/MVAc2MODwtgGU5RHR0gIPD3h5QVLS6mOl0iwfTt69ZI21UKB8Dm4ZSU3HJRm+rasULUqVq6Em1ueevFCITZuxJw5+PiRaz2LF2PJEtmCwlWr0LcvHw4WT0oKXFzg5ycX2QETEjB2LLZtKyBD05IlMDPDgAEMWBk8GGIxzp1joCvuKY+OcPNmfP0qwy6uCxeQkQFv71IZ7d4d8fFl+BcT90RGKs646A/c3GBpiUWL8jQ2aoShQzF9OtdiWrWCqakMM4Vfv2LDBnkMB0Uikbv7H3XqtF6wYFHxR7MPIXB3h4MDhg2jLQUAMHYs3N0LKCsYHY0NGxhLCiMUom5dBAYy0xvHlDtHGB+P+fMRHCztoCgAX18YGMDWtlR2hUIMG4Zt20rVSbni1i0FWSmTj6AgbN+eP/HeokW4ehUXLnAtZvFiLF4sbVC4ahWcnORxsLpbt4E7d2Y8f75s8eKLc+b40JaDZcsQHw8/P9o6AAA7duDJE8yfn6fx5MmTvr5rhgx55ueHfCnWSkO/fvjnH8Z64xQOdnJwTxH7CLt3J0uXytBVUhJRUSErVzKg6t9/iaEhychgoCtacLaPUCIhurrkv/9+tpTdfYS/c+gQMTcnaWl5Go8fJ+bmJDOTazHt2pEdO/I3/r5R7PNnoq9PXrzgSJVMqKo2ABIAAtzQ0OizfTun37J89+rsWVKlCnn7ljsBRfDvv8TIiNy7l6dxxoxFOjruwA5VVZvY2NgSdFvY1smEBCIQkNevS9AlW/D7CAtg+3bEx2PGDBlOyc2yPXIkA9ZNTWFjg8OHGehK4Xn6FJUqwciItg52+N//0KxZ/sF2R0c0aABfX67FLFggVVC4ahWcneUrHPzyBaNHo2JFZGc7APuANCBEWdluzBhoaEBXF23awMcHr19zJ+n1a7i5ISSEyTCrxBCCESMweTKsrfO0h4aeSE4OBtxycmYcOBBWyNklwcAABgb46y8Gu+SIcuQI37+Hlxe2boWKigxnrV2Ltm1RqRIzGvglM1KiSBsnCiQgAMeO4cyZ/I2BgXjyhFMlHTqgalXs3VvUMZ8/Y+NGeHlxpak4Tp5EixYwNMTBg/DwwKtXs+vU2auq2sDe/sXnz94ZGXj/HkuXQlsbQUGoWROqqqhdG0OH4tgxFneqZGXB2RkzZ6JdO7ZMyERAANLSMHNm/nZNTWMgCiBaWlctLBj+adOpE44fZ7ZLTuAgOOWeAodGnZzIggWy9RMXR5SVyZUrTOkiWVnE2LgM56bibGh04kTi55enRZGGRnM5d46YmpIvX/I0rl5N2rWTKucfg1y4QOrUIb++sfnGvmbNImPHciqpQBITyZQpRE+PKCkRe3tpv5hZWeToUTJkCKlVi1SoQIRCYmxM2rUjy5eTxERy4sQJFRVTobCmrm69xMTEEqj6ca9GjiROTly/d4Xx7BkxMiK/D3yeOkV0dV+Ym3evWtV2zBgvSYnkFpFV7uZNIhCQrKwS9MoKUg6NlhdHuGcPadRI5rdn8GBiYMCkMELI9OnE25vhPjmDM0dob0/++SdPi+I5QkLImDHE3T1Pi0hEbGzIrl1cK2nbluze/fOfvz7pPn0i+vrk1SuuJf1KaChp2JAIBMTYmMyfX6rn7M2bZMIEYm1NNDQIQAQCB+AhQIA/+/QZWoIOc+/V5s3E0pKkpJRcGIOIRMTenqxbl7/9xAlibEyuXy9t/0WnV1VXJ0FBpTXBFLwj/OkIExJI5crk1i3ZOhGLiYYG804rNpaYmBAp3hp5hBtHmJ1NtLTyZ6NWSEeYmkpq1ybHjuVpvHWLVK5MPn/mVMn586Ru3Z9B4a9POi8vMm4cp2J+8OEDGTaMaGoSZWXSoQO5fZv5/pWU6gEpAAEOAT6qqqRyZWJvTzw8SHCwVO4/OTk5MpIYGZEnTxiWV2KWLSMODvlj09BQYmREIiIY6L9oR9i6NbG3Z8AKI/CO8Kcj7N+fzJolcycHDhAlJVaeR23akCNHmO+WA7hxhJGRxNo6f6NCOkJCyLVrxMQkz/pYQsi4cWTUKK6VtGlD9uz59vrHk+7TJ2JgQN68Yd16Tk5O//7uhoY2Li7DxWLx8ePE3p4IhURPj0yZwuIq0CFDxggETQWChUKh6alTFy9dIsuXEycnUr8+0dMjQiERCIiGBqlVi3TvTqZMIUeP/oxHz549q65eWyisoao6Zf9+MVsSZeTePWJomH/pZkgIqVKFREczY6JoRxgcTCpUYMZQ6eEd4TdHePQoqVevJF8kGxvSrBnDwnLZuZP06MFKz2zDjSMMDCzADSiqIySEzJhB+vfP05KURExM3g0evGj+fN/PXMWGZ8+SunWJSETIL0+6mTPJ+PFcWO/de4hAMBx4AAxXVl6trEwcHMidO1yYPnr06IwZMx49elTgX+/cIX5+ZMAAYm1N9PWJsjIBiLo6MTUlQmFT4B4gATzmzJnLhdbiyMwkVlb598Ns3kyqVCEPHzJmpWhHmJNDhEJy6RJj5koD7wj3EkI+fSKVK5dktUtSElFSImFhzGsjhKSnE319+dptIyXcOEJ3d7JpU/5GBXaEGRmkQYOf0RghJCMjw8DAFtinrLzVwqIVZ0ratCEhIYR8f9IlJHAUDhJCdHWtgDiAALHq6q04rHopMwkJJCSETJhABIL6gAQgwLbGjV3kYYLQ25s4OeVpCQoiNWqQp0+ZtFJsCca6dfP/tqMFv48QACZOxKBBaN1a5hOXLoWGBnr2ZEEToK6OAQOwfTsrnSsACr93Ih9qati1C5Mn482bby0xMTECQRPAVSQanpioHx8fz42SuXOxePHPDQYrVmDgQFSrxoXpmjUbA5uAt0BAu3ZVlZW5MFoyDAwwcCDWrkXjxtUEgiHANoFgYaVK42rUgJsbzp+nJuzGDWzfnidlmp8fli9HeLi0xeaYom9fXLzIqcXSwoFP5p7ciDAsjNSqVcIK4IaGZORIpmX9QnQ0MTX9NgxVhuAgIkxOJtra5HcjChwR5uLjQzp3/rbA4cuXL0ZGdsAH4KWmZmNuVurm0ro12buXJCcn54aD3GRIuXOHKCuLTEw8tbQsOnb8X5b8rL4vjoULF3bv/r8rV64QQt6/J76+pE4dYmlJfH3Jp0+cKklLI+bm5ODBny2+vqRePVYC+mIjwg8fiEBA4uOZNy0r5X1oNDj4oKkpuXy5JKffvk2EQtbfRTs7cvo0uyYYhwNHeOECadOmgHaFd4Q5OaRZM7Jx47d/nj59rmHDjtbW3SwsImbM4E7G6dOkfn3y9WvytGlk0iQuLP77L1FXJz17cmGLDX73CrdvE09PoqtLnJ3JuXMc7SwcO5YMG/bznwsWkPr1ybt3rNgq1hESQvT1yfTprFiXiXLtCB0cBqip2Tg5lcgNEtKlC7GwYFZRAWzcKC/D6NLDgSP8808ybVoB7QrvCAkhMTHEwCD/dM7Hj6RRI5lzQZSY7OxsPT03bW0bZeWO16+znlr082eio0OaN2fbDosU5hUSE8nGjcTKipibE1/f/AuDmeXX5AwSCZk6lTRpQhIS2DInjSPs35/UrcuWAOkp13OEX792yMw8devW1BKcK5Hg4kUuauIMGoQLF/DhA+uGyhaRkYpZdEIaLCwwaxbc3SEW/2w0NMSFCzh0KH/xJpbYtGlbamrdlJQokWjZ7Nm/pedilPR0WFhAXx/Xr7Nqhw6VKsHTE9HRCAnBixeoXx8uLggLy/PmMkJSEjw8sHkzdHVBCCZPxuXLOHcOBgYMG5KJKVPw/Lls1S4popiOMCNDHzDKyRFIZE8sGBAAoRAjRrChKw9aWnBywq5drBsqW5S3lTL5mDwZysrw98/TmOsLQ0OxZAnrAl69+pCdbQUAaBgf/x97hiQSWFtDKMTjxxAq5nPoG7a22LgRL1/CwQELF8LMDN7eePlS7OY2sVo1u9at+75//740/Y8fjz590LUrJBKMHIk7dxAeDn19puSXkJYtoapadpYEchCcck+nTn9oaIx0cRkt01lpaWne3t6VKm3u1Sut+KOZ4Pp1Ym4uL8kJpYHtodH4eGJoWPCfysPQaC7PnxNDQ/L4cf72//4jlpZkyRJ2rT969MjIyEZJaZ2eXs+1azezZ8jGhujocL2ihA2kGSf8ldu3yejRREMjWElpFkAEgrPdu5cktVsuP0p6iURk2DDSrh2RUU5JkPKS7e1JK+72/hSMlEOjdBYpZ2Zm5q6zatu2rZqaWoHHfPjw4d27d9bW1sq/rKROTk6+cuWKtrZ2q1atlAovrauj837ixHrLli2TSZWRUcO0NCfg1fnzjUSiWGX2V3AbGj57/XqEnl5a06Z1jx/fXtitKD+U83Awl1q1sGQJ3Nxw6tTXChWUtLW1c9uNjBAejo4dIRBg9my2rFtaWkZGHvn777/t7BY0ZW2QunNnPH6MJ0/oBy7cY2sLW1soKb1av741AEJanz69yMICNWvCzAxmZj9fFFGGTCKRXL16NTVVZdw4+yNHBCoqGDAAqak4dQrq6txdS9EMG4apJZmeogH7Ljk/nz59ql+/frt27Tp06GBhYZHw25Tu69evq1atmusV/vtlivnJkycmJia9e/du0qRJhw4dilhjXURh3sK4efOmUNgeIAARCnueOXNGptNLRtu2/YDbAKlQYeWqVQEcWCwlbEeEc+YUuiqk/ESEuZiaztDUbGVo2GzevDxVoT98IPXrkz//ZNe6rFGOTAweTFRUyN277FnglJLdq6ioKH39ZsCuihVdlixZ+/w5OXeObNxIvLyIszNp1YrUqkXU1EitWsTBgXh6El9fEhpKbt8m794RiUTSunVvPb2xqqojatcenJVF+vYlPXtyV45YykvOTTHDYPWeEiC/EWFgYGCtWrVOnDgBoG/fvgEBAQsXLvz1AAMDg9OnTxsZGRkbG//avmjRokGDBq1atSo7O7tp06aHDh0aOHAgU6qqV69OSDyQCQiB5zU5qUCakPAJMAeQlVXv7dvbHFiUc27dwuTJtEXIAS9fvkxPf5iWdijXAkEAACAASURBVDUtjQQF2c+YMfpHXGhsjLNn0aEDBAI5KhAoPd7e2LcPp0+jcWPaUqhiY2Nz+fLWEyfO2tqOcnBwAFCrVv5jUlLw6hVevvz2/5s3v72WSJ5mZKiJROsBJCd379r1PyMj49BQ2SqtcoCyMszMEBhYkpQmXMOBT85HkyZN9nxPJ7Vv3z7r3/MrE0II+fTpE/JGhBoaGre/55/38fFxdnYuzEQJIsL4eKKktEUgqCMUmrq6yja5WGI2btyhr99FKPTT1m4SExPDjdHSwGpEKJEQXd1CV5mXq4gwJibG0NA5d3xCQ6NdfHz+dKNv3pA6dYi/P1sCWIoIV60iQiHZt4+NvqnBavRcII8evdPSag1IAJGycpP+/ZM4zkgn/SVPmVLorD83yG9E+ObNG1NT09zX1atXf/v2rTRnJSYmpqen/3rimXzlvX8hKSnpxo0bAoEg958VKlRwdHQsonOJBE2bKjVs6BEV5Z7bImZ8jXNBeHgMbtbMav/+x0ePhtWta8KN0dIgFouFQuGPG8sssbHQ1VXS1y/0NojFYvm/RYxQp04ddfWXQmEPIENNrWabNnobN0ratyc/DqhcGefOwcFBCZCMG0eK6KpksHGrDx4UzJgh/OsvSf/+RJHeRu4/lpUrqxPyEWgMSASC7B07VAUCTiVIf8mTJ2PNGqX4eHHe0T3uyK08XOxhFBxhdnb2j3UoKioqWVlZ0pyVe9ivJ2ZmZhZ2cGJi4o0bN364WG1t7U6dOhWxuKZfvwpfv+LOnXTptDCJubn53LnmoaHqERFZjRvLvNmDY7KysiQSCUvfuRs3lG1tUdjngRCSnZ0t5aelrPPkyb/JyboSiS+gIhS6+fh8HT68YufOkj//zNHS+vatNjLC338LundXE4tFo0czvF2L8Vt98aJw8GC1SZNyRo7MUbD3kPuP5bNnzzQ0WqalbQKEurrOb968/hEhcIP0l2xggEqV1NeskSxcmMO2qgKRSCRFPPl/QMERVq5cOXfYE0BCQkLlypWlOcvQ0FBJSenTp096enq5J1apUqWwg83MzJycnAYMGCBNz5s34+xZXL4MIyMNaY5ng+HDsWePWsuWtOxLi0AgUFFRYWk9bXQ07O2hoVHwuyAQCNTV1Qv7q9ySnp5etOacHLx4gbi4b/89fYq4OHz+TEQifcAKQFaWdrdu4l69BDNnKtvYKAcGok+fb+fWq4dLl9Chg2qFCqp//MGkbLFYzOCtvnsXfftiyBCsXq0CyNlEVqlh9l5Jg5WVFfAQOKqklK2j897c3JylQZrCkOmS27XD8eMqK1bQed+l/OFOYSNrq1atwsPDc1+Hh4e3atUq93XRAaySklKLFi0KPLE0PHyIsWOxeDHl6dzhw7FvH9LTaWqgjoLtnUhOTm7SxKFmzW41atg9evQotzExEVevYtMmeHvDxQV2dtDRQYcOWLsWcXGoVQtTpuCff5CaatG6dba+vruurquqqrmHh75IhI0bsXcvZs6Eiwu+/5JE9eoID4efX56aA3LFmzdo1QrdupWdvdVyj1hcQV39mLPzo5kzX12/fpxjLygrkyfj6VPIntqEW9ieq/yd6OhobW3t5cuXr1y5UkdH5969e7ntmpqaZ8+ezX09a9asiRMnApgwYYK3t3du4/Hjx/X19YOCgry8vAwNDT9+/FiYCSkXy2RlEX190qFDqS+JCXr2JLt20RZRHOwtlsnKIpqaRZUKKXOLZRYtWqmsvBEgwAMDg/9ZWxN1dVKlCmnfnnh6kpUrybFj5MkTUtgmIIlEcuvWrXv37mVlkenTiakpCQ8nhJD0dOLlRSpXzlN89dUrYmZGNmxgTDxTC0ASEsp8KtFi4X6xzNSpefJrc4+sl6yqSnbuZElLMcjvYhkrK6t//vln586dhJBLly5ZW1vnti9ZssTc3Dz3ta6ubsWKFX19fX890dHRMTQ09MiRI1paWhEREYaGhqVU0q4dCMHp06Xshhk8PLB2LYYMoa2DEtHRqFsXmpq0dTBHYmKKSFQXAGCkqZmydSvq1sX3TRDFIxAIfuxnX7kSnTtj6FC4ucHHB76+cHWFhwdCQ7FhA6pVQ40aOHv22157T09WLqcEpKfD0hIGBoqZSpQW9+5hzx48eEBbhyxYWyM4GEOH0tZRBOy7ZApIExHOn0+UlUlsLDeKiicnh1SpQuLiaOsoEvYiwoAA4ulZ1AFlLiLcsOG5klJjHR1vA4MWR46Elb7Djx+JoyOxs/v2IcnOJr6+xMiIbNz4LUtfXBwxNc0TKZaY0kc5YjGpXZsYG3O3y5sWXEaEYjFp0YIEB3NmsGBkveS1a4mGBktaiqFcV58oluvXsWQJ1q3D9xCUPsrKGDIEW7fS1kEJBSs6cfEi5s+vdfHixYMHO965s79v316l79PQEMeOYdQotGmDTZugogIvL1y4gK1b0a4d4uJQty7On8ecOdi9u/TWSkvTpvj4EQ8fotznDWSS9euhrIzhw2nrkJHRo5GZichI2jqKgAOfzD1FR4SJiURDg7i4cKlIKmJjiYkJkeLnCzXYiwjr1yfR0UUdUIYiwuhoYmJSwqLQ0vDoEbGyIv37f6s/JxaTjRuJoSHx9SUiEXnyhBgb3zE0bGZo2MTRcZhIJCqBiZJFOYmJicbG1kpK1ZWVW6uqxr9gvZqhXMBZRPj+PTEyIo8ecWOtKEpwyWZmZMgQNrQUAx8RFkqLFjAwwP79tHX8hrk56tbFyZO0dXBOSgrevoWlJW0dTPD8Obp3R0AA2rZly4SlJSIiUKUKmjTB1asQCuHpiRs3cOYMWrWCSARt7WkJCfsTEu6Eh1cJDT3Alo7f+N//hv/332Cx+LVINE1fvz8naQrLERMnYtSosvo16dUL587RFlE45c4Rjh+PFy/kd/bewwPBwbRFcM6tW2jSBOxX+2Cdjx/Rowd8fNCvH7uG1NTg7w9/f/Tvj4ULIRajdm1cuAAPD3TogP/+SwKqAUhLqx0bm8CulF+Ii0sEcjfD2iclfSrmaB5ZOHMGd+9izhzaOkrK9On4+BGfP9PWUQjlyxEePoygIISGompV2lIKwcUFN24gPp62Dm5RjB2Eycno1g3u7hg5kiOLffogOhoREWjdGi9fQiDAqFG4exc5ORKgN7BMIPDz8+u9fDnzVdHzIZFgxAi8f78AmARsFwhcBwzoyq7J8kRGBv74A/7+clRiSVZq1EDFili3jraOQihHjvD9ewwciFGjfibmkEPU1dG/f7nbeqwAK2UyM+HoiBYtMGsWp3aNjXHyJAYMQPPm2LsXAKpWhZGRCrAIqAdMdHW9dOEC7OwQEcGWhlOnoK+PAwewd2+Ho0fnOTmd2rRpaHDwWrbslT98fNC8OXr0oK2jdLRujd27X3/8+JG2kILgYLqSe35fLCMWk2rViJUVLUUyEBlJataU07L1LC2WqVqVvHxZzDHyvFhGJCL9+pEBA4hYTE1DZCSpW5cMHUpSU8n06T4VK7oJhZt0dJoYG79wdib+/sTEhHh6kqQkqXqTcjVEYiJp144IhcTJqdDMAAoP24tlHj4kBgbk3TtWjchGCS45JydHW7sJ0EkorDlkyB9sqCoQfrFMHpyckJiIK1do65ACOztUrIhLl2jr4Ir375GdDTMz2jpKCiEYOxYpKdixA0J63yc7O0RFQSCAnR0GD57v7V2/f//LZ86sf/asprU1liyBgwNycmBpiZ07mbH4118wNsbTp4iMxOHDUFVlplueX5FIMHo0Fi9G4ZmVywb+/v6pqY2B8xJJTEhIGG05+SkXjnDLFpw4gZMnoaNDW4p0jBhRjpbM3LyJ5s1piygF8+YhOhqHDtH3BNra2LED8+ahVavFPj4xBw+27dNnbELCyzlzEBcHMzMcP46WLbF6NTp2RGxsyQ3FxcHcHDNnYuZMvHsHGxvmroEnL8HBkEjkKFtQicnKygJyE0cpAwKJnOUeVXxH+PgxxozBwoUsLmdnnKFDcfIkEhNp6+CEMj1BuHYtDh7E339DS4u2lO8MGoRKlcIyM7dJJJ4JCVO2bAkDUKkSFi9GbCxq1cLr18jKQuvWWLgQspYPyg1Q6teHhgbevsXixaxcAk8uHz9i7lxs2EBzpIEpJk6cqK5+QiBwATp26NBcKGeXJF9qGEckQtu2aNUK8+bRliILlSqhe3fs2UNbBydERpbVJaN798LPD2fOwMCAtpS8mJjoAY8BoqwcGRhoamWFWbNw9SoqVYKvL+Li0KYNxGLs24cGDXD+vLTdXrgAQ0Ps3o3t23HvHmiVWi0/TJuGYcNgZUVbBxNoaWklJT3dvbufULhj7txQ2nLyo+COsF07SCRyvZGzMDw8sHkzbRHsQwiiomBnR1uH7ISHY/JkhIWhRg3aUn4jJMS/UaPJlSvbuburff7cd9s2VKoEb28YG8PFBUeOYMoUREeja1d8/Ih+/eDkhKKX8qWno0cPdO6M5s2RmCjf2ZMVhUuXcOUK5s+nrYM5lJWVBw1yrVfPbPVq2lJ+Q5EdoY8Pbt3C9ev0J29KQIcOSE9HVBRtHSwTGws9PbmLqIolMhIDBuDQIXwvnSJf1KtX7/798/HxUZs2rVRSEtjawssLV6/i8WP06oXz52FhAScn6OoiJAT9++PcOdSqhRUrCi4aFxgIPT3cuYMbN3DyZJn8NpU5srIwZgzWrZOjIXemGDQIly/TFvEbiukIExOTlyzZvWiRZN06WFjQVlMiBAK4uyv+kpmyuJX+2TP07YuNGykXcy4BRkZwc0NoKN6/h68vMjMxfTrCw9GzJywsMHcuatVCVBT27NnTvr3DtGnTXr6EpSUmTcKECfjwQao1TeHhFydPnn/o0FH2r0aRWboUjRrB0ZG2DhaYOhUpKbh3j7aOfHCwk4N7LCz8gIlKSvJRcrekvH9P9PVJWhptHb/A+D7C8ePJ6tVSHSkn+wjfvSM1a9Kvg8Mgz5+TNWuIgwPR0iIVKxKB4C7QENgAdBcI/K2tSXy8tF2FhZ3S1e0CnNTRGbJq1Xo2VcsRjO8jjI0lhobk7Vtme2WSUl5yjRrcJeCWch+hgBBC2xczj7p6cmamDlA7KytGtSwP5fTujf794eZGW8d3MjIyVFRUlJnLCmpvDz8/qUIrY2Pj+/fvG1NdoZGYiLZtMXQoZs6kqIItvnzBmTNwd2+ZnT0LcAQygIYODs8JwdevxZybnY20NHz8OCE1dTBgDyTY2AyLiioX+eNTUlK0pS+4XByEoGNHODlh4kSmumSeUl7y1KnYtQsJnCTBlUgkYrFYRUWl6MPKfp7jgqhS5faLF/WBjPv3VcviQowfeHhg1So5coTMkp2Nhw/LzC60jAw4OqJrV8X0ggD09DBwIFavVr59+ybgCEQpKYm9vCAQoFKlYs5VUYGWFoKD669eHZaZ2VQgOJqRYZmZyRcjlJmdO5GainHjaOtgEy8vrFmDN29gakpbyncU0xFqaAQDl3v1cnN0ROfO+Osv6OvT1lQievbEH3/gyZOyOtNZNPfuwdwcGhq0dUiBWIwhQ2BmhhUraEthmYsXTxoZWWdk7BQKcehQgIODDOf6+HgmJs49fbplkyY2ysorbG2xd6+CrP7nhi9fMGsWwsKgpERbCpsYG8PICCtWyFEObsVcLGNpKdq71y8szDcmBrq6sLbGzp0oi2PAysoYOlRhc3DL/0qZT58+DR06sVkzx06d9icnY+tWRdjaXDRaWlrp6c+Tkx+LxW/6yJifXllZef163xcvIg4dCtq/X9vLC506wd+/TH71qDB9OlxdYWtLWwf79OyJo/K0oErBv9aVKsHfH8ePIyAAHTogJoa2INkZNQo7diAnh7YOFpD/nDLOzmNDQtpERgZdu7Zz9uybZXm6mQJubrh1C/v2wclJfgvRyQ9XruDMGSxcSFsHJ3h749274ueeOUPBHWEuNja4eRODBqFNG3h7IzOTtiBZqF0b9erh779p62AB+Y8Inz59JZE4A9XEYufo6Nu05ZQ9atbE5cuoUwdNmpSjPPIlIDsbo0cjIAAVK9KWwgl166JiRaxaRVvHdwp2hJ8+fTp06NCCBQvGjRs3adKkpUuXnj9/Pj09nWNxDCIUwtMTDx4gPh4NG+L0adqCZEEhy9YnJeHtW9SvT1tHkVhYNBQI1gK39PS2dexYdpLVyhOqqvDzw7ZtGDIEkyYp5thG6VmxAmZmcHKirYNDOnbE/v20RXwn/2KZc+fOBQQE/P3332KxGICmpqZIJMrKysp9PXDgwEmTJjVs2JCCUiaoXBk7dyI8HH/8gbp1sX69HC1bKoL+/TFlinwtsio9kZGwtQVzGzGYRyJBVlaAg8O6ihW3jR3rY8Wv+igFnTrh7l0MH45WrRASgjp1aAuSJ549w5o1iIykrYNbvL1hbw85WVr8MyJ89epVly5dunfvnpmZGRgYGB0dnZmZmZqampmZmZycfPXq1fnz50dFRTVu3HjUqFFpaWkURZeSjh1x9y5sbdG4MZYvh1hMW1BxqKvD1ZWxMnJygvyPi/r6QijUPH3a+8CBoI4d21NWU/YxNERYGAYPRsuW2L2bthp5YtIkzJqFmjVp6+CWpk2hpoYNG2jrAPCrI4yIiKhdu/bLly/PnDkzevRoKyurChUq5P5JW1u7VatWM2fOvHPnTkRExIcPH969e0dJMDOoq2PhQkREIDwcdna4eZO2oOLIHR2VsxpepULOV8rcvQt/f+zapfjLRLlEIMCkSbh4EStWwM0Nqam0BdHm8+fPe/ZI3ryR6+3z7NGypbwsif/5LXd1dQ0KCjItbvTN1tY2LCzM3NycZWFcUKcOzpzB7NlwcoKbG+7efeXhMc3dfcrz589pS8uPjQ0qVcLFi7R1MIc8R4SZmXBzg78/qlenLUURadAAERHQ1YWVFW7coK2GEklJSY0ata9ff9CwYXZTptwpLvOJYjJxIh4+lIvf9/zPXTg7IyYGFSuKmjbtv21b5507u7dv75Ila8VS9lGkJTNv30IkkscCRrnMmAFrawwYQFuH4qKuDn9/LF+Ovn3h51ceNxoGBGx58mRwQsIZsTg0KEiBii3JgqMjlJQQEkJbhzSO8Nq1a0OGDOn8HQ40cU+lSpg69U3FirUI6UZIl+zshs+ePaMtKj9DhuDUKXz6RFsHE9y6JVUpAyqcPo2wMAQE0NZRDnB2xq1bOHoU3brhwweIxeL4+Hix/E/aM0FaWqZYnLtVomJGRpna0cUoNjZyMU1YzKI9kUg0atQoX1/fqlWrciOIFlWrVq1Q4SlwD1DKynpQU/5mritWRK9eCAlRhOkEuZ0g/PQJI0di167is2vyMEKNGrh8GX5+aNz4hUTiIhRWqVAh/tKlA3L4BWSWkSOH+fn1VlO7pqYWsXRpOY0IAXh6YuxY2iKKjQgzMjLs7Ox69+5t+x1uZHGPqqrq6dM7HBxWNG++VCgMfv1aHjNgenhgyxbaIphAbicIR4/G4MHo0IG2jvKEkhK8vNCggV9Cgt9//x3/99+Vs2f70RbFOlu3VuvZ88rJk8737x/t3bsHbTnUGDYMIhFOnaIso5iIUFtbW1dX9/Xr1zXkdj6HOaysrM6dCwGwZQtcXBARIXf5oNu1Q2am/IZTUkII7t6FHFYF2boVL15g717aOsolWlrZQO73TePz52zKalgmOhqbNyM6WtPEpKwVd2YaoRCWlvD3R/fuVGUU8bfatWvr6elt27bNzMxM7zucKaPIyJGwsZHHEUiBAMOHl/klM0+eQF9f7uqBvHwJb2/s3Ak+oSgVFi2aYGIy2sBgXMWKYyMjJ+7bR1sQa4jFGDkSvr4wMaEtRT5wc8PVq5Q1FOUInz9//uXLl9xixF++w5kyuqxfjxs3sGMHbR2/MXw4QkPL9gYsORwXlUjg7o65c9GoEW0p5RVra+u4uH9OnHB78+byhQuN5s2DmxvKclbHQvnrL2hrw92dtg65Yfx4ZGQgIoKmBhm2T4jF4oyMDPakyBWamggNxYwZePyYtpS8mJigTRscPEhbRymQw6HdpUtRoQImTKCto3yjra3dvHlzbW1tGxvcuYOcHDRrhkePaMtilFevsGIFNm+GQEBbitygpoaaNbF6NU0NeRyhvb19YGBg7mtCyODBg2/+knNl//79GvI2acYmDRpg+XK4uMjdz9KyvqFQ3iLCqCgEBCA4mH82yRHa2ti7FzNnomNHudhnxgiEwNMTXl6oXZu2FDnD2Rnnz9MUkMcRvn//Pjk5Ofc1ISQkJOTVq1cURMkNw4ejaVOMGkVbR1569MCLF2WytiKArCw8fowmTWjr+E56OgYPxtq1CpXQXGFwc0N4OJYuVZBh0m3bkJCASZNo65A/ZsxAYiLi4qgJ4DPLFMP69Xj0SF4S4uWirIxhw7BtG20dJeLePdSrB3V12jq+M20a7O3h6kpbB08hNGiAmzchFsPODg8f0lZTCj58wKxZ2LpVriuu0EJPD5Ur0yxPyDvCYlBXR2goZs5EdDRtKb+Qu+m7LJZ2k6tx0TNncPYs1q6lrYOnSLS1sWcPvL3Rvn0Z3kc7fjxGjZKjsRB5o3dvhIVRs847wuIxN4e/P1xckJJCW8p3atYk6uprGjTo4+29JDu7LG26kp+VMgkJGDECW7dCR4e2FB4pcHPDlSvw94ebG8pcFbgTJ/DgAebOpa1DjvH2xocP+PiRjnXeEUrFwIFo3RqenrR1fCcoaOv790+fPg1au1Y8b94K2nJkQH4iwrFj4e6Odu1o6+CRmvr1cesWKlSAnR0ePKCtRmqSk/HHH9iyRS4q0MotNWpATw9+lHIK5XeEc+bMUVJSUlJSUlFRATBo0CCl7wwZMoSGQnkhMBBPnmDrVto6AADh4ZGZmSOAKhkZnpcu3aItR1qSkhAfj/r1aesANm3Cy5dYsIC2Dh4ZUVfH5s2YNQsdO8Lfn7Ya6Zg+HY6OaNOGtg65p3NnahvD8szburi4fFKM6gYsoKaG0FC0bg0bGzRuTFlMr15tz51bk5w8U0kpuGfPMhPU3LoFW1soKVGW8eIF5s5FeDifRKas4uaGZs3g4oKoKKxfDy0t2oIK5/JlnDxZtpf5cMasWWjcGKmpFN7QPI5w5cqVXNsvU9Sti7Vr4eKC27cpTyy5uw8Si8WHDq2OjGzesqWcbe8oHHkYFxWJMHgwFixAw4aUlfCUBgsLRETA2xt2dggNReXKCXFxcVZWVtra2rSl/SQrC2PGIDCQL2YiFVZW0NREYCC8vLg2zc8RyoarK9q1k4vJQg+PoSdPblu3bszMmUryUOJZGuRhpcySJdDRwR9/UJbBU3pyq/vOmYO2bS/XqtWjd+/QevXaylUl0fnz0bgx+vShraPs0LYtdu2iYPenI3z27FlCQoI058TExCQmJrImSd5Ztw5xcdi8mbYOAICrK9TUsH8/bR3Fcft2VJ06rcLCbA4cmEro1SOPiMDGjdi+nU8iozgMHYqaNdekpu7/8sX//XvfFSs20Vb0jTt3sGNHmZnIlBOmTEFMDEQiru3+dITR0dG1atWaMWPGo0Ky+xFCLl++PGjQoMaNGyclJXGlUO7InSycOxd379KWAggE8PXF7NnIyqItpUiGDJny/HmIRHLnzJn006dPU9GQmoqhQ7F+PSpXpmKfhy309dWA3CdS0u3bav/9R1kPAJEII0dixQoYGdGWUqZwcICqKoUEJj/nCPv166epqTlz5kw/Pz8LCwt7e3tzc3M9PT2RSPTly5fo6OgbN27Ex8d369YtKirKzMyMa6XyRJ06WLcOrq70JwsBtG0La2sEBGDaNMpKiiAlJR2oDiAtrWF8/AcqGqZMQdu2cHKiYpyHRdaunde169CcHCNt7eTWrY82aIDBgzF7NoyNqUlatQq6uhg6lJqAskvTptiyBSNHcmuV5EUikYSHhw8bNqx69eo/jhEIBA0bNpw8efLDhw8JE+zZs6d69epaWlq9e/f+/Pnz7wfExMQ0b95cQ0OjYcOG169fz20MDQ2t9Qv3798vrH8XF5e9e/cyIrUIRo8mLi5sG5GKJ0+IoSH59Il1Q+np6Tk5OSU4cfLk+UpKg5WVV1Wt2uTDhw8l6MHIyKhkJ169eq1t2/42NkNNTWOSk0vQQXkkuazdKbFY/PHjx9zXb96QiROJvj7x8iKJiayb/v1excURfX3y/DnrpmnB6scjJISoqDDWm1gszs7OLvaw/I7wV75+/RobG/vq1au0tDTGdBHy6tUrLS2ta9euZWRkDBgwwMPD4/djbG1tfXx8RCLRtm3bKleunHslwcHBnTt3fv6dzMzMwkxw4wgzM4mNDQkKYtuOVIwdS6ZOZd1KiR1heDgxM7u0Y8fOhISEkpkumSP8/PmzkZEtEAPcNDKylUgkJbNe3ihzjvB3Xr0inp7EyIgsWECSklg0lO9eSSSkUyeyZg2LFqnD9sdDWZkcOcJMVww4QpZYtGhR3759c18/ePBAQ0MjPT391wPu3bunqan5o7FWrVrHjh0jhAQHB/84sWi4cYSEkKdPiZERiYriwFQxfPxIDAzIs2fsWimxI+zXj6xfXyrTJXOEN2/e1NcfBxCAGBt3/xE08BSNAjjCXGJiiIsLqVyZrFlDCv/lXCry3auNG0nz5kQkYsWWnMD2x6NJE9KxIzNdSekIKWyfiIuLa/S9ELilpWVWVtbbt29/PeDp06e1a9dW/16hoFGjRnHf63Ncvny5evXqdnZ2AQEBhN7iwx/UqYOAALi64t27lPj4eIpKDA0xaZKcJjOMj8fFixg8mIJpS0vLzMwbwEWB4Ji29ldDQ0MKInjoYWGB/ftx7hyuXUPduvD3Z3dZ2fv3mDcPW7bQzxpRphk+HDducGqRQkWQxMTEH45QKBRqamp++fLl1wO+fPmi9UtqAR0dndwD7O3tz5w5Y2pqeu/evREjRigpKY0dO7ZAE/fv3w8NDR04cGDuPw0MDOLi4pTZT7M5XgAAIABJREFUKX/SrRsCAo7Wru2vo2NQr16FsLAdSpS+BJ6esLXVunAho1kzMUsmMjIyVFRUZL2T/v4VnJ0hEGSVJms5ISQ1NVXW0tA7d6oYGu7q3Hmdnp66l9f2FPnJmy7fpKam0pbAJNWrIzgYUVFKy5errlolnD49e+jQHKaeB7/eq9Gj1T08JDVqlOqjLv+w/fEYMgSTJ2ufOZPesmVpH2USiURFRSU3Y2hRMBN/yoKbm9vcuXNzX4tEIqFQ+DzvtPKhQ4caNWr045+Ojo6rV6/O18maNWvat29fmAnOhkZzqVy5CZABkIoVp546dYozu78THEzatmWx/xIMjWZlERMTUvpVViUYGr18mZiYkKdPS2u6HKIwQ6O/c+0a6diRWFiQHTuYGcD8ca8OHCD16pGMDAb6lHM4+HiYmxPp5sGKQX6HRuvVqxf9vbjfw4cPNTQ0qlSpku+AZ8+epX+vSP3gwYN69erl60RJSYnIwdBoLoQAUAKQk1Mhh2qRQHd3pKTg2DGKEvJz8CAaNkSDBlzbffkSAwZgzx7UqcO1aR55pmVLXLiAdesQEABra6xb96xBgw6VK9s5Og4rzZc3KQlTpvAlJhhj4EBcvMihPQZ8roy8efNGS0vrzJkzSUlJTk5OY8aMyW1ftmxZSEhI7mt7e3tvb++0tLSAgABTU9PcKOTAgQMxMTFJSUkXL140NTVdU/jCLI4jwpUr1xsatqlYcYiKisPr11mc2S2QU6dIvXpEit9AJaEEEWHLlswsAJMpIvz6lVhYkA0bGLBbPlHgiPAHEgk5epRoajoBtwGiprZo/frNJegn916NGEEmTGBaorzCwccjKYkIBAyMJMlvRFitWrXdu3dPmTKlTp066urqy5cvz21PSEj4kbBmz549t2/fNjU13b1797Fjx3InpR48eNCrV6+aNWtOnDhx0qRJEyZM4F58gUyfPvbevX1Xrsz08jozbJgq9/mBfqVbN1SvLi+FvO/dw9u36NWLU6NiMQYMQOfOGD2aU7s8ZQuBAH36oEqVj4AlgMxMq+nT4x0cMGkSNm/GzZsyFOK+cAEXLmDpUhbVljd0dGBqCs7KQAiI3AwwMoirq6uTk9OAAQM4tiuRoFcvNGyIFVRr5UZHo1s3xMYyn/VG1sUyI0eidm3MmsWAaWNj4/v37xtLkSxkyhQ8eoSTJ8HO6qhyQUpKilyVcWAPf/9NixYdT07uoqe3PSwsJDvbIioKjx/j0SNER0NbGw0awNLy2/9tbfF9MTsASCSSkJD9UVGxhw65BAZaOjrSuwxu4ebjMXky9u5FKRPmSSQSsVhc7GIZ/lHBJEIh9uxB06Zo2hTOztRkWFuja1esXInFi6lpAPD1Kw4fRkwMp0a3b8fJk7h5k/eCPFIxaZJny5ZNHj9+3KnT8WrVqgFo3frnX+PjkesXr17Fpk2IjYWR0U+/eOKE96lTOamprdXUhtWrtwcwp3YZisisWVi7Fu/eoWpV1m3xTwuG0dXFoUPo0uXbt4UWS5fC2hqenjA1paYhOBi9enGa7/Gff+DlhX/+ga4ud0Z5yjpNmzZtWkh5sCpVUKUKfoR62dl48gSPHuH+fRw8iFOnLopEkQCystKPHz89fTrvCJnE2BgGBli5EmvWsG7r5xxhbGzsgQMHWDdYDrC2xqpV+N//kJxMTUPVqvD0xMKF1AQQgk2bMG4cdxZfvoSrK3bvxm9LjHl4mEFVFVZWGDgQf/6J48dha2sCXAdEOjrnGzbkvSDzdO+OI0e4MPTTEV6/fn3h9weniYnJ9evXubCvoAwZgvbt4eYGijOws2bh1ClqtaJOnYKmJpo358hcair69MHs2ejcmSOLPDyhoYGtW680NW0xfnz9bt260ZajgHh54c0bfP3KuqGfjlBXV/fLly+SslLsXO4JCMCnT1i1ipoAbW3Mno3Zs+lYDwzExIkc2ZJIMHgwmjaF3Kwj5ikXVK9e/cqVI48ehS9Z4k1bi2JiaQkdHS6KG/+cI7Szs/v69WuXLl1q166dnJy8YsWKAlfobdy4kXVRCoGKCvbuRbNmsLamFqaMGYP163H2LLp04dTu8+eIiABnA+2zZiExkTtzPDw8nNG+PfbuxYIF7Fr56QirVat24MCBFStWhIWFZWVlhYeHF7hKnneE0mNqij17MHQoIiJQrRoFAcrKWLIEM2bAwQFCDreMBgXBwwMypgUtIbt24cABRERAVZULczw8PFwyfTratUN2Nrtf8DxPx169ev3zzz/x8fGGhoanT5/+UhAsalFEOnbEhAno35/dnPdF8L//QUcHu3ZxZzEjAzt3YswYLmzduIFp03D8OPiqEjw8Cknr1qhQAZs3s2ul4DBh8+bNv6f3/JVnz55t2LCBHUmKhpcXTE0xdSo1AX5+mDsX31O3sk5ICOztUbMm64bi4+Hqiq1b0bAh67Z4eHhoYW+PrVvZNVGwI3R0dNTX1y/itNjY2DUcbO5QCAQCbN2KixexbRsdAc2bo0ULLiaccwkK4mLXREYG+vTBpElc52/j4eHhmB49ou/enW5l1TEyMpIlExRyjZZDtLVx5Ai8vXHnDh0By5dj9erSJiuShuvXkZTE+uIgQjB8OBo0wLRp7Bri4eGhS2Zm5qxZPQlp9eDBhJYt/5fOztAW7wg5ol49rF2Lfv3w+TMF6zVrYvBgLFnCuqHAQIwbx/rCnPnz8e4dNm1i1woPDw91bty4IZE0ApwAJ4nE5tKlS2xY4R0hd7i6wskJAwdCzFYB+aJYsAChoexm/kxIwKlTcHNj0QSAgwexezcOHeKXifLwKD5NmjQBHgFxwFNCHjdr1owNK7wj5JQVK5CdzUVk9ju6upg2DXPnsmhi40Y4O0NPj0UTd+5g3DgcOwYjIxat8PDwyAmVKlXasmWxtnZvNTUXYHtcnAEbVnhHyCnKyti/H8HBOHmSgvWJExEVhatXWelcJGIxuWhycnJWVtbDhx+dnLBpE6ysWLHCw8MjhwwfPiw5+UlGxt3OnVv16IHsbOZN8I6Qa4yNERqKESPw4gXXptXUsGQJpk9nJQPqsWMwM2PFRUVGRurp1U9KynRwcLW1PdqnD/MmeHh45J+//4ayMit5skroCIVCofTVWXnyYW+P2bPRvz8yMrg2PWgQef7cS0/PpmHDjo8fP2aw59xlMmzg6TlLLA4AKgHbzp9nosgvDw9PGURZGZcv4+pVrF7NcM8yOML09PSUlJTc1927d3/48CHDWsoTEyeifn388QfXdk+e/Ds9Penr16hHjwIHDpzEVLePH+PJEzg5MdVfHkQiJSAHACAWCASs2ODh4SkLNGiAP//EzJmIi2Oy24Id4YABA7Zs2ZKv8eDBg7Vr1yYUCwspFlu24N69bCcnfycnz8OHj3Fj9NWrt5mZdoAAsPj8mbHqJuvXY9QoVpZxHj+O+PiVQuFUIFEo7B8QQKmaBg8Pj3wwYwb+z96dx1OZ/XEA/9xruZaUspMKk1IkSUVKWlRK0aaFmhaU0jbatIyWmbapaBet2heJtOmHNGmRlNapVBShHdnd5/fHncro4i7P3Tjv1/zBc89zzvfeka/nec75Hhsb2NuDxq2SuCTC8vLy06dPGxgYAMjMzDx//jzneL9+/d69e/fq1SvaBm/YlJXRvn1gZOTHyMipkyeHXbwYK4ZBhwxx1tLaxmTuUlCYpKrqXFFBQ58FBThyBN7eNHRVVUUFFi7EzJk4e7bD16/P1dVVHj9O8PT0oHkYgiBkTXw8SkowbBhtHXJJhO/fvy8rKzM0NARw9erVWbP+vYGmra3NYDDev39P2+AN3s2bV4DFQKfPn/2ioxPEMGKLFi1u347esoWKjBzTuvXyceNoWNS4fz/69oWBAR3xfZORge7d8c8/uHsXtrZQUlJSVFRs0qQJnWMQBCGblJQQG4voaNpqkHJJhIqKigA4jwM/fPjwfceJvLw8iqJqr0FK8KVbN2sFhV1ANoMR3rSpjXgGbd68ua+vj7Oz04kT+PABkycLe4chJITmaTLR0ejSBcOGISIC6up09kwQRP3QtSvmzsXUqcjKoqE3LolQQ0NDX1//r7/+un//flhYmLy8/PHjxwEEBQVpaGg0l8jGevXUzp2rp0zJsrHxmTOn5549rgcPinV0ZWVEReH1a0yZIngu/N//QFHo0YOekCoqEBgIPz+cPo0FC0BmxhAEUZP169G6NWxtaeiK+2SZVatWRUZGdujQQVdXd9WqVe7u7s2aNVu3bt3ChQsVSWEr+qiqqm7fvvrWregNG7wuX8bSpQgMFGsAKiqIjsarV/DyEjAXbtsGPz96Mtbr1+jVC3fu4M4d2NnR0CFBEPXbtWt4/x6TJgnbD/e1gBMnTrS1tc3MzHR0dFRQUFBXV09NTbWzsxs0aJCwAxI1MDNDUhIGD8abN9i5E2JbpamigrNnMWgQfHywaxd/Ke31ayQm0rPr79mzmDIFU6di2TKR1+wmCKJ+UFfHsWNwdcXw4RAmOzHq5XIId3d3Nze30aNHSzoQvhUWYtQoKCjgyBGoqIhv3K9f4eyMNm0QElJbLiwuLlZQUPheS2HxYnz9CiE3pqyowKpV2LsXhw+je/cam+no6KSlpeno6Ag1GMGbgoICNTU1SUchGxrgZyVtb3nCBBw/jqwsLoWO2Wx2ZWWlgoJC7T2Qv72lS6NGiIqCjg4cHZGXJ75xVVURHY379zF7Nq8F2MrKsGcPpk0Tatw3b+DoiGvXcOtWbVmQIAiiJvv3Q1cXDg6C90ASodSRl0dICAYOhJ0dzdUTate4MS5cwM2bmDuXp/bHj6NDB7RpI/iI//sfunVDnz64eBHkSo8gCIHduIGnTzFzpoCnk0QojRgMBAYiIAA9eyIpSXzjNmmCixdx7RpPuVCY4qKVlQgMxKRJOHoUgYHkoSBBEELR0UFYGLZtE3B3HfIbSHpNmoQDBzB8OKKjxTdokya4dAlXr8Lfv7ZmqanIzhbw6XRuLgYMwN9/49Yt2NsLFiZBEMR/eHpi0CAMGICiIr7PJYlQqjk5ISoKPj7YsUN8g6qrIzYWCQmYN6/GNlu2wNcXcnK89vnkyRMrK6cWLWzGjv3T2hrdu+PSJXI7lCAIOkVGQk0N/frxfWIdifDr169hYWEfPnwQMC5CaDY2uHoVmzZh6VKR7CPIFScXxsVxX9f46RMiIzFxIh8djhw5/e7dTa9f3zh27NG8efHkdihBELRjMhEfj5s3sXYtnyfW/vKHDx+8vLwyMzMFD40QmokJrl1DbCwmTkR5uZgGbdoUFy7g1CmsWFH9pbAwDBkCbe26O6moQFISVqzA06dfgPaAHJNpr6Ii9i2JCYJoGNq2xcaNCAjA3bt8nEX+LJcNWlpISEBREQYOxJcv4hv0f//D8eNYterHQTYbO3fWMU3mxQvs2oVRo6ClBU9PZGXBzq6HquosBmO/hsauAQNEsMM0QRAEAGDmTHTvjt69wfvuOmSXeZmhpIQjRzBzJnr0wLlzEE/NV21txMXB0RFMJgICAOD8eYamJmx+qhCem4vERFy+jPPnUVEBe3v07YugIOjrAwCbveHUqYiXL9+4u5/mbGxCEAQhInFx0NbGkCE4e5an9nUkQnl5eT09vTqX5RPiISeHbdsQHAxbW2rcuNCUlCu9enVesMBPXpQF2b7nwry8lw8fzn340G7pUm+gCYDCQty4gcuXcfkyMjLg6Iju3eHtDWvr6p0wmcyRI0eILkiCIIjv5OVx8SK6dSu3tw/466+xdnZWdbSv/WV9ff3s7Gz6wiNoMGsW7t8PX7fuJkUtS0raU1KyYeXKBSIdUUcHR4++s7QcBswA3v/2W+f8/GeXLyM5GTY26NsXISGwsiLzXwiCkBY2NlBWdrx+vfelS9p1FvEnv7pk0ufPSRTlC7QpKpp9+PC1z59FPuLFiycYjC7AZGBBSUmjt2/fBgQgNxexsViwANbWJAsSBCFdiouzgRX79unX2ZL89pJJTk62amo7gXQWa7OCgl2rVnBxwYEDoD0j5uRg507064cVK6wpKgn4CmQxme82btRxdASLRfNwBEEQdJGTqwSezplT9w5zJBHKJC+v8cuXd3JwCAgI0HjwwP/tW3h74/JltGoFe3sEB+P9e6H6f/0au3bBxQVt2iAmBp6eePu265QpfeXlzVgs+61bA5nkApAgCOl2+PAmZeWBLNb/6mxJtmGqV4qLcfkyTpxAVBTMzTF+PNzd0aQJr6e/eoUzZ3DiBJ48gbMzRo6Ek9N/LvuqbcMkTmQbJnGStn12pFkD/Kxk6C3zuA0TWT5Rrygrw8UFLi4/MuKiRbCzw8iRGDr034xYXFycnJxsaGhoZGTEOevFC0RH48QJ/PMPBg7EggUYMABkpjBBEA0ET4nwzZs3aWlplpaWBgYGog6IoMXPGXHOHNjZwdn549q1AwoKbJnMNG9vDwWFySdOoKAAbm4IDESvXpDExR5BEIQkcf+15+7ubmZmFhgYCCAhIcHZ2bm4uJjFYh09etTV1VWsARLC+Z4RCwoQHY11605nZHgCfkDpunU9f/tt8r596Ny5tl3pCYIg6jcuUx4qKirOnDlj923lxaJFi1q3bn316tVx48bNmjWrsrJSvBES9FBTw9ixWLBAVVmZU0K9oFUr+TVrYGNDsiBBEA0al0T48ePH0tJSExMTAHl5ebdu3VqwYIG9vf2qVasyMzNJAW6ZNmLE8M6d72prO+jq9t258w9Jh0MQBCF5XG6NcibYlJaWAjh//jxFUX369AHQrFkzAB8+fPg+yYKQOQoKComJkQUFBaqqqmQJBEEQBLheETZt2tTAwCAsLOzLly9hYWFWVlacOesZGRkAtLS0hB81KSnJy8vLy8vr2rVrXBvk5OQsXLhw3LhxYWFhbPaP5ZBRUVETJkzw8/N7/Pix8GE0WGpqaiQLEgRBcHD/bbhq1arg4GB1dfWkpKTFixdzDp49e1ZLS6tFixZCDnnnzp0BAwZYWVlZW1s7OzunpKRUa1BeXu7g4PDp0ydXV9fg4OA//vj3Dt7Jkye9vLz69eunra1tb2+fk5MjZCQEQRAEwX3W6K+//mptbZ2amtqxY8cOHTpwDurr62/evJkh9MyK4OBgHx8fX19fABkZGUFBQeHh4VUbREVFMRiMnTt3MhgMQ0NDFxeX+fPns1isv/76648//vDw8ACQmpoaFha2ZMkSIYMhCIIgGrgaV41ZWFhYWFhUPUJXoZYbN25s3LiR87WDg4Ofn9/PDXr27MnJuF26dCksLHz+/Hnbtm1v37598ODB7yfGxsbSEg9BEATRkNWYCN+8eRMWFvbw4cPi4uKzZ88COHv2rJqamoODg5BD5uTkaGhocL7W0tJ6+/btzw2af9t2lslkamhovH37VkNDo7KysvYTv7ugeSHiQcT4peM53zIYDBUVFSHDJiTry6Qvbfa0Ef6GBMELiqLIR82jBvhZydZbttSxTJiYUHsb7okwOTnZycmJc2fywwfOsjM8ePBg//79ws9SUVJSKisr43xdUlLyc4pSVlYuLy///i2njbKyMoDaT/zOqNBIXk3+ezZVVFQkM11pUVFRwWQyJTLRZuvWrR4TPVRVVcU/dANUVlamqKgo6ShkQwP8rGToLVMU1bpZ6zqbcU+EU6dOtbKyOn36dGpqKueZHIDBgwcvWrQoLy9PW1tbmMgMDQ2/L0bMzMz8uWybgYHBgwcPOF8XFBR8/vy5efPmTZo0adSoUWZmJmcKK9cTv2tT0satY0Msui1qEiy6vc9j39JdS0nRbfGQoarKEtcAPysZesucott1NuPyp/2nT5/u3LmzcuXKJk2aVL3+bdmyJYCsrCwhIxs2bFh4eDhFURRFHTx4cPjw4ZzjJ0+ezM3NBTB8+PDY2FjOpNDDhw936tSJM1V1+PDhBw4cAFBSUnLixInvJxIEQRCEwLj8ac9ZSv9zwv/06RMA4a8GfH19T5482aVLFwaDUV5ePn36dM7xX3/99fTp0/369TM3N584caKNjY25uXlKSsrJkyc5DZYuXero6PjgwYOcnBwTExM3NzchIyEIgiAILllNW1tbS0vr3LlzHTp0qHpFePToUVVVVVNTUyGHVFdXT05OTk5OBmBjYyMnJ8c5/vjx4++r9Tdt2jRt2rSsrCwrKyt1dXXOQRMTk6dPnyYnJ6upqVlaWsrQ01qCIAhCanFJhEwmc86cOYGBgZWVlQYGBmw2+9GjR8eOHVuzZs2cOXNYVfdpFZScnFy3bt2qHTQ0NKz6ramp6c9JV0lJqUePHsIHQBAEQRAc3O9zLliw4P3794GBgRUVFQDat2/PYDAmTJiwYsUK8YZHEARBEKLFPREymcwNGzbMmjUrPj4+JydHXV29Z8+eZmZmYg6OIAiCIESttpkvLVq0mDBhgthCIQiCIAjx454IMzIyalp7YWxsLMp4CIIgCEKsuCfCrl27cpb0/YyiKFHGQxAEQRBixT0RhoaGlpSUfP+2oKDgypUrZ86c+fPPP8UVmFC+fv2akJBAKssQhABevHixceNGV1fXvn37SjoWghAHBu9XeGvWrLlw4UJCQoIo46GHufmihw8T9PVLs7LuSDqWekWCJdZ0dHTS0tJIiTVRi4mJGTx4CjAYSHRz6xwRcUjSEUk7Gao3RhcZesucEmsKCgq1N+OjevK4ceOuXLny4sUL4QITh/x8F+BadvaH70W6CYKoXVERHjzAr79uBYKAUCApMjJR0kERhDjw8af9ly9fABQVFYksGNqw2XJAKcCQlRLpBEGjsrKyrVvDHjxInzJlhJ2dbbVX2WxkZeHFC7x4gZcvf3zx+TOMjFBcbAF8BAB8pqg2Hh7YuBHCldknCGnH06zR8vLy9PT033//XV1dXfgSa2LQpMm+rKwLwKrVq7FokaSjIQjxmjTpt4gIzeLiwVFR8w4c2K6o2IGT7Tj/PXkCRUUYG//7X8+e+PVXGBujVSswmfD3992woScQBmS7uKy8cgW6umjbFn/9BWdnSb8xghANPmaN6uvrh4eHy8Q1lpnZu8mTVy9fPnrxYnz5gjVrJB0QQYhRXNyN4uJkAB8++E6e/L9u3ToYGf2b84yMYGQEJSXuJx48iKNHWzVtmnn06PXffrP291fs2RM3b2LuXLi4QEMDM2ciIACS2I+SIESIp1mj8vLyzZs3t7CwUKrpH5CUkZOT09fH+vVYtw7r1+PzZ+zcKemYCEL0Xr3C2rV4/74Vg3GWonqqq0ceOjS1d2+ezo2Jwbx5cHKCri5sbc29vRV37EDPnujaFdeu4dMn/PYb/vgDK1di6FBs3UrulxL1CFUfjRo16siRI2w25eRETZpEyclREyZIOqZ6oaioqLy8XCJDa2tr5+TkSGRomZCeTs2cSWlqUjNnUvfv57q6TjYz6xUUtJPH02/coLS0qIsXqaZNqcxMKj8/Pz+f0tCgfv7Iw8MpY2OKwaDMzKizZ2l+F7IoPz9f0iGImwy95crKyrKysjqb1ed7HAwGdu1CVBT27MHhwyDbFxL10sOHGD8ednZo2hRPnyI4GObm2qdPhz16FD9rlg8vPTx6hKFDceAAbt2Cmxs428CoqWHYMOzZU72xhwfS05GSAm1tDBkCDQ0sXQo2m+53RRBi9OPW6Pnz59fw8DDtypUrooyHZi1b4vffsXMn4uPh6AhHR8THSzomgqBJWhr++guxsfDxwT//oEkTQTrJyoKzM9auhaMjJk1CbOyPl6ZPx5AhmD8f3/YM/cHKCgkJyM/HvHnYsAFr16JvX2zZUjFp0ojk5EetWxteuXLq+06iBCHlflwRKigoNOKBBGMVjK8vFBWRkoI7d3DjBmyrTyYnCNmTlAQXF/Tvj/btkZ6OwEABs+CXL3B2xsyZmDAB+/bBxgbt2/941dIS+vo4f77G0xs3RkgIioqwZw/++QetW+9JTNQvLr6WltbXwWGEIAERhCT8uCLs27dvvayoxGRi71507QpnZzx6BAsLWFggNRWSqI5CEMKKj8eqVXj1CgsW4NQpCDOJu6gIzs4YMABz54LNxoYNXG6ETpuGHTsweHAdXXl4wMMDenqHcnLWAVqA17NnewWPjCDEqz4/I/zOyAgLF8LLC61a4fFjZGTAwgIVFZIOiyD48fff6NMHU6bA3R3//ANvb6GyYGUlPDxgZPTv4qLTp9G0Keztqzdzd0dKCl6+5KlPd3dbBmM5kAwsp6hheXmCh0cQ4lTbZVFeXt6LFy8KCwurHpTRq8bZs3H6NEJCMHUqnj9Hu3YwMcHjx1BRkXRkBMFNSkrKokUbWSzF9esXPnvWZtUqFBVh3jyMG8fliR2/KAre3igrw/HjYDAAYMMGLFzIpSWLBQ8P7NqF1avr7jYoaE1lpX9kpF/XrhZpaZtatsSFC3BwEDZaghA5rnNJc3Nz+/Tpw3t7acNZPlHt4JMnlKYmlZ5OURT15Qulp0fp6FCfPkkgPNlFlk+IR35+vq5uJ+AucF1BoZOlZeXx4xSbTVv/CxZQXbpQhYX/fpuQQJmaUpWV/wng+9fPn1Pa2lRxMd+jeHhQTCa1Zo2QwUo7GVpLQBcZestCLZ/w9vZ++PDhoUOHBg8ePGnSpAsXLsycOVNdXT08PFyUSVm02rTBb7/B2xsUhcaN8fw5FBRgbIwaNl4kCIl5/vx5ebk1YAl0a9TIMCbm7ciR/166CW/7dkREIDoaqqr/Hlm/Hv7+NdaLMTFBx46IiOB7oPBwBAdj8WIMHgyyjSkhzbj87LPZ7EuXLm3YsGHs2LGampq6urr9+/cPDg4ODAxcs2YNJcs/0f7+yM/H3r0AoKKC9HRoaaFNG2RlSToygqiiSRPT/PwUIJHBuNikSba+vj5dPR89ijVrEBv7oy7MkydISYGnZ21ncabMCGDGDNy8iYQEtGoF8siQkFpcEuG7d++Ki4utra0BsFis/Px8znFPT8+ipNmwAAAgAElEQVSHDx8+f/5crAHSSl4e+/dj0SK8eQMAiop4/BhGRjA1xdOnkg6OIAAAL19iwABVb++Do0Yd9fQ8n5BwgkHTxWBcHGbNQnQ0Wrb8cXDtWsycWWP1UQ4XF7x+jbQ0QQa1tkZ2NpSU0KoVEsm2ToRU4pII1dXVGQwGZ9MlAwODJ0+ecI4XFxdDRrZhqoWZGaZPx9Sp/37LZCI1FdbW6NABqakSjYwggPv34eAAf39s3Wp27Nj2/fuDWlbNWkK4fRujR+PUKVha/jiYlYUzZ+DtXce5cnKYNEnwgr2NG+OffzB8OBwdsXatgJ0QhOhwSYQsFqt9+/bJyckABg4cGB8fv27dusuXL0+ZMqVJkyatW7cWe5A0CwjA27c4ePDHkcRE9O4NG5tQObmWcnKGo0Z5SS46ouFKSECfPti4se7MxK/0dAwdip07qy+QCArCxInQ0Ki7B29vHDuGb7eHBBEejs2bsXgx2c6JkD5cp9AcPHgwKCiI8/XMmTM5d2ZUVVWPHj1K43we0eE6a7Squ3cpXV3q7dsfR7Kzs4F2QAlQxmS2ffr0qcijlEFk1qjonD5NaWlRsbH095ydTRkZUWFh1Y9/+UJpalIZGVxO4TotcMQIascOYYO5dYtSVaVatKDevRO2KykhQ1Mo6SJDb1moWaPjxo2bNWsW5+vg4OC8vLxbt25lZWW5u7uLL0WLkqUlJk+GT5WKxJmZmUymNsACFADjlzwuISYIOuzdixkzcOkSaF+mm5+PQYPg5YXJk6u/tGMHnJ3RogWvXU2bhq1bhY3HxgbZ2WCx0KIFeWRISAvuiTDrv9MoNTU1bWxsmghWzVBaLVuGFy9w/Pi/33bt2lVF5TUwG1jEZL7pzeMebgQhtLVrsWoVEhLQsSPNPZeVYfhw2Npi0aLqL5WWYvNmzJnDR2+cfxPXrgkbVePGePr030eG69YJ2xtBCI97IrS2tu7cufOuXbsKCgrEHJDYKCpi927MnPljVve7dw+WLFEbMkSNzb799i0pRUqIHEVh7lwcOoSrV/HLL7R1m52d7eQ09pdf7Dp12qSpiS1buLQ5eBCWlnynXi8vAddR/IzzyDAggDwyJCSPeyJcsWIFg8Hw8fHR09ObMGFCfHw8uz5uONalC8aPx7d7wFBSUlq5cuWZMwGmpgpk80JC1MrKMGYMUlPx99+gb6EgALi7T798eXJ6evw//9z28Lj080p5isLGjZg3j++eJ0zAuXO0rQicPh3XryMxEa1a4f17evokCAHUWFkmOTn58ePH/v7+V65c6d27d8uWLRcuXCjTiwi5WrkSaWnVq2acPYvUVEFKaRAEjwoLMWQISktx/jwaN6a584yMNxTVB2BVVjrfv//w5wZRUVBRgaMj3z2rq8PN7d+SFLSwscHr15CXR4sWJUZGA1isVubmjtXqGxOEqNW2+0Tbtm0DAwPT09PPnTtnb28fHBxsamoqtsjEg8VCWBhmzsTHjz8OmphgxAhMmkT23SZE4sMH9OsHfX2cOFHHSnbB9O3bQ05uIRCjobHVxcXp5wbr12P+fAE7nzEDO3agslKoCKtq2hTPn6NRI79Xr+zLyh49fNhz0CAP2nonCB7UvQ2TnJycoaGhgYFB48aNKVmur1YTW1uMGFF91kB4OMrK8NtvEoqJqL9evoStLfr3x549otoUc9asdWpqbfz8bl66tL191Z12AQA3b+LtWwwbJmDnVlbQ1salS8IGWU1l5S1gFKACjHvw4BXNvRNErWpLhO/fv9+8ebO1tbWFhUVoaKiLi8vVq1fFFpk4/fknkpJw5syPI4qK2LQJW7aQktwEnR4+hIMDpk1DYKAIR9m2Td7ff+LmzSusrKx+fnX1avj7C7WXk8ClR2sxcKAtg7EASASW5Of7/O9/NPdPELXhurowOjra1dVVUVGRyWT26dMnPDz869evNC90FKU6F9T/LD6eMjCgnj17l5ub+/2gsTFlZ0d3cLKMLKgXxvXrlJ4eJeqiFB8/Us2aUTV9VE+eULq6VFFRHZ3UvmK6qIjS1KRevhQwwprMnj3PxMR+zpyFnP2bPDxo7l9EZGh1OV1k6C3zuKCe+62ZKVOmsFisOXPmeHt7Gxsbizk3S0SvXmjWbKmFxRU1Naara5ddu9YBOHMGlpa4cAEDBkg6PkLGRUdj8mSEh6N/f9EOFBICV1fo6HB/9a+/4OsLZWWhhlBWhqcnQkPxxx9C9VPNpk3rNm369+sRIzB6NBISkJgIIyM6RyGIn3G/NRoTE/Pq1as1a9Y0kCwI4P3792/fXikpSXz3LiEy8t6bN28AmJtj0CB4kCf3hEAqKiqePn369evXAwfg5YWoKJFnwYoK7NiBmTO5v5qbi4gITJtGw0C+vggLQ2kpDV1xNXQo3r6FhgZMTbFtm6hGIQiOGhfU07Xzi6yorKxkMBS/fadQUVHB+er4cXz9iiVLJBUXIavy8vJMTe169AjU0+uxcGFSYiK6dRP5oBERMDb+z/4SVQUFwcMDmpo0DPTLL7CwQGQkDV3VRF0dd+9ixQrMmgVbW8j4tjeEVKtx1lpMTExkZGRWVlZ5eXnV47GxsaKPSgJ0dHQGDmx/7tygDx+YdnbNW7VqxTmupIRVq7BoEfz9oa4u0RAJmRIUFJqRMYfNHgO8MjaeYWp6VgyDbt5c41TnggKEhuLmTdrGmjYNW7ZA1OWHFy2Ciwt694aODqKj0auXaIcjGibuV4Rz584dPHjw+fPnZX33Qb7s3x988+YWZ+egQYP+s/Hab79BR0fw6eZEw0RRFEVxbqsw5OXFsSL1zh28eYMhQ7i/umsX+veHiQltww0divR0PHhAW4c1MTdHTg6cndGnD5fS4QQhPC5XhJWVlTt37vT19d28ebOcMJOsZZCxsfHkydi2DV7/3ZEwIgLduiExET17SigyQtZ06+bFZA5WV4+Rl3+4ceNmMYwYFAQ/P+7rIsrLsXkzTp+mczh5eUyejJAQ7rVM6cVk4tgxREfD3R2XLpEZNATNuFwRvn//vri4eNKkSQ0tC3I4OyM1FdnZ/zloY4PevTFmjIRiImTNq1eYNk3nwoXriYmLnj9P7NnTvu5zhJOXh5gYTJzI/dXDh9GmDTp1onlQHx8cPgyxVeZ3cUF29r8zaMSQfYmGg0si1NTU1NHRqbYTU8PBYsHFBadOVT9++jTev8fKlZKIiZApnDqiixejb1/5du3aNWrUSAyDbt8Od3c0a8blJYrCpk2ClNiuk54eHBxw5Aj9Pdfk+wyaOXPIDBqCNlwSoZyc3KZNm5YuXfrq1SuxxyMV3N1x9Gj1g40aITAQK1ciP18SMREygqIweTI6dcL06eIbtKwMu3bBz4/7q+fOAaB/y1+OadOwfbtIeq7FokVIS8PLl9DWBqlBQwiP+6zR06dP5+TktGnTxszMTEtLq+pL9XXWaFX9+uHXX/HqFb5NHf3XokUIDsaYMYiJkUxghPRbsQKvXyM+XqyDHjuGDh1gZsb91fXrsWABRLQeqm9fFBfjxg1xLA6pql07ZGdjwgQ4OWHMGHbjxr8lJKR6eAwICFgo1jiIeqHG5RMdOnQQZxxSRV4ebm44fpxLhf5jx9C7N1JSYG0ticgI6XbmDEJDcesWWCyxjrtlC5Yv5/5ScjIyMjBypKiGZjDg44MdO8SdCAEwmQgPx6hRcHVdw2Z/BpYuWbK8pKRixQqy7JfgD/dEePz4cZGOSlHUjRs3cnNz7ezstLW1ubZ58eLF3bt327Zt265dO86R9+/fZ2RkfG9gZmamoqIiogjd3eHvzyUROjige3eMGIGXL0U0MiGrHj6EtzdiYmjeZbdOSUn4/LnGmjXr1mHuXFFtc8ExaRKMjfHuHf5780hMXFzAYu0vLr4GaFKU6t69/iQREvyqexsm2lEUNWrUqMmTJx86dKh9+/bXrl37uc3evXu7det28uTJfv36rVq1inMwKiqqb9++Pt9UTYq0c3BAbi4eP+byUmQk3rxBUJDoBidkz8ePcHXFxo3o3FncQ2/eDD8//LwNPYAXL3DlSo1TSemiro6hQ7F/v2hHqYWhoTZwAigFTmZljbG1RWqqxIIhZFJN1bivX7/u7u5ubm7erl07zpGgoKDdu3cLXw48Pj7ewMCAU788KCjIwcGhWoOSkhJtbe24uDiKop49e6asrMzZEWL37t2urq68DCHA7hM/mzWLWr6c+0v+/hSLRcnUhhz0ILtPcFVRQQ0cSM2fL4Gh37yhNDSoL1+4vzp1KrVsmSDd8ru9wM2blIkJVVkpyFjCy8nJadWqq4JCi06d+sXFlXbrRjGZlLExdeqUOEaXoa0Y6CJDb5nH3Se4XxFGR0f36NHj0aNHxsbGX7584RxksViBgYGU0HvzRkZGDho0SE1NDcDo0aOvXLnyser28EBSUhKTyezVqxeAX375xdLS8hxn3htQVFSUlJSUnp4ufBh1cnevcV74+vVQU4Onp6hDIGSDvz8qKvDnnxIYeutWjB+Pxo2rH6co6tWrr8eO0VNiu05dukBdHZcvi2Osn+no6Lx8eaOsLCMl5ZKjo+L168jMhKUlRo1CkyZYuhRscRT2IWQY90cHs2fPHj169P79+xMTEz2+bb7Qu3fvadOmZWdnGxgYCDNkVlaW5beqwDo6OiwWKysrq1mVBVBv3rxp3rz596rfhoaGnL0gADx79mzRokVPnjxp06ZNZGRkM67LpoBPnz5dvnz58+fPnG8bN27szn9JRBsblJTI3b3LtrDgknQPHoSzs9ydO5U1FTiulyorK5lMpqQKsldWVlZWVkpk6FocPMg4d46ZlFQJQMzRFRdjzx65q1erfypXr14bM8YvP79xkyYaTZocqaxUrKGDGgnwUXt7M7ZvZ/TpIxU5R1cXJ06gsBD+/sz16xl//QUPD2rjRrYoJhVI54+lSMnQW2az2bxcNXFJhHl5eS9evDhx4kS1X3mc/JeTk1NnInz+/PlkbjUBt27damFhUVZWJl/l2b2CgkLpf3dzKS8v59pgzJgxkyZNAlBcXDx48OBly5Zt3bqVawAFBQVVrxpZLJaLi4sAhXKGDVM4fBiBgeU/v9SjB6yslNzdGffulfDbrewqLS1ls9kS+TdAUVRZWVmp6Db+Ecjdu8x581jnz5coK7PFH9r+/fI2NmjevPqn4uOzOCfnIqBbWbkyPPywhwffJZEE+Kjd3LBokXJ6elnz5iK/W8MjBQUEByM4GLt2ya9erbBnj1z37uyQkLKWLenM1lL4YylqMvSW2Ww2L7/5uSRCTvL7OYtmZ2cD4GWipp6e3p/c7hNxtnTQ09N79+4d50hRUVFhYaH+f6fZ6erqvn///vu37969s7e3B6D8bTtRZWVlDw+PHTt21BRAixYt3NzcRo8eXWeotfPwwIgRWLtWgeslUHQ0DA1x+LDKlClCjiMzGAyGgoKCvEjnINY8tLKysujmCQsgJwdjxmD3blhbK0kkgF27sGkTl3+SZWXlgDqAigqd4uJSAT60yspKfs9SUYGHBw4eVF6xgt/RRG72bMyejWPHEBDAbN9eqUMH7NxJ23oPAT4rWSdDb5nHP9y5PCPU0tIyMjIKDw/Ht6TIsW3bNi0tLVNT0zo7VVVV7c4N57lg9+7d4+PjOYk2Li7OxMRET0+PEzGbzQZgY2Pz+vVrzqTQoqKi69evcxJhVY8ePeKcJVJWVmCxkJzM/VU9PXh7Y9YslJWJOhBC6pSWws0NPj5wcZFMAHFxqKyEoyOXl+bP92GxBqqo+DdvvmvsWJEtIfzJtGkIC0M5lxsoUsHdHenpuHMH6uro3h2Ghti/H2w2e86c+d26DTx8+LCkAyQkh+sUmgMHDgAYN25cYGCgtrb2yZMnR44cCSA4OFj4aTzFxcUmJiaTJ08ODQ1t0aJFSEgI5/igQYMWL17M+XratGldunTZu3dv//79Bw0axDk4ffr05cuXh4SETJs2TUVF5dq1azUNQcusUY7ff6fmzq2tgbo65eFBy1AygMwa/W7KFGrYMIrNllgAQ4ZQu3Zxf6migmraNP3YsbjCwkLBOhdsWmBFRUXz5kt1de08PPwEHlo83ryh3NwoOTlKTm4hMAE4zWC0P3DggABdydAUSrrI0FvmcdZojcsnQkNDqxZXa9So0Zo1a9g0/bvPzc0NDAycPn36mTNnvh+MiIj4+++/OV9XVFSEhob6+Phs2rSpqKiIczAuLi4gIMDHx2fVqlWcR4A1oTERPn5M6evXNi88MpJiMql//qFlNGlHEiFHUBDVsSMlwV/16emUllaNC3iuXqWsrITqX7DfdJs27WCx5gMlCgrbvL0lsZqETwUFFJPZFigDKOBso0aTZ8ygDh2icnP56ESGsgJdZOgtC5sIKYoqLS29fv16REREQkJCQUEBfbGJHI2JkKIoS0vqypXaGrRpk8lijdPUtNywYRNdg0onkggpirp8mdLTo169kmQMc+dSCxbU+OqiRdS3eysCEuw33ZgxM4DrAAXkWVkNECoCcdHSsgDOApXAVD29dS1aUMrKFEDJy1MaGlSHDpS7O7V+PXXrFvfT169fP2bMuFs1vVxP1b9EWNusB0VFxW7iLyAofdzdcexYbVvyvn7dt7T0j9LS+fPmjbGysnTk+tyGqBdevYKnJw4fRsuWEouhsBAHDuDOnRobREcjLEyMAX0zalT/CxdWffr0m7LyQXd3ZwlEwL+4uMMODqO/fJnaoUPb27e3fS/Qk5qK2Fjcv4/79xEbi4AAlJeDxYKWFpo3R4cO6NEDmza5pKYqUlTvY8eGX7161M7OTqJvhRAc90SYmJhYxm0GSOPGjVu1alVTddD6avRodO2K4OAaCzaWlJQDIwCw2eNPnz5NEmF9VVgIFxcsWYJevSQZxv79cHSEoSH3VzMzkZcHGxvxxgQAcHUdzGIprFp1Tl6+//z54pukIwxzc/MPHx78fNzKClZW/zmSlYXLl3H9Ou7fR0wMDhxASUkOcBNgstkq69ZtjowkiVBWcf/VPmrUqNzc3JrO6d69+/79+01MTEQWlXQxMoKxMeLi4OTEvYGystzXr8eBtsCxoUPXizc6QkwoCuPHw84Ovr4SDmPLltou+KKiMGgQ99KjYjBwYH9j4/5OTqLa9UmCDAwwYQImTPhxREXlY3HxE6AdkPLiRZ+MDEneJyCEwf2fy44dOzQ0NHx9fc+fP5+cnBwVFTV27FgDA4MzZ85s37795cuXLi4uslJZgBacu6M1uX07uk2bLU2aeDEY65jMPmKMixCf5cuRm4vNmyUcxoULUFXFT+uJfoiJwaBBYgzoJ23aoLKyQWzPcuLEZhZrEIPRQknp/pgxkzt3ho8Par6CIKQY1yeHtra269atq3Zw2rRpY8eOpSgqOTkZwNWrV4V/kiki9E6WoSgqO5vS0KBKSupo5uJCNW0qsdLDYtAAJ8sUFBRkZmZGRLANDKisLPGPX13//lQtk/wLCyk1NerzZ2FHEXI2xNixVFiYsDHIik+f8s3MqNhY6t07asECSkODWrCgxjLo9UP9myzD5Yrw48eP169fHzJkSLXjQ4YMiYmJAdC5c2c9Pb2XDeFPvm/09GBujosX62gWEYGKCpHvekOITUREtImJg5XVDHd3p6NHi8W80eDPnj7F3bu17bIbG4tu3dCkiRhj4sbREfHxEo5BbOTksGwZFi2ChgbWrMGdO/j0CaamWLsWMlKGjOB2a5SiKADPnz+vdvz58+fUt7prioqKSkqSKSslKbXfHeWQl8eRIzh4sMZiNIRs+e23P/Ly4j98OAP0Tk8/KelwEBSEqVNRy788id8X5ejdG3Fxkg5CjNzdUV6Os2cBoEULhITgwgUkJKB9exw9CtHvlEMIi0si1NDQ6Nq16/Tp0xMTEzlHKIqKiopaunSps7MzgHfv3r1588bIyEiskUrayJGIicHXr3U0GzQIPXpIrOwWQS82m+JMKGMyFSX+UPzzZxw7hqlTa2xAUTh/XioSobExFBXx9Kmk4xAXBgOBgVi8+Md+Tx074vx5hIUhKAjW1nXfTCIki/tkmQMHDsjJyTk4ODRu3NjExERVVXXo0KHGxsbBwcEAnj175u3tbVVtcnF9p6mJbt0QE1N3y7Nn8eUL5swRfUyEiHXoMFdBobeW1nhj48hRo0ZINpjduzFoEHR1a2yQkgI1NfzyixhjqlmvXg3rotDVFSoqOHXqPwd79cKNG1i3DvPnw94ef/8toeCIunBfPmFqanr//v2TJ0+mpaXl5OQYGhpaW1u7urpyth2ws7NrmEtHOXdHR42qo1mjRti+HVOmwMcHbduKJTJCBM6cwYMH7mlpvcvLc9q1ayfANl40qqzE9u04erS2NmfPYvBgcQVUF0dHnDtX2/Vr/bN8OWbOhJtb9QXHffsiNRWnTsHTE7/8go0bYWEhoRCJmoh60o5E0D5rlOPLF0pdndcpedbWlLEx7SFIWMOZNfrkCaWjQyUni23AOkREUN2719Gmc2cqPp6e4YSfFpiZSWlrS7IoudhU/awcHGqb01taSoWEULq61MiR1MuXVGLiVT+/xQcOHK6UtYnmDWLWKFGTxo3h4IAzZ3hqfOECXr/GkiUijokQgcJCDBuGP/9E586SDuWbzZsxc2ZtDd6+xfPn6N5dXAHVxdAQamp49EjScYjXqlVYtqzGfdkUFeHtjadPYW2Njh2vOTkt2bKl+/TpiUuWrBVvmER1PxLhsWPH9PX1OXu+W1pa6tdAcqFKBV7mjnJoamLdOqxZg4wMEcdE0IpTQaZXL0yaJOlQvnnwAM+ewc2ttjYxMRgwAAoK4oqJB46ODesxIQB7e7RujX37amujpoYFCzB8+IWSkoXAwIKCoJMneZh6QIjSj5vZRkZGI0eObNu2LYAhQ4bk5+dLLirpNXQopk3Du3eoskVVjWbPRkgIBg/G/fuij4ygyerVyMur42mcmG3ahBkz6khyMTEYPlxcAfHG0REnT8LPT9JxiNeff8LVFZ6eUFaurZmDg9nx42cLC/syGJHt25uJKzqCux+JsEuXLl26dOF8vXLlSgnFI+1UVDBgAE6fhrc3T+0vX0arVtiypcH9OpBRly9j+3bcvAlFRUmH8s3794iMrGMpQmkp4uMRGiqumHjj6Ag/P7DZEit8KhGdO6NTJ4SG1nEr29NzzP37zyMi7L9+bd+2LSlQLGEN6SeUJrzfHQVgYIDFi+Hvj48fRRkTQYeMDIwfjyNHYGAg6VCqCAnB8OHQ0KitTXw8LCygqSmumHijpwdtbaSlSToOsfvjD/z5JwoLa2vDYDDWr1+Wnn7jzp3de/Y0u31bXMER3NSYCM+cOdOjR49mzZo1b96cc2TdunVBQUHiCkx6OTvj7l1kZ/PaPjAQurpSscyZqEVxMYYPx5Il6NFD0qFUUVGBkBDMmFFHMykpKPOzBlVr7TsLCzg4YNs2nhrr62PDBkyYgJISEYdF1Ix7Ity/f7+rq6uSktLQoUO/H9TV1V29erXE62tIHIuFIUNwkp96W7GxuHUL4eEii4kQmq8vzMwkvMXSz06eROvW6NChjmYxMVK0grCqhpkIAaxciU2bwONECw8PmJmBPI+SIO61RgMCAmbNmhUbG/vrr79+P969e/e8vLysrCzxRSet+Lo7CsDUFN7e8PZGUZHIYiKEEByMu3cREiLpOKo4d+6SubnjlCn9+/e/UXvLBw9QWQlzc/HExZ/evXH1KhrgH8+mpujfH5s28dp+xw7s24ebN0UZE1EzLokwNzc3Ozt74k97KOjq6gLIy8sTR1zSrW9fpKfj1Ss+TtmxA+rqqHKBTUiLpCSsWYOICKioSDqUbz5+/Pjrr4sfPjz29euujRunVlRU1NL47Fn8tFWMtNDQgKEh7tyRdBySsHw5tm3jdXKAlha2bcOvv6K4WMRhEdxwSYSKiooAin/6H/Lq1SsATSS+xYsUkJeHmxuOH+fvrHPnEBfH63p8QjzevoW7O/buhVTVkM/IyAA6AtpAS6Bl7X99Su0DQo4Ge3e0VSsMG4b1PE8IdXWFpSWWLhVlTEQNuCTCZs2atWvXbvv27RRFMRgMzkGKotauXdu8efNfpKSmr6S5u/O91MzKCqNGwcOjxsIThJiVl2PUKPj5YcAASYfyX2ZmZixWKnBUXn5P06Yf9PT0amr58SPu30evXmIMjk8NNhECWLYMYWF87Fm/bRuOHMG3XX8I8eE+WWbNmjWHDh3q27fviRMniouLt2zZ4uDgEB4evnr16u+psYHr2RN5eXj8mL+zDh2CoiLGjBFNTASf/PygqYl58yQdx0+UlJQmTIjq0OH5kiXvr12LquUf3blzcHSsbYdCievVC9euNdA//vT1MW4c1qzhtb2GBnbuxMSJdSy9IOhXUxHS6OjotlW2TmjevHl4eDh9pVBFS0RFt6uZPZsKDOT7rEuXKCaTSkgQQUCiV5+Kbh84QLVtS335QmOXtGGzqdateSr5PXo0FRpKfwD0VlXu1Im6do3G/qRL7Z9VXh6lqUllZvLR4fjx1IwZwkYlUg2o6PbgwYMfP378+vXr27dvP336NDMz08PDQyypWWYIcHcUQL9+GDAAbm4/9vAkxC81Ff7+iIhA48aSDoWbCxfQqFHdJb8rKxEbi4EDxRKTEBry3VEtLUyejD//5OOUzZsRFYVLl0QWE/GTOirLNG/e3NraunXr1uSO6M+6dkVpKe7d4/vE06dRVgYvLxHERPDgwwcMH47t22EmrSUed+zgqSbf1aswMpKuOjhcNeRECGDBAkRE4MULXts3aYLdu+Hjw+syREJ4pMSa4BgMvhcUcigqYt8+7NuHlBQRhEXUqrISY8fC01PqSlR/l5mJpCS4u9fdUsrni37Xsydu3my4lVOaNsW0aVixgo9T+vZFv37w9xdZTMR/kUQoFM7dUYri+8QRI2BnJxu/xeqZ+fNBUVi2TFXgx1YAACAASURBVNJx1CwkBOPH87SoUaq2pK+Fmhrat2/Qq8XnzsW5c/zNrdu0CXFxOH9eZDERVZBEKJSOHaGsjFu3BDn3/Hl8/gw/v4KkpKSyhjmpTuwiInD6NI4cgZycpEOpQVkZ9uyBj0/dLV+8wJcvsLYWfUx0aOB3Rxs3xty5WL6cj1NUVREaCi8vfPoksrCIb0giFNbIkYLcHQXQqBFGjgzfutW+R4/f1dRMnj17RndoxH/cvw8fH5w8WcdODpJ16hTMzdGmTd0to6IwaBBk5dl9A0+EAPz8kJiIu3f5OMXREW5umD1bZDER35BEKKwxY3D0qIDVFCMiVgDRbHZsWdlKLy/yQEAkCgoKdu/eu2fPMTe38g0b0KmTpAOq1Y4dmDaNp5ay8oCQw94ed+7g61dJxyE5qqpYsACBgfydtXYtrl9HRIRIQiK+I4lQWG3aQEcHf/8tyLkURQGcfcdZFRUNrzKx6JWVlVlbO/n6fvL2flRcPGz8eEkHVKtHj5CeDheXulsWFuLWLfTpI/qYaKKigo4dcf26pOOQqKlTcecObtRRRP0/VFSwfz+mTwep8SxSJBHSQLC5owACA32ZzO5M5khgvZnZRrrjIpCWlvbpk2VZ2dzKyuUVFaUfpXt/5G3b4O0NBYW6W168CDs7qKmJPib6kLujLBYWL8bvv/N3lq0txo2rY797QkgkEdJg7FgcOZJ45MjxfD4X/syfP/fVqyvHjrmHhl7ds8dUsGxK1EJXV7e4+BFQCnwGchtL5/p5AEBhIY4exZQpPDWWrfuiHCQRApg0CenpSEjg76xVq5CWxneVf4J38pIOoD7444/5hYU5v/7aVk+v9717/+Nrgw5DQ0NDQ0MAT59i3DgYG8PGRmSBNjynTjVXVZ3YuLEdi6WwZcs6eXnp/YE/eBC9e/O0Op7NxvnzsrdNga0tHjxAfr6UVvMRDwUFLF2KJUv4e5iipISDBzFoEHr2hK6uyIJrwMgVIQ3OnPlfRcWBsrKAd+/cLl++LFgn69bByQkODnzUqidqt3s3Nm3CrVsTs7NTXr68MXhwf0lHVJtdu3idJpOcDE1N6do3ihdKSujcGdeuSToOSfPwwMePfFdQ69QJkybxtK6GEABJhDRQUpID8gCwWP9oa2sL3M+5czAygoVFAy3VT6/9+7FiBeLj0bKlpEPhwbVrKCyEoyNPjWNiZGMd/c/I3VEAcnIIDMTixXwX4vj9d7x6hYMHRRNWw0YSIQ327v2rRYtBCgpW9vb6PXr0EKar5GSw2bCzoyu0BurkSQQE4OJFmbls2rED06fzuijw7FnZe0DIQRIhx8iRqKhAVBR/Zykq4sABzJ2L169FE1YDRhIhDRwde2ZkJC9blmpiwvPOYzVQUcG9e3j4kKdSkwRXkZGYOROXLqHKNmJS7f17nDsHT0+eGmdnIzMT3bqJOCbR6NIFT5+SUilgMLB8OZYsQUUFf3vQWFpixgxMnixIWUeiFiQR0sbVFadP0/ADamCA2FicPMn3NGsCwKVLmDoV0dFo317SofAsLAzDhqFZM54aR0dj4EBI8aSf2igqols3XL0q6TikQN++RRkZQzQ0upmY2D558oT3EwMC8Pkzdu8WXWgNEUmEtDE3B4uFtDQaurK3x549WLUKJ0/S0FvD8b//wdMTUVEyU4ETAJuNXbswdSqv7WVx4URV5O4ox/btu0tK+uTn33rxYpeXVwDvJ8rLY/9+LFxYuHr1nv37w0sa7KYetCKJkE5DhiAykp6uJkzArFkYPZps1cSrpCSMHYvjx9Gli6RD4ceFC9DSqnsPXo7iYiQmor9Uz36tA0mEHHl5n8rLOfO4DD984O9mcevWFUzmgCVLPk2d+sbWVjbnTUkZkgjpNHQozpyhrbeNG9GrF3r2JNWV6nbzJtzccPgwHBwkHQqftm/nddUEgLg4WFmhaVNRBiRi1tZ49Qrv30s6Dknz9h6rq7ucxVrBYrkGBPD8EwAA+OeffwATNvu3kpJFb982evPmjYiCbDhIIqSTnR3evsXLl7R1ePky9PRgZYWKCtr6rH/u3YOrK3bvlqXamxyZmbh5E6NG8dpe1u+LApCXh709rlyRdByS9ssvv9y/H7t7t6WiYqiLC88/AQAAXV1d4AlQBOSz2Rka0rydiowgiZBOTCYGDeJ7VnTt0tJQXAzhFmXUZ//8A2dnBAfL5NK6HTswYQJPe/BynDsnk2+zGnJ3lENTU3PcuKF9+rQ+fZq/EzU0NNavn9uiRS9VVaeBA1cqKyuLJsAGhCRCmtF7dxSAigqSk3HnDsaOpbPb+uH5c/Tti/Xr+biokh5lZdi3D97evLa/dw9ycjKzJqQWJBFW5eEhyBr5CRPcMzJu3blz4/z5wWQ5ivBIIqSZkxNSU2l+BGJigkuXcPw41q2js1tZl5kJJycsWyarfyKcPIkOHWBqymv7s2cxZIgoAxKXjh2Rm4ucHEnHIR1cXJCWJuAaeVNTuLggKIjumBoeCSRCiqLmzZvXrFmzpk2b+vv7s9nVl5R+/Phx1qxZ3bt3NzEx+fDhw/fjxcXF48aNa9KkiY6OTnBwsHij5hWLhT59EBNDc7cODti8GYsW0TYrVda9eYPeveHvDy8vSYciKN734OWoBw8IOZhM9OjB9w4M9ZWiItzccOSIgKf//ju2bSPz6YQlgUR46NChqKiox48fP336NCYm5uBP9wXKy8s1NDR8fX1fvHhRWWXr99WrV799+zY7OzsxMXHlypU3b94Ub+C8ov3uKIevL/z8MGoUHjygv3PZkpcHJyd4e8PXV9KhCOrRI2Rk8PHA7907PHpUfx4Vk7ujVXl6Yv9+Ac9t0QLjxmHtWloDangkkAj37t3r6+uro6OjpaU1Y8aMvXv3Vmugo6OzbNmyAQMGVDu+Z8+e+fPnq6qqtmnTZuzYsT+fKCVcXBAXh6Ii+nsOCoKtLbp1a9Czz9+/R58+GDcO8+dLOhQhbN0KLy8+CsScO4d+/cBiiTImMSKJsKru3VFSgrt3BTx9yRKEh5MCpEKRQCJ8+vSphYUF52tzc/Nnz57xclZxcXFWVhaPJ1IU9fXr10/fFBcXCx8279TV0bkzYmNF0nl8PLS1YWWFn+4oNwhfvmDAALi6YvFiSYcihMJCHDuGyZP5OKXe3BflMDfHly/IzJR0HNKBwcC4cQgPF/B0LS14e2PlSlpjamBEUrIwLS0tMTGx2kEmk+nr6wvg06dPjRo14hxUU1Or+hSwFh8/fgTA44kPHjw4d+6cv78/59tGjRrdv39fTk6Oz/chuAEDFE+eZPbuLZLqR9euwcxMVVMzoKIi2sLCKCJitwrvE/CFU1xcrKCgIP7tbR89elRYWBgTExcWNqpLl8r580sLCsQcAp1CQxV69ZJXUyvm8V2UlyM2ttHq1V8LCsRUa7mwsFDUQ3TvrnzxYsXo0eWiHkjUaPmshg9nOjmpLFlSKNi/rWnTGB07qvr4FJmaiuMPZDH8eNCFzWYrKCgoKCjU3kwkv9E+f/78/Pnzage/5yFNTc38/HzO11++fOFxAz9NTU0Gg5Gfn8/Z/732Ey0sLJYtWzZ69GhBoqfDmDFYuxbKygqiSBlqaujZ0zc6uhI4l5S0w81tclLSWfqH4UZeXl78ifDMmTPDhs1ksxmTJ//VoUPJtm0TGQxFcQZAu337sHkz1NTUeGz/v/+hbVuYmDQSaVTV8B6eYPr1w/Xr8l5eSiIdRTyE/6wsLWFkhFu31ASrn6emhrlzsWGD6uHDQgbC+4ii/fGgC5vNrjrRpCYi+Y3Ws2fPnj171vRq27Zt792717dvXwBpaWlt2rThpU8Wi9WyZct79+4ZGhrydaJEGBigZUskJaHmj0Eod+6kANuBlsCstLTeIhlDagQEBLHZocB4YN/z58MZjImSjkgoV6+iooK/UnD17L4oh6Mj1gi7a1m9wllQKHAh2ZkzYWqKe/dgaUlrWA2DBJ4RTpkyZevWrenp6S9fvty8efOUKVM4x8eMGXPv3j3O13fu3OF8fe/evZRvZaenTJny559/vnv3Ljk5+ejRo5MmTRJ/8LwT0dxRDicnGwZjHfAY+Ku0dHhqqqgGkgaVlWoAp2zda2Vlmb+A4Kya4HEPXo6zZ+tDQZlq2rZFRQWd9Qhl3ZgxiImBwDcdVVUREIAlS2iNqeGgJGHVqlUtWrQwNDRcuXLl94O9evW6efMm5+suXbpYf9O5c2fOwdLS0unTp+vq6rZu3Xrv3r219D9q1KgjR46ILHyepKVRrVqJsP/Ro6doaFg6Obk7OZUxmdTixSIc67uioqLy8nJxjERRFEVlZ1OurlS7dtnKysaAnIKCflxcnNhGF4W8PKppU+rjRz5OefKE0ten2GyRxcRNfn6+GEYZO5bavVsM44gWjZ+Viwt14IDgp5eVUcbG1JUrdIVTI/H8eNCisrKyrKyszmaSSYSiJg2JkKKo1q2pe/fEMdCePZSCAmVsTOXminYgcSbC48cpHR1qwQKqtJSiKEpLSysnJ0c8Q4vOH39QXl78nfLXX9TUqaKJpmbi+U0XGkp5eIhhHNGi8bM6fpzq10+oHvbupXr0oCmamtW/REhKrImQi4uYCsFMnPjvLSZDQ8FLVEiPjAw4OeGPP3D+PNasgaIiADD4upkoldhshIbCx4e/s+rlA0IOR0fExUk6CGni4oKUFAizq5KnJz58wKVL9MXUMJBEKEIifUxYjYEB0tMxfTo8PODsLKvbNlEUdu2CjQ3s7JCcDCsrSQdEq5gY6OrC2pqPU758QUoKetfT6VAmJlBQwNOnko5DaigpYdgwHD0qeA9yclixAgEBoMS00KaeIIlQhOztkZUl1ukAGzfi6lVcuwYdHcjcDJqXL9GvH/btw5UrCAxEXSt/ZA+/xUUBXLiAnj352KdJ5vTqRUrM/IenJ/btE6qHYcMgL4+ICHriaSBIIhQhzvaE0dFiHdTODrm5sLBA5874/XexDi0wzoVg167o1w9Xr8LMTNIBiUBGBlJS+N4uqh7fF+Ugtdaq6dEDRUVISxO8BwYDy5dj8WJZvS0kESQRipY4745+p6SEhASEhGD1arRrJ+2FSdPT0bs3DhzA1atYsABirP8jVtu3Y8IEKPG8+uPt27eLFv1x6tR6O7v6vN1c796Ijyf38X5gMDBmjCA7FFbVvz/09XHoEE0xNQAkEYpWv35ISZFMKpoyBS9forgYzZvjxAkJBFAnNhu7dsHWFgMGIDERUlwgQVilpdi/n48do0pKSrp1c1m71qi4uJm7e73YhLAGhoZo1AiPHkk6DmkyYQIOHQIP5VBqs2YNAgNRWkpTTPUdSYSipayMPn1w7pxkRjcwwMuXmDABo0dj2DDpqtP98CFsbXHwIK5dw4IFYNbrn8QTJ2BlhdateW3/5MmT4mIrihpLUZM/fWqWnZ0tyugkjHNRSHxnagp9fWE/ky5dYG6O0FCaYqrv6vWvH+kgkbujVYWE4MoVXL4MbW3Bt3qhUUUF1q6FoyPGjkVCAh/pQXbxO02mRYsWwF0gD8hkMDJ4rMcro8h8mZ95egq+GcV3f/6J1avx9SsdAdV3JBGK3ODBuHxZJNsT8s7eHnl5MDeHtTWloeEsJ9eSxWoVIcaJZWlp99u166Wv32ncuKW2toiPx+3bmDWrnl8IcqSl4fVr/ua8NGvWbPnyPxQU3Dt29Dp1aqf4t/sQp969ceWKdN2ukLgxYxAVJXi5NQ4LC/TsiS1baIqpXmsAv4ckrVkzdO6M//1PwmFwZtD07r3+40cjNjujrCx23Dh/sY0+erTf48ehb9+mHD2a1atX7IULaNFCbINL2Pbt8PHhexJQXp6Tn198aupFO7tuoolLWujpQVNTqHmS9Y+WFrp3p+FO0ooV2LgRn+rzdCt6kEQoDhK/O/pd06ZpAKc6vUlJibqZGaZNQ3KySMaiKDx+jN27MXEinj0rBH4BGAxGZ1PThrIf6+PHj1ev3nLkyJWJ/G+YERGBYcNEEJNUIosofubpKezcUQCtW8PVFRs20BFQvUYSoThwEqGQ08BoMW/eLCZzJbCJwXAzMzPq3BmXLsHWFnJyMDSEp6ewCbu8HCkpCA7GqFHQ1kafPrh4ER07YsiQPo0b+8jJhWhphQ4aNJCmdyPVbt9O6dHj14CAxsXFm0+c2MXXuc+fIy8PtrYiCk3qkET4s6FDkZyMnBxh+wkMREgIcnPpiKn+IolQHFq2hKEhkpIkHQdgY2OTmhozefLDoKA+jx6dCg9HejoqKhAXB2dnpKRgxIh/k+KwYTh48D9Pbths9rFjx2JjY6v1mZuL6GgsXAh7ezRtivHj8egRBg/G7dvIzsbx45g1CydPrjl61G3LFsadO+f09fXF+p4lZN++iA8fVgITysv3h4TwVwH21CkMG9YgHqByODri6lWp+EtReigpYehQCL/Rrr4+PD2xejUdMcmgv//+++3bt3U2q88P4flSWVmZn58vuv6dnHDsGMzNRdW/oqKiqqoqLy07dOgQFhZW7aCDw4+tYlNTER6OCxcwZQomTICWFjp1wujRFX5+7QoKLIB37duvP3Mm7u+/ce0a/v4bb96gSxd0747AQNjZca8HxmAwBg5sEBeC35matlBUvF5W5sRgJBkb8/dE9NSphrVprZYWmjdHaio6d5Z0KNLEwwO//Ya5c4XtJyAAZmaYNQtGRnSEJTv69BmVkPBiw4aw2bPr+AdIEuG/du3aNXfuXGVlZRH1T1EoL6fh77saOqcUFBTy8vJo6c3KClZW2LgRAG7cwL59iI/HpEm3KivtgT0AHjzo2r9/kb29SvfumDULZmb8bTPbQEydOun33/0aNerctq1BaOhO3k988wYvX6JnT9GFJo04O1GQRFhVr1749An3/9/encY1cbVtAL+SYCCsCtYAssoPsBhcUFFUBAQtahFRq7WiiEtbrEpFse5IqdraumLdah/Bra6IW21lkU0J7oCCWBXr8goiIAjIlsz7YSqligLJkAnh/D8xk8k9F7ZwM5k552TC3l6uOh07YtYsfPcdfv2VoWStREJCmlR6v6Ki8V9PpBH+49WrV7NmzVrXOm8rl5WVGRkZtUTl/v3Rvz8AREc/GzPmKUUBqOVyC9LTVXkmaEacOdPO0nL71avN/ivh6FGMGgWVHjHRADc37NyJhQvZzqFMOBx8+in272fg44GgINjYIDtbNSfybVBlJShKAFRLJOqNHtxm7kIQ8hk9erSFxQsu157DsTExGaxJ2uB7URRWr8bKlbJcK0dFYezYFsik3FxccPEiamrYzqFkpkzBvn0M3D3V08P8+Vi5koFIrUJiIoRC8PkLOBw7HZ3URo8njZBoqvv3U7Ozj+XmitXVI44eZTuNcjt5ErW18PJq9hvz85GZqbILEL6Hvj6srHDlCts5lMyHH0IoRGIiA6XmzkVqKq5dY6CUMpNKMXkyhgyBuzsqKma8eJHu52fX6LtIIySawdTUtHNn/chIzJnDwIPdKmzVKoSEyHI5ePw4Ro5sxiIVKoOiKKHwUEBA4P79h9nOolwYGVAIQEMDixdj2TIGSimtjAwYGSE6GmfPIioKXC60tbW1tbUbfSNphESzOTlh+nT4+5PVcxp2+jRevcIomRaNaFPj6OvbsSMiKelMevr4r746uWNHJNtxlMhnn+HECWbmaJwxA3fuICGBgVJKKCgIvXqhSxfk52PYsOa9lzRCQhYrV6KoCG+NwiAAICwMK1fKMgqwqAiXLsHTswUyKb2jR2MqKpYBA0tKVhw9eo7tOEqkUyc4OuLkSQZKtWuHkBAVvCh89Ah2dti6FZGRSE2V5Tk+0ggJWaipITISS5bgzh22oyiZ339HeTl8fGR578mT8PBoo4/jOjn1EAh+A8r4/AMDB/ZkO45yYWQxCtqkSSgrw9mzzFRTBmvXoksXAHj4EL6+MhYhjZCQUdeuWLECkyaRJ/3+Y/VqhITIOCnMsWNt8XlR2vLl82bOrO3Y0atHD+nSpV+zHUe5+PhALGbmrjyXi8WLK6ZO/e6jj6ZER59ioCJ7XryAkxOWLMGyZcjKgjyLlZFG2DqkpKTo1DNWOX5fzp6NDz7ADz+wnUNp/PkniotlbGYvXyI5uXmrNakSPp+/aVPYkSPnOZxv27Vrx3Yc5SIQwMsLhw4xU+3YsaCCAp1z5xZOm7YlNbXxoQXK6ehRGBnh8WPk5CAkRN5qpBG2DrW1tRoaGgWvHTx4kO1EAMDhYNcubNmCS5fYjqIcvvsOK1bIeDl4+jQGDYKuLtOZWpUBA3D7NoqK2M6hfBj8dFQsvkpRgYCouHj6uXPJzBRVoMpKjBiBCRMwdSoePYKVFQM1SSNURikpKf379zcwMLC0tNy2bVvdfo3XlOdPZmNjbN0KPz+WVx5WBjExKCzEJ5/I+Pa2OY7+DXw+BgxQ2cca5eHmhmfPcPMmA6UcHOzV1HYBuQLBXje3VrbaZVIShEJcugSxGPV+Ncqrjc3j1BznH5z/q/AvBZzoww8+dDZzrr8nMDBw+vTpAQEBpaWlxa9X1aysrNyyZQv99UcffWRtba2AbE0xZgyOH8c337T1tbDDwrBsmYyXg69eISaGyR/s1svdHXFxbXQMyXtwuZg4Eb/9hlWr5C21Z8/GRYvWXL36Z3a2r6lpK5jTdvr0ryMjjwHo0mXGvXshH32EkycZnoOQNMJ3mn5ieu6LXAWcyKK9RW7gf07k6OiYnp4uFovNzc0tLCzonRKJ5O7du/TXAwcOVECwptuyBT16YPhwjBjBdhSWxMUhLw/jx8v49j/+QN++6NiR0Uytk4cHJkxgO4RS8vPD8OEIC5N3fS5dXd2tW9cAWLUKQUE4fpyZeC3kyZMnERGnpdK/APz118DISP8pU5q3lktTkEb4Tte/vF5YUaiAExlqG76xx93dffny5Xl5eRMnTvz000/pnVpaWhs3blRAHhno6WHvXnz2GW7cgIEB22nYEBaGFStk/yuVXoCQANCjB0pK8OABXv8FSPzDzg4dOiA5+d8V0+QUHIzu3XH6ND7+mJmCLSE3NxewBDQAcLlCI6McgDRCBdJT19NT11P8eWtra6dMmXL79m0zM+b/e7ccZ2dMmICZMxEVxXYUhUtIwOPHeP0XS7PV1ODsWaxdy2imVovD+WdJpmnT2I6ifOhHZphqhHw+wsMREAAPD+Wd1c/GZgBFVQHLAKm6epabm1tLnIU8LKN01NTUDA0No6Oji4qKnj59msjIhLsKsWoV7t1j7Nm2ViQ0FMuXy345GBsLOzsYGzOaqTWjbxMSb/vsM0RF4dUrxgoOHYpevZR3BFRFBUQirqXluXnzKufPr3727KZay6xPRhqhMjp16lRaWpqrq6u3t3dSUhKAjh07jlD6+2/q6jhwAMHB+PtvtqMo0IULePQIkybJXqEtj6Nv0NChiI0lM9k2wMgIffviFKPj4DdsQHi4Mk4RJZWie3fweLh1S2P9+p9++umnpkyfLRvy0agy6tat2/79++vvEYlEkZGtYCbibt0QFITJk3H+PHg8ttMoREgIli6V/XJQIsGpU1i6lNFMrZy5OXR1kZmJ7t3ZjqJ86E9HZX4s620mJggOxty5+OMPxmoyom9fFBYiN1cRH9uSK0KCYQsWgMfDhg1s51CI1FTcvSvX5WBiIszMYGnJXCaV4OGB2Fi2QyglHx+kpCA/n8ma8+bh4UNm5vVmypAhyMrCjRto314RpyONkGAYl4u9e7FuHTIy2I7S8ui5/Pl82SuQcfQNIrcJ30VLC15eOMzooo18PrZvx9y5KC9nsqzMfH2RkoK0NJibK+iMpBESzDMxwfff47PPUFnJdpSWJBYjJwdTpshegaJw4oSMS1WoNnd3pKSguprtHErJ15eZpXrrGzwYTk74/nuGy8pgyRIcPIjYWIV+ME4aIdEi/PzQrRtWrGA7R0tauRJLl8p1OZiaig4dYGvLXCZV0aEDbGyQlsZ2DqXk4YH/+z/cvs1w2Q0bsGMHcnIYLtvcDD/8gP37MVixM96QRki0lJ9/xoEDOH+e7Rwt4+pVZGdj6lS5ipBx9O/h4UE+HW0Yl4sJE/Dfx+kYYGiIb77BnDkMl226I0ewYAE2b2ZhaiHSCImW0rEjdu+Gnx9ez5aqUlaswOLFcl0OAoiOJjcI38ndnTwv8058/u41axxMTHrv23eEwbKBgcjPZ2dOjLg4TJyI+fPx1VcsnJ0Mn/hXfn7+1atX2U4hiwplXfph6FB4eeHrr9Eahn40w5UryMyU9/fFtWvg8WBvz1AmlTNoEDIyUFra1pemeltBQcGuXb9IJOInTyRBQQPHjBmpqanJSGU1Nfz8MyZNwrBhaLExew24fh2envD1ZW1+JdII/2FtbX3gwIEvvvhCkSd99Oj/amqMAA5QoadXZWDQQeZSgxX8mXqT/fgjHBxw+DCTI59Y9+23+OYbqKvLVYSMo38/DQ3064fERHh5sR1FyRQUFHA4VgAfAIdjUlxczFQjBDBoEJydsXo1Vq9mqmQjHj3CwIHw9EREhILO+DbSCP/h5eXlpfAfuM2bfwkNPVlcPFxX93+nTh0WiUQKDqAAmprYvx8jRsDJCZ06VanL2T2UwPXruHqVgeXCo6KwZw8TgVQXPYiCNMI32NraCoWPy8pCqqqqNTSknTt3Zrb+unWwt8fkyfjwQ2YLN6CoCPb26N6d4elymovcI2TT3Lkzf/992bhxGt7eR1WyC9J694anZ4K1dU9TU9eRIydLJBK2E8klNBSLFkEgkKvIrVsoK0OfPgxlUlFkWH2DeDze5ct/REba79jRr7o6Oj2d4fpCIZYuVcRTM9XVsLeHvj4uXmzxc70faYQs69evX1jYtPPnLVR7ZsXk5MVVVbEFBanJyYanT59mO47sbtzApUuYMUPeOlFRGDcOHA4TmVSXgwPy8/HkCds5lI+6uvq4ceOmTx+9RXM1rwAAEVtJREFUdq1aSwzYnT0bhYUMD9t/g1QKkQgSCbKy5F1hUX5sn58AbG0hEODGDbZztKTq6lqgPYDKSuOiohdsx5Hdt99i4UJ5LwdBbhA2DZcLV1fEx7OdQ4lNngyRiPkBuzwetmzBvHkoLWW4cp1+/fD0KTIylGIFKNIIlcKoUSx/RN7SAgImGRiMUldfzOcfGD16FNtxZHTrFsRifP65vHVyc5GfDycnJjKpOjLXWqPoAbsJCQyXHTgQQ4fiu+8YLkvz9MTNm8jIQKdOLVK/uVSzERYVFZWUlLCdohm8vFpHI8zMzHz06JEMb1y69OukpJ+OHnU1Nk68eFGWh2Orq6ulUqkMb2RQSAgWLID8D+gdOQIfH6VenSMuLq6mpobtFADg4YGYGLZDvFtFRQXrK4Z27Ij//Q/TpjF/9bZ2LXbvbmDS4Pj4+Go5pr/78kvExyM1VRFzzRcWForF4saPo1SRqalpYGAg2ymaobaW+uAD6tEjtnM0ZubMmT///LM8Ff78k7Kyol69avYbuVxuZmamPKeW061bVKdO1MuXDJTq14+KiWGgTsvp0qXLX3/9xXaKf1haUllZbId4h7S0tN69e7OdgqIo6osvKH9/5stu2UI5O1NS6X922tjYZGdnN7eUvf0QLteMy7XicM6ePctYwvc7fvz4qFGjGj1MNa8IAVCt6uETHg+enmgVD5HI+Q87bBgcHJRibt/mCg1FcDADo4yfPMG9e3B1ZSBSG0GmmGmKdetw4QKOHmW4bEAAqqrw22/y1tm1a9fNm1pS6QOp9BKHM8vTk4lwTdDE31cq2whbndby6aj8Nm3C1q0sz+3bXNnZSEjAl18yUOrYMXh5yb6QbxtEbhM2hZYWIiIwZw7y8pgsy+ViyxYEB0POe02HDhVTVE+AA+gDFOu3Od5AGqGyGD4cFy7g5Uu2c7Q8IyMsWsTm3L4y+PZbBAUxM+kUmWi7uTw8kJAA5bhlqdScnDB9OgMPc72hb18MH47QUFneK5Vi+XLo6CAx8XMO5xCHE8LhfGJl1ZnL+oCJ/+K0ro8Qm8jAwEBbW9vGxobtIM2TnT3P2Picnt4ttoO8U3Z2to6OjomJiZx1KIqXnr6yS5d9urpNvTCMj4/v378/g1NJNVFtbW15Of/u3bU9ey7n8V7JWU0iUb9+/XsHh2+4XKVeau/ixYu9evUSyD9MhCEZGcutrCK0tGR5UKtFlZaW3r5929HRke0g/6AotYyMFebmR9q3z2SwbE2NXnp6SPfu3/L5LwCkpqb26NGj0R/G6mrdy5fDa2t1jY3PWlntqq2tfPDggUAgMDU1ZTDb+xUUFAC40djoNNVshIcOHeJyuR06yD51J9GggoICgUCgrcjpeF/Lzc21VMBDZgQA4MGDB+bm5hwy4L8xEonkyZMnZmZmbAdRqFb0v0dVVZWmpqabm9v7D1PNRkgQBEEQTaRcH9QSBEEQhIKRRkgQBEG0aaQREgRBEG0aaYQEQRBEm6ZqjfDy5cuJ9eTm5rKdSBU8e/Ys4b9z+iYlJeUxNHA3Nja2vLy8bjMnJyc7O/uNY0pKSgICAt5+b2Zm5v379+s2r169euXKFUZStU0SiSTxv54/f852KCX19OnTRh/KVxmXL19+8no1LKlUGhsbm5+fT29WV1fHxsZWVVWxl64BEokkNja2/ozTt27dunXr3SPTWnSeN8WzsbGxt7cf8tovv/zCdiJVUF5ebm1tvW/fPnrz4MGDFhYWpaWl8lemJ5i4efNm3Z7AwMAZM2a8cVheXp6xsfHbb//kk09CQkLor7du3WpkZHTjxg35U7VZpaWlAAYMGFD3E5SYmMh2KCW1Y8eOAQMGsJ1CQfz8/IKCguivr1+/zuFwQkND6c2EhAR9fX2JRMJeuob5+/v7+vrSXz9+/NjAwCAlJeVdB6vgRE8hISFjyVJvjNLU1IyMjPz4449dXFzU1dUDAwP37t2ro6PDdq5//fDDD9u2bUtISGh1sygooQMHDpibm7OdglAirq6uW7Zsob9OSEgYMWJE3ZobCQkJLi4uyjZTDIANGzbY29tHRUX5+PjMmDFj5syZAwcOfNfBKtgIiZbg5OQ0ZcqUWbNmqampjR07dujQoWwn+ldISMihQ4eSk5MVOWMFQbQdQ4YMmTFjxosXL9q3b5+YmDh37typU6dWVVWpq6snJiZ6e3uzHbABenp6O3funDZt2p07dx4+fBgdHf2eg0kjJJpq1apVXbt25fF4e/bsYTvLv7Zt26arq3vhwoVOSrLEJ0GoHDMzM1NT05SUlJEjR6alpe3bt8/BweHSpUuOjo5isXjjxo1sB2yYp6enm5vb0qVLU1NT1dXV33Ok0l3PEkorOzu7tLT05cuXpYwvACoHkUiUn5+fnJzMdhCCUGWurq6JiYkZGRmWlpZaWlrOzs6JiYlpaWkCgUAkErGdrmHl5eVpaWlaWlqPHz9+/5GkERJNUlVV5e/v/+OPP06ePPlz5ua353A42traL+stulFaWtqsu4/Ozs7Hjh2bOnXqwYMHmUpFEMQbXFxc6AeJBw8eDGDw4MH0pqurqxLeIKQtWrSoW7duUVFRAQEB9Ozb76Kk3wChbEJDQz/44IMZM2asWrXqzp07+/fvZ6qynZ2dWCymv6YoSiwWN/cPzKFDh0ZFRc2cOZP0QoJoIUOGDLl27dqJEydcXFwA9OnTJz09/dy5c67Kusb0hQsXDh8+vHPnTg8Pj5EjRwYFBb3nYHKPkGjc9evXt23bdvXqVQ6Ho6mpGRERMXr0aA8PD6FQKH/xsLCwiRMn0lP4//777+rq6pMmTWpuEboXjhs3Tl1d3cfHR/5UBNGo+/fv1w1vVVNTCw8PZzdPizIzMzMzM0tKSjp+/DiAdu3a9ejRIzY2duvWrWxHa0B5efnUqVM3bdpE/47asGGDSCSKjo4ePXp0g8fzVq5cqdCALczQ0LBPnz7t27dnO4hKycrK8vX17dmzJ71pampqY2PD5XKNjIzkL25lZTV27NinT58+f/7c1dV148aNb9/WLi8v37Fjx/z589/Y36FDhx49ehgbG9N13NzcHj9+LBKJeDye/MHaIC6Xa2pq6ujoyOfz2c6i7LS0tKysrIxeMzY27tWrF9uhWlbXrl09PT379OlDb1pbW/fr18/T01MJ12PKzc21tbWdOHEivamhoeHi4lJeXm5ra9vg8WQZJqIVyM/Pd3BwqJvbgiAIgkHkHiFBEATRppErQqIVqK2tzcnJ6datG9tBCIJQQaQREgRBEG0a+WiUIAiCaNNIIyQIgiDaNNIICYIgiDaNNEKCIAiiTSONkCCIZti9e/fly5fZTkEQTCKNkCCIZggMDHz/0m4E0eqQRkgQLKAo6tmzZzU1NbK9vaCgoKioSP4YEokkLy+v/uof9dEhm77q1rNnz4qLi99/QFlZWbNTEkQLI42QIBgQHh5uaGhYWVlJby5evFhfX3/v3r30ZmZmpr6+fnx8PIDs7Oxhw4YJBAKhUKipqdm7d+8LFy7Qh8XExOjr6yclJdWvvHbtWqFQWNf2tm/fbmZm1qlTJwMDA5FIlJCQ0GCebdu2GRgYvLEM2xdffGFnZyeVSgFIJJKQkBChUGhkZKSnpzdo0KCsrKy6IyUSSVhYmKGhoVAo1NPTMzExOXjwYHV1tb6+fllZ2fr16/X19fX19etmKo6MjDQzMxMKhfr6+t27d6e/U9rs2bP79u176tQpCwsLoVC4YMECmf6BCaIlUQRByI2+bRYbG0tvdu/enc/nT5o0id7csGEDn89/+fIlRVHJyckLFiyIi4vLzs6OiYkZOHCgnp5efn4+RVE1NTWGhob+/v71K9va2np7e9Nfr127lsvlLlq06MqVK2Kx2NvbWyAQZGVlvZ0nLy9PTU1tzZo1dXvKysq0tbWDg4PpzZkzZ2pqaq5bty49Pf38+fP9+vUzNDQsLCykX/3888+5XG5QUJBYLL527drOnTsjIiIkEklMTIxAIJg4cWJMTExMTMydO3coivrtt98AjB8/XiwWx8fH9+/fn8/nX7t2jS7l7++vp6dnZmb266+/Xrx4USwWM/MvThDMIY2QIBhQW1urr6+/ZMkSiqLy8/M5HM6sWbOEQqFUKqUoauTIkc7Ozg2+sbCwkMfj/frrr/RmUFAQvVIxvUlfLEZFRVEUVVJSoq2tPXv27Lr3VlZWmpubBwQENFh55MiRNjY2dZu7d+8GkJGRQVHUrVu3OBxOeHh43at5eXkCgWDDhg0URWVlZXE4nKCgoAbL6ujo0N9mHTs7Ozs7O4lEUvcd6ejojBs3jt709/cHEBcX12A1glAGZD1CgmAAj8dzcXGJjY1dtWpVfHy8np7eggULtm7dmpWVZWtrm5ycXH9d0GfPnh06dCg3N7e8vByAhobG3bt36Zf8/f3Xr18fHR3t6+sLIDIy0sDAYMSIEQAuXLhQVlZmYmISGxtbV8rc3PzmzZsNRvLz8xs/fvylS5ccHR3pUn379rW3twdw7tw5iqI6dOhQv5SxsTFdim5a06dPb8o3XllZefv27RUrVtQtU66vrz9s2LD6H/Dq6OgMGTKkKdUIghWkERIEM9zd3QMDA4uLi+Pi4tzd3S0tLa2trePi4kpKSkpLS93d3enDzpw5M27cOHNzc2dn5w4dOnC5XB6PV/dAikgkcnBwiIyM9PX1raysPHz4sJ+fH71AY35+PoA1a9bUtRyapaVlg3m8vb07duwYGRnp6Oj4999/JyUlbd68mX6JLjVnzpw33kI/NfP8+XMAJiYmTfmuHz16JJVK31iZ0tjYuLCwsG6TkQWcCaLlkEZIEMxwd3eXSCSJiYlxcXHBwcH0HroRamlp0ZdlAFavXu3o6BgfH0+vHiyVSjdt2lS/jp+f37x58x4+fHjx4sUXL174+fnR+/X09AAcP37czc2tKXn4fP6ECRMOHDiwfv36yMhINTW1CRMm1C918+ZNek3jN9DrWufn5+vq6jZ6Fi0tLQAFBQX1dxYUFNCnoL3RuQlC2ZD/QQmCGV27djUxMfnll19yc3Pp6z93d/fz58//+eefgwcPrlvzPTc3t2fPnnQXBBAXF1f3rClt0qRJ7dq127dvX2RkpEgkqlv33MnJqV27dkeOHGl6JD8/v+Li4pMnT+7Zs4e+QKT3u7i4AHhXqcGDBwM4fPhwg69qa2u/evWqbtPY2NjMzOzMmTN1eyoqKuLi4vr379/0nATBMrZvUhKE6pg8eTIAU1NTevP58+f0xdCPP/5Yd8zw4cONjIyuXLlSWVkZExPTpUsXDQ2Nr776qn6dMWPGmJiY8Hi8devW1d8fHBzM5XKXL1+em5tbUVFx+/btzZs37969+z2RRCKRhYUFgNOnT9ff7+3traWlFR4e/uTJk7KysoyMjLCwsD/++IN+1cfHRyAQhIeHP336tLi4OCYm5sSJE/RLHh4e1tbWZ8+evXLlyuPHjymKCg8PB7Bw4cK8vLz79+/7+PhwOJzz58/Tx/v7+9d/ZocglBBphATBmIiICAD1xz84ODgAuH79et2enJycrl270n+G6urq7tmzx8zM7I1GSE/doqam9vTp0/r76eF99T+xNDc3P3z48HsirV27FoBQKKypqam/v6KiYtasWfTdR1q3bt1SUlLoV8vLy6dNm6am9s+tE4FAsH37dvql9PT0QYMG0Z+ILly4kKIoqVQaGhoqEAjogw0MDPbu3Vt3ItIICeVHFuYlCEWrra29d+9eRUVF165d6/pH09XU1GRnZ1dVVRkbG3fu3FmeJBUVFTk5ORRFmZiYdOrU6Y1XS0pKcnJyNDU1LSwstLW131+qvLw8KyuLz+fb2dm1a9dOnlQEoWCkERIEQRBtGnlYhiAIgmjTSCMkCIIg2jTSCAmCIIg2jTRCgiAIok0jjZAgCIJo00gjJAiCINq0/wetZTvzMAsTwwAAAABJRU5ErkJggg==", "text/html": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ], "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "metadata": {}, "execution_count": 7 } ], "cell_type": "code", "source": [ "plot_bandstructure(scfres; kline_density=10)" ], "metadata": {}, "execution_count": 7 }, { "cell_type": "markdown", "source": [ "or get the cartesian forces (in Hartree / Bohr)" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "2-element Vector{StaticArraysCore.SVector{3, Float64}}:\n [-3.5336251758888616e-16, -1.7807673505274671e-16, -4.694043518787138e-16]\n [-4.890859093161763e-16, 2.5304492174484446e-17, -1.1342162044248493e-16]" }, "metadata": {}, "execution_count": 8 } ], "cell_type": "code", "source": [ "compute_forces_cart(scfres)" ], "metadata": {}, "execution_count": 8 }, { "cell_type": "markdown", "source": [ "As expected, they are numerically zero in this highly symmetric configuration." ], "metadata": {} } ], "nbformat_minor": 3, "metadata": { "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.8.5" }, "kernelspec": { "name": "julia-1.8", "display_name": "Julia 1.8.5", "language": "julia" } }, "nbformat": 4 }