{ "cells": [ { "cell_type": "markdown", "source": [ "# Custom potential" ], "metadata": {} }, { "cell_type": "markdown", "source": [ "We solve the 1D Gross-Pitaevskii equation with a custom potential.\n", "This is similar to Gross-Pitaevskii equation in 1D example\n", "and we show how to define local potentials attached to atoms, which allows for\n", "instance to compute forces.\n", "The custom potential is actually already defined as `ElementGaussian` in DFTK, and could\n", "be used as is." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "using DFTK\n", "using LinearAlgebra" ], "metadata": {}, "execution_count": 1 }, { "cell_type": "markdown", "source": [ "First, we define a new element which represents a nucleus generating\n", "a Gaussian potential." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "struct CustomPotential <: DFTK.Element\n", " α # Prefactor\n", " L # Width of the Gaussian nucleus\n", "end" ], "metadata": {}, "execution_count": 2 }, { "cell_type": "markdown", "source": [ "Some default values" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "CustomPotential() = CustomPotential(1.0, 0.5);" ], "metadata": {}, "execution_count": 3 }, { "cell_type": "markdown", "source": [ "We extend the two methods providing access to the real and Fourier\n", "representation of the potential to DFTK." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "function DFTK.local_potential_real(el::CustomPotential, r::Real)\n", " -el.α / (√(2π) * el.L) * exp(- (r / el.L)^2 / 2)\n", "end\n", "function DFTK.local_potential_fourier(el::CustomPotential, q::Real)\n", " # = ∫ V(r) exp(-ix⋅q) dx\n", " -el.α * exp(- (q * el.L)^2 / 2)\n", "end" ], "metadata": {}, "execution_count": 4 }, { "cell_type": "markdown", "source": [ "We set up the lattice. For a 1D case we supply two zero lattice vectors" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "a = 10\n", "lattice = a .* [[1 0 0.]; [0 0 0]; [0 0 0]];" ], "metadata": {}, "execution_count": 5 }, { "cell_type": "markdown", "source": [ "In this example, we want to generate two Gaussian potentials generated by\n", "two \"nuclei\" localized at positions $x_1$ and $x_2$, that are expressed in\n", "$[0,1)$ in fractional coordinates. $|x_1 - x_2|$ should be different from\n", "$0.5$ to break symmetry and get nonzero forces." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "x1 = 0.2\n", "x2 = 0.8\n", "positions = [[x1, 0, 0], [x2, 0, 0]]\n", "gauss = CustomPotential()\n", "atoms = [gauss, gauss];" ], "metadata": {}, "execution_count": 6 }, { "cell_type": "markdown", "source": [ "We setup a Gross-Pitaevskii model" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "C = 1.0\n", "α = 2;\n", "n_electrons = 1 # Increase this for fun\n", "terms = [Kinetic(),\n", " AtomicLocal(),\n", " LocalNonlinearity(ρ -> C * ρ^α)]\n", "model = Model(lattice, atoms, positions; n_electrons, terms,\n", " spin_polarization=:spinless); # use \"spinless electrons\"" ], "metadata": {}, "execution_count": 7 }, { "cell_type": "markdown", "source": [ "We discretize using a moderate Ecut and run a SCF algorithm to compute forces\n", "afterwards. As there is no ionic charge associated to `gauss` we have to specify\n", "a starting density and we choose to start from a zero density." ], "metadata": {} }, { "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "n Energy log10(ΔE) log10(Δρ) Diag Δtime\n", "--- --------------- --------- --------- ---- ------\n", " 1 -0.143534078866 -0.42 8.0 \n", " 2 -0.156029488060 -1.90 -1.10 1.0 1.36ms\n", " 3 -0.156768935198 -3.13 -1.56 1.0 1.11ms\n", " 4 -0.157046205361 -3.56 -2.33 1.0 1.08ms\n", " 5 -0.157055305404 -5.04 -2.93 1.0 1.01ms\n", " 6 -0.157056367778 -5.97 -3.56 1.0 1.08ms\n", " 7 -0.157056405556 -7.42 -4.20 1.0 1.09ms\n", " 8 -0.157056406830 -8.90 -4.82 1.0 985μs\n", " 9 -0.157056406917 -10.06 -5.73 1.0 1.01ms\n" ] }, { "output_type": "execute_result", "data": { "text/plain": "Energy breakdown (in Ha):\n Kinetic 0.0380295 \n AtomicLocal -0.3163465\n LocalNonlinearity 0.1212607 \n\n total -0.157056406917" }, "metadata": {}, "execution_count": 8 } ], "cell_type": "code", "source": [ "basis = PlaneWaveBasis(model; Ecut=500, kgrid=(1, 1, 1))\n", "ρ = zeros(eltype(basis), basis.fft_size..., 1)\n", "scfres = self_consistent_field(basis; tol=1e-5, ρ)\n", "scfres.energies" ], "metadata": {}, "execution_count": 8 }, { "cell_type": "markdown", "source": [ "Computing the forces can then be done as usual:" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "2-element Vector{StaticArraysCore.SVector{3, Float64}}:\n [-0.05568113281378982, 0.0, 0.0]\n [0.05568335021094699, 0.0, 0.0]" }, "metadata": {}, "execution_count": 9 } ], "cell_type": "code", "source": [ "compute_forces(scfres)" ], "metadata": {}, "execution_count": 9 }, { "cell_type": "markdown", "source": [ "Extract the converged total local potential" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "tot_local_pot = DFTK.total_local_potential(scfres.ham)[:, 1, 1]; # use only dimension 1" ], "metadata": {}, "execution_count": 10 }, { "cell_type": "markdown", "source": [ "Extract other quantities before plotting them" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "Plot{Plots.GRBackend() n=4}", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nOzddXhUV94H8HNl3CfuxBNCCBLcXZoCLVAKRUpbli41WurKbre7NZa2VLbuXloo7logeEKIu9u4z1x5/0heSoFASCa5I7/P02efhL3MHM7Mvd97j2IsyyIAAADAX+FcFwAAAADgEgQhAAAAvwZBCAAAwK9BEAIAAPBrEIQAAAD8GgQhAAAAvwZBCAAAwK9BEAIAAPBrEIQAAAD8GgQhAAAAv+ZxQfjcc8/ZbLZOHkxRFCwRxyGXy8V1Efwa1D+3oP655cb697gg/PLLL7VabScPdrlcDMP0aHnAddjtdq6L4Neg/rkF9c8tN9a/xwUhAAAA0JsgCAEAAPg1CEIAAAB+DYIQAACAX4MgBAAA4NcgCAEAAPg1CEIAAAB+DYIQAACAXyO5LgDoWS4GFRnYYgNbYUIUg5wMUvBRsBAlKrD+aowHN0IAcIdmUaGeLTKw9VakcyA+jjAMRUtQkgJLUWJiuDz3Fqhp39RoQz+UMbvrmD8a2TAxlqLE4mSITyA+jsqM6FgTKtQzpUZ2QAA2pw8+rw/WR4ZxXWQA/EWLHf1WyWysYE40s2FiLE2FhYuRSoAsFKIZdLoFFRuYMhM7MACbEoHfGYclKuD07FkQhL5mTx37zkX6jyb2thj8vmT82/G4SnDtI60UOtbEbqxkhv1O91djj/YjZkRhcMIB0HNONLPr85g9dcyMSPz+VPzHSbiSf+0jrRQ60sjurGXGbKUTFdgDqfj8OJyA87NnYJ62aHVkZGR2dnZERERnDrbZbHw+nyCIni6VV9hRw649S5td6In++Pw4XNLpmxwXg34qZ/6bx7gY9PpQYnpkZ882k8kkk8m6WFzQbVD/3Lqp+j/byj5xkq40oUf64cuTcBmvs+/iYtD2GmbdBabRhp7NwJcm4jjEIULIrd9/CEJfUKhnHz1BV5rRK5n4nJiunyebq5inTzHREvTeKCJBfuNXgQsxt6D+udXJ+m+1ozXZ9N469qVB+L3JXX+qO9zIPneatlLo7eHE6FAIQ3d+/2GwhHdzMeif55ixW6lpkXju7eTtfbp1tzg7Br9wOzkjCh/xO7XuAkN71j0SAN7nhzKm/6+uYCEqmk/+LaVbbZtjQ7HDWeQT6fiiA/TfjtJG2APKfSAIvViejh26mcpuZs7dRq7uh7tlCCiJo9X98OzZ5LZqZuI2qs4CYQhAV5hcaNEB+pXzzOYp5BvDCGmn20KvA0Pozng8bx6JIdR/I7WvHk5P93BPEH7zzTeJiYnBwcH33nvvdbbV3b59e2Zm5tatW93ypn7ukyJm4jbq4TR82zQyQuLmdpI4GbZ3Jjk1Es/cRO2qhZMNgJuTo2UzN1FSHjo5mxwS5ObTU85DH44mPh5DLDtEv3CGhpab7nNDEBYUFDzwwANffPFFUVFRVVXVyy+/fM3DDAbDU089pdVqNRpN99/Un1kptPggveEicziLXJ7UU8/0OIaeG4D/OIm89wi97gLsfgxAZ/1ayUzdQa0dhH80mhD12MD8KRHY2TnkyWZ20naqqcOnD9ApbriMfvbZZ3PmzBk1apRKpXr++ec//fTTaw7AWbNmzerVq4OCgrr/jv6s2syO2UoRGDoxi0xR9niH+dhQ7Pgs4ttS5p7DtBPSEIAb+dc55tETzM7p5ML4Hu94ChahHdPJ8WHYsM3UmVZ4MOw6N3xUhYWFGRkZbT9nZGQ0NzfrdLorjtm/f39ZWdk999zT/bfzZ8ea2OG/U4sT8C/H9eCd5hWiJNiRW0mdA2XtokzQPw9AB2gW/e0o/Xs1kz2bHBjQS6M6cQytHUSsH47P3EX9VA73ql3khqupVqu9NIZVLpcjhFpbW9Vq9aUDLBbLww8/vGnTJqwT07VbWloiIyMv/ZqVlfX99993dLBfTZ/YVIM/epr8aDg1JcxhNvf2u38xHD12mhy3hdk41hUkbL/3tFgsnflMQQ+B+ufW5fVvo9Dy4zwHg7aMoyS0vZfP0CmBaPM47I4j/GKtc3UK1avvzZ1Ofv+FQiFJ3iDp3BCEarXaaDS2/dz2Q0BAwOUHrF27dsiQIQaD4cyZMxaLpbKysqqqKiYm5pqvFhQU1Pl5hARB+EkQrs9j1l9g9t5CZKg7WIii5306Aa09S888iO+dQbQNz2FZViqVclUeAPXPrUv1b3ahBYeocAn2+ViCh3ewklMPGy5Fx2exWbuxegfv3ZGEP0y6d+P33w1No4mJiRcvXmz7OS8vLyAg4PLHQYQQwzAXLlxYuXLlypUrKysrv/jii6+//rr77+snWISeOUV/WsQcm0VkqDn+dq8dRNybjI/bRleZoUMCAIQQMjjR9J1UggL7ahzB7Sr2ERLscBZZYmTv3E87aC5L4nXcsLLMhQsXxowZc/DgwaSkpPnz5/ft2/eNN95ACL388suZmZkzZsy4/OBhw4atWrVq2bJlHb0arCxzOZpFK4/S+Tp26zRSzc2N5jW8m8+su8AcmEkEIDOsbMIhWFmGWyaTiRHIpu6ghgVjb4/wlHVAnQxacpDW2NlNU0i3TF70WJ61skx6evobb7yRlZUVHh4ul8tfeumltj8vKipqamq64uCIiAg4dTvJyaCF++kaM7tnpgelIELowb74E+n4xO10rdVDzn0AOGCmsBk7qeGelIIIIT6Ovp9AJMixyTsonYPr0ngJWGvUQ9koNG8fxcexHyYSAo/89224yKy/QB2+lRfp7un8oJPgiZBDFgpN3uoYFEy+O9KDUvASFqEns+k9deyuGWSIiOvS9AzPeiIEbmehUNZuSiXAfp7koSmIEHooDb8vgZ6yg26GybzAz9hpNGcPlSRnPTMFEUIYQm8MI26Pxcdvo+qtnvW044EgCD2OhUK37qJiZdhX4wjSsz+fh5KpBXHYlB2UFlpggN+gWbTkIC3nYW9nujwzBS95cSB+bzI+egtdYYIsvB7PvtD6H6MLTdlOpSixj8d4xwDotYOISeHYrN2U1V8mLwG/xiJ09yHaTrM/TCQ89GHwrx5Pxx9JwydtpyshCzsGQehBDE40bQfVX429N8orTrF264YTKUps9h4K1mADPu/xbLrMyP4wkeR2psRNeaQf/kwGPn4bXWaELLw27/kwfZ3OgabsoIYHYx+M9qYURAhhCP1vFCEisHsPwzr4wJe9cp7ZV8dun05KemuBQ3dZkYI/MwCftB2y8NogCD1CWwqOCcXWD/eyFGxD4ujHiUSlmX0iG+bxAt/0aRHzWRGzcwap5Gxxp25ZmYI/NwCfuJ0uhSy8CgQh97QONHkHNSEMWzfMU0eIdoKIRL9PIXfWsrBnE/A9W6qZF87QO6cTod48FWFFCv7iQHziNrrEAFn4F972hO9zNA40ZTs1JQJ7bagXp2AblQDtnE6M2kKHidGint+DBoDecbyZve8IvXUqmajwxvaav7g3GccwNGk7vXcmkeT9/xx3gSDkks6Bpu+gJvtECraJlGA7phMTt1GBQmxqBJxmwOuVGtl5e+kvx7l/o3mu3JOEiwg0YRu9ewaRpvKRf1Q3wW07Z5ptaPw2alok9rqvpGCbvkrsp0nkkoNUjhaaX4B3a7CiqTvoVzLx6ZE+FRgL4/FXh+LTdtL5ejhJEYIg5EqjDU3YRt3WB/tXpk+lYJuxodh7I4msXXQ1bFIBvJbJhW7ZRd2bjN+d5IPXySUJ+OtD8Snb6QtwwwpNo5yos7CTttNLE/FnB/jgCdZmXixeZ0EzdtJHbyVVnrRiOACd4WLQ/H3U0CDsOd89SRfF4ySGpu2ktk0jBwb41CPvzfLZz9hjVZrYcdvo+1J8OQXbPNIPnxGFzdlD2WFKBfAqLEIrjtB8HHtvlA822Fzujjj8vZHEjJ1UdrNfPxf6+LXY05QY2PHb6Ef74Y+n+0XNvz6UCJdgSw7SjF+fZcDLvHCaLjR4zSJq3XRbH/yzseSsPdSRRv89S/3icuwhLmjZCdvpFwfhD/T1l2rHMfTFWKLVzj4GE+2Bl/hfAfNTBbtlKin2m46jmVHY9xPIefuonbV+moX+ckXmXHYzO2UH9d9h+D2+2PF+HQIC/TaF3F/PvgkT7YHH+72Kefkcs3M6ESTkuii9a2I4tnkKefch6pcKfzxP/euizJU9deysPdTnY8k74vyxwpV8tGMa8e5F5rsyfzzHgLc41sSuOEr/PpWIk/lBk+hVhgdju2eQq08wHxf63Xnqj9flXvZbJbP0ILVxEjkjyh/PrjYREmz7dOKxE/Quf217AR6u1MjO30d/NY4cHOi/52l/NXYki3jjAvP0Kf/qy4Ag7Flv5zGPHGf2zCRHh/rv2dWmrxLbOJlceog63QpZCDxLrYWdvJ1+bSg+zbcmzndBrAw7nEXuqmUfOkb7z24yEIQ9hUVo7Vn6o0LmyK1EP1jHCCGE0KgQ7NMx5K27qCJY8xd4DIMTZe2iH0rDFyfA9RAhhEJF6EgWWWpk5+6l/WTDbfjge4SNQgv20Qfq2aO3kjFSSME/ZUVjrwwhZuyk662QhYB7VgrdsouaEoGt8Y8ZTZ0k5aHNU0k5D03aTjXZuC5Nz4PP3v0abWj8NkpAoN0zYFGVa7gnCV+Zgk/bQescXBcF+Le25WMS5Njr3rwDWg/h4+jL8cS0SGzE79RFnY/ftkIQutk5DTt8M3VLNP7VeEIAJ1cHnsrAp0dit+yiLP7R8AI8EMOiuw/RJI59MsYf5s13BYbQ2kHEy4PxidupbTW+nIUQhO70UzkzfSf15jD8xYE4nFrX9/owIlWJ3b6HcvjX8DTgKR46TtdZ2R8nEiRcBa/rrgR88xTy/qP0azk+u0IUfAXcg2LQUyfpp08xe2aQ82KhVm8MQ+ijMYScjy06QFN+N20JcOyZU/SpFvb3qaQQmm06YXgwlj2b+LWSuXM/bXZxXZoeAJdsN2iyoak7qFwte2oO2V8Nj4KdRWDo2wmEjWbvOQyLkYLe8+/zzNZqdsd0Us7juijeI1yMHc4iVXw0dDNV4HO7GEIQdte+enbwJmpcGL5tGhkAQ2NuEh9HGyeRNRb2/j/8Z84S4NL6PObLEmb3DDhbb5qAQP8bTTyVgY/fRn1R7FPNOBCEXedi0HOn6WWH6K/GES8NwqFXsGtEJNoylczXsY8ch95C0LPez2fevcjsm0mEibkuitdalojvn0muu8AsPkgbfaWZFIKwiwr17MjfqRwNe+42cmI4ZGC3SHlo2zQyu5l97ARkIegpHxUyr+cy+2YSkRI4YbslTYWdnE3KeGjAr9Rhn9i8CYLwptEsWp/HjN1KrUjBt04j/W2V+h6i4KNdM8gjjewa2LAJ9ICPC5lXzjP7ZhJ9/HJBbbcTkeiDUcS7I4lFB+jHs2mbl8+DgiC8ORd17Ogt1JYq5vgs8m8pUHvupOSjPTPJww3smmzoLwTu9HEh86/zzP6ZRLwcUtCdZkZh528jG22o/6/UoQYvPmvhUt5ZVgo9e4qesI1anoTvu4WEM6onKPlo9wzyaCP78DHIQuAe7+Yz/86BFOwpgUL0zXhi/XBi6SF6+WG6xc51gboEgrBTfqlg0jZS1RaUO5f3txQYFtODVAK0ZyZ5TsOuPApzKkB3vXmBeSuPOXgLpGDPyorGLs4lAwSo30bXe/mM180MhiC8gVMt7Lit1Cvnmc/HEt+MJ0JFXBfID8h5aOd0stTALjlIu7ztjAKe46Uz9KdFzKFbCFj4vhdIeejNYcT+meTmKibjV2q7Vy3JBkHYoXw9O3cvffteemkifmYOOT4MzqXeI+Wh7dNJkwvN3UvbYfQMuEksQqtP0Fuq2cNZZASMEe1FaSps9wzytaHEkyfpsVupo14yphSC8Bpyteyd++mJ26iRIVjxfPLeZJgjyAEhgTZOJmQ8NH0nZXByXRrgPVwMWnaQPtPK7r8FBnVzIysay72dvC8ZX3aInrSdOuDx42ggCP9ibx07cxc1Yyc9JAgrW8Bbk46LSK7L5Md4OPp6PJGhxsZtpRqsXJcGeAMLhWbvofROtHs6qeRzXRo/hmNoaSJeNJ9cloj//Sg9/Hfqp3LP7TuEIEQIIb0TbbjI9NtIPXqCnheLly8g16TjEohAD4Bj6O0RxIJ4fNQWH1zhELhXsw1N2k6FibBfJxNwC+sJSBwtTcTz55HPZODv5TNxP1L/PMd44Kbcfv1lcTFoTx37TSmzo4aZHoW/O5IYF4ZBI6gHeiYDj5KgiduoHyeRY0PhIwLXUGRgb9lFL0nAXxwEXRmeBcfQ7Bh8dgyeo2U/LGDSN9LDg7HFCfjsGFzsGRHkGaXoXRYK7a1jNlWxW6qYZCW2OAF/dyRPDSvwerbFCXiEBLtjH/XGMGJJArRkgL842MAu3E+9OpRYlgjfDc+VocbeH0W8OYzYVMV8Xcrcf5SeEonfFoNNj8K5XQPdX4LQQaNTLeyhRvZgA5PdzA4Lxm6Nxl8eTMKqg15kQhh24Bby1t10gY79VyYBt/2gzadFzPOn6e8mkhNgaLc3EJNoUTy+KB7XOtDvVcwvFewDx1ypSmxSODY2DB8ZjEl7fXssjGU9q7k2MjIyOzs7IiKiMwfbbDY+n08QV+6tSbOozsKWGlG+ns3Tsmda2Xw9m6bCxoZi48PwCeEY9P+5hclkkslkvfymGge6fQ+lFmBfjSdk/r2fHCf171EoBj1xkt5ew26dSiQqejsFof7dxcmgP5rYg/XMwQb2rIaNkWKZgVg/NdZPhSXIUYwU413rOd+N9e/dgXCgETvUwmIYbXEhowsZnajZzjZYUb2VDRJiCXKUqsT6q7FlSfjAAAy2ovYNAQK0Zyb58DF6xO/UpilEAqwY4q80DrRgH8XDUfZsGCDq3fg4mhCGTQgjEEIUg/J07JlW9qKO3V3LlBlRnZUNEmLhYhQiwhR8JOOhB/riaSp3nvjeHYQCHCl5CMexCDGS85GCj4KFeJgYhYsxAcSe7+Lj6H+jiQ8LmdFbqE/GkFnRkIV+52wrO38fPS8W+/cQgoDP34eQOBoQgA0I+PNDpRjUaGPrLKjZzuqdyOJCbm879e4gHBnMjo/ECAK6x/3RyhQ8Q40t2E9nt2BrB8HV0I98WsQ8e5p+fyQxNxbOfd9H4ihSgkVKEEI9dZLD1wh4seHB2Ok55Ilmdsp2ygMnJwG3M7vQ0oP0W3nM4SwSUhC4C3yTgHcLEqJd08mJ4XjmJmpnLWShL8vRskM2U3wCZc8mk3t9aAzwYRCEwOvhGHp+IP7jRPL+o/SjJ2CRbh/EIrQ+j5m6g3pxIP7JGMJDZmEDnwFBCHzEmFDs3G1kvRUN3UzlaOHR0HfUWNhpO6iNFUz2LHJhPFyygPvBtwr4DpUA/TiReKI/PnUH9cp5z13hF3TelyVM5iZqfBh+KIvsI4PmUNAjoIkB+JolCfiEMOy+I/SmSubTsUR/NVw9vVKthb3/KF1nRXtmkPAhgh4FT4TAB0VKsB3TyVV98Sk7qBfOQK+hl2FY9EEBM+g3angwfnI2pCDocRCEwDdhCC1Pws/fxivUo/SN1N466DX0DjladtQW6rsy5mAW+fxA/JprawHgXtA0CnxZmBj9PInYVsP+7Sg9JAhbNwyHZdY9lsGJXjxD/1DOvJJJ3JsMa6qD3gO3W8D33RKFXZxLpijQwN+of59nbBTXBQJ/xbDo0yIm9ReXnUb583j3QQqC3gVBCPyCiET/GEycnE2e17Cpv1DflTHQVOoh9tezmZuor0qYLVPJD0cT3O5LB/wTNI0CPxIrw36aRBxtZNdk0+suMK8NISZHwLMHZ85r2GdO0aVG9J8h+DxYLw1wB4IQ+J3RodiJ2eTGCuaBY3SEGP0rkxgZAnHYqwr17NqzzKEG5vmBxN9SYEQM4Bh8AYE/whCaF4tfnEsuTsTvOkjP2Ekda4K20t5QqGeXHqTHbaMGBGClC3gP9IUUBNyD7yDwXySO7knCi+aTt/fBFx+kJ22n9sAsix5ztpW9cz89fhuVosRK7uA9nYFLoEEKeAb4JgJ/x8fRihR8eRL+XRnz2AmaT6DH0/F5sfCk4h4sQrtr2XUX6EI9ejQd/2QMz+27qgLQTRCEACCEEImjpYn4kkR8WzW7Po9+8iTzQF/83mQ8SMh1ybyWhULflDLv5DF8Aj3aD18YD/cWwENBEALwJwyhrGgsK5rM0bIbLjLJP7uyovCVqfgoGE1zMy7q2I8KmW9LmXFh+LujiAlhUHvAo0EQAnANGWrskzHEG0OJz4uZFUdoDKF7kvG7EvBQEdcl82BGF/qpnPm8mKkyo3uSsHO3k1Gwjg/wBhCEAHRIJUCPpeOPpeNHGtnPi5m+v7hGBmOLEvDZMTDQ408uBu2uY78rZbbXMBPD8acz8JlROAEJCLwHnM0A3NiYUGxMKLGBIn6rZL4tZR74g54Wic+PxWZE4X67WzrFoAMN7M8VzKZKJkmBLYrH3xnJg3VhgDfy15MYgJsnIdHiBHxxAq5xoN8qmY8KmXuP0BPC8Fkx2MwoPMQ/Wk0NTrS7jvm9it1RwyQqsPmx+Ok5ZLQUHgCBF4MgBOCmBQjQfcn4fcm4zoG21zC/V7Nrsl3xMmx6FDY1Ah8RjJG+NTySYdF5Dbunjt1Vy5xuZUeHYrOi8deGkuFiyD/gC9wZhA6HQyC4XsuIy+Xi8WAOEfAdKgG6KwG/KwFRDPFHE7u7jnnsBF1sYEeGYOPC8NEh2JAgTEBwXcouoVmUq2WPNLKHGthDDUywCJsaiT3enxgfhvltazDwVe75RtfW1i5cuDA3N5fP569fv37x4sVXHPDSSy99/vnnLS0tEonkwQcfXLt2rVveFwAPQeJoXBg2Lox4JRPpHOhwI3OwgX30BFOgZzMCsGFB2NAgbFAgliDHPHmHoWoze1bDnm5hTzSzp1rYSAk2OhSbG4u9O5IXJua6cAD0GPcE4erVq/v163fo0KEzZ85MnDhxwoQJERERlx8QFBS0Z8+e5OTkwsLCsWPHpqWlzZ8/3y1vDYCnUQnQ7Bh8dgxCCFko1JYrP1ewz5xmtHY2XY2lq7H+aixFiaUqMQ7nY2gdqMjA5uvYizo2V8vmaFkSQ4MCsSFB2Jp0YlgwpoaRL8A/YCzb3cUV9Xp9UFBQSUlJnz59EEJZWVljx4598sknOzp+3rx5aWlp//jHP675/0ZGRmZnZ1+Rox2x2Wx8Pp8gvLPtyfuZTCaZTMZ1KbyJzoFytWyulr2oYwv0bIGetdMoQY7FyrBYGYqWYtFSFC7GwsUoSHjjNtXO1D/FoGY722RDdRZUY2GrzWylGVWY2BIDS7MoWYGlqbBUJdZfjfVXY/DYd1Pg+88tN9a/G54Iq6qq+Hx+WwoihPr27VtWVtbRwVqt9ujRo/fff39HB7AsazAYxOL2M1IoFIpE/jEaD/gBlaCtBfXP5lG9E5Ua2QoTW2lCpUZ2fz1qsDL1VtRsYwUEChRiagFS8pGCj4lIJCGRiETC/w9Ip5Pk8+m2n10MMruQnUY2CumcrMGJ9E7UamdNLhQkRCEiLFKCIiVYpAS7NRrFyvAEOQarxwHQxg1BqNfrJRLJpV9lMllFRcU1j3S5XEuXLp02bdrkyZM7ejWNRjNixAgcbx91N2XKlM8++6yjg+GJkFsWiwXDPLjLyxuQCKUIUYoQoaAr/y+jC9M4kN6J9E7M6EI2GtlozEYhB9Ne5yTrFKD//5lg+wgRH2fFJFLwkIKPFDxWLUBqfgdNPhQym3vqH+Un4PvPrU7Wv1AoJMkbJJ0bgjAwMNBgMFz6VafThYSEXH0YTdNLlixBCH388cfXf7XON40SBAFByCGWZaVSKdel8FlShMKve4DJ5JLJ4LGOM/D955Yb698N05369OlDkuTFixfbfj137lxqauoVxzAMs3z5cq1W+8svv/D5/O6/KQAAAOAWbghCiUSyaNGi5557rqGh4fvvvz937tzChQsRQrm5uVlZWW3HrFq1at++fQ899NDRo0f37t1bUFDQ/fcFAAAAus890yfWrVv36KOPjhw5MiwsbPPmzWq1GiHEsixFUW0HmM3mvn37vvPOO22/zpw58+qnRgAAAKD3uWH6hHvB9AkvAsPHuQX1zy2of265sf59a0lEAAAA4CbBooE+S2fX57UU5LUUVhlrGsxNGpvOxVB2yi7nyxRCebA4MF4Vm6iKHRw6IECk4rqwAPgXs9Nytim3SFNSpq+sNzcZ7Eaj00jiPAHBVwjkYdKQKHl438Dk9KC+oZJgrgvr+yAIfU2rVbOn8tCRmhPVxtp+QSlpgSmzEqeHSUICxGo+zhOSQoPDaHCYGi1NZbrKP2pPbjjzSbA4cGzUyOlxE0MkV81lAwC4j9Fh2lN56EDVkXJ9Vf/gvqkByVkJ0yJl4QqBTCGQuxjKSTt1dkOjuanSWHO0Jvv9s58rBPKxUcMnxYyNUURxXXyfBX2EvuNMY85vxdtymi6OjR4xPnrUwJD+JH7jmmFYpkBTvKfi0P6qI8kBCQtTbx8U2r+T7wh9JNyC+ufWTdV/qa7i+/xfT9SfHhE+ZGrc+IHB6TzixlvxMCxbqCk+VH1sT+XBKFnE7KQZ46NH4Rh0aSHk1u8/BKEvOFZ38uu8n6wu2/yUWZP6jBORXZlk7aJde6sO/5D/q5AUrshYkhk24IZ/BS7E3IL651Yn679EV/7x+a/L9JXzU2bdmjBNwuvKiq4UQx+tPbGxaKvWplucNm9q3AQC8/frHgRhOwjC/Nai989+ZqPsy9LvHB05HO/2gk8sYg9XH/8k5+tgSdADg+6NU8Zc52C4EHML6p9bN6z/Fmvrh+e+OtuUsyz9zpnxU3i4G7qicpovfnnhhxarZuXAZaMjh3X/Bb0XBGE7fw5CrU33/tnPc5rz7s1YPDV2Qvcj8HI0S28t3f157ndTYycs77+oo0dMuBBzC+qfW9epf5qlfync8u3FX+YkzVzY9/auNdJcx8n6sx+c+1wpVDw65P5oeaR7X9xbQBC28wdGV+8AACAASURBVM8gZFh2U/H2Ly/8kJU4dUnaHUKyp3aN0zsM/zv7xdmmC08Nf2hwaMbVB8CFmFtQ/9zqqP7L9JX/OfaWSqhcPWRlhCysh96dYZnfird/deHHrIQpy9Lv5BN+t3QlBGE7PwzCamPd6yfewTH88WEPRss7VUvddKYx57UTG4aHD141aLnwrze2cCHmFtQ/t66uf4Zlvr34y8aiLfcPXD49bmIvlEFr020480mprvzJ4Q+nB/nXcl0QhO3cFYS0g3EaXC4LTTsY2kEjhHAeTvBxgYonUPEw3CN2WmER+0vhlm/yfr67/52zE2e6ty30+iwu64bTH19sLXxh1ONJ6vhLfw4XYm5B/XPrivpvsrT869h/+QTvmeGPBIoDerMkR2pOvHX6w0kxY+7LWOw5j4YOvcuucTIutu2iSvBxnI/zpCRfSpISNzy9QBC260IQusyUtclha3HaWx22Fqdd43ToXCzN8hUkT0rifJwUEggh2skwTsaudbpMlDhcqEqSqlJl8ljONvButrb+5/hbLpp6duTqcGkoJ2XYV3Vkw+mPFqfdMTclC0MYggsx16D+uXV5/R+pObHu5HsLUm9bkHpbb96kXmJ0mNadfL/KWPv8yEcTVHG9X4A21ga7Nt+kKzKbqm2kABcG8gkBTgjaL6q0g3GZKZeJYmhWqOIJA/jCQL4oSCAK4ouCBAIlD91MzUEQtrtOEDIU6zS4HDqXXee0a1z2Vqet1WFvdiICiYIE4mCBKIgvDBQIA/lCFY8UdxilDMWaqqz6YrMm14gRWPiYgOBMJUb06hf9YPUfb5363/yU2Qv73s7tFKIGc9NLR18Lk4Q8OfwhCU8MF2JuQf1zq63+KYb+8PyXh6uPrR3zZGpAErdF2l1x8P2zny5Kmzc/ZRZ2U6nSTSzS5Bnrj2jsrc7ADLkySSqPkxCCDi9WtINxaJ02jdPe6rS1OG0tDluzg7IzokC+KIgvDGi/MgtUPL6SR/Cv/ToQhO0az2tMxXYMw1iGpe0M7WQoG01ZaZeJYlwMX85ra9tsv+8I5IsCBV1/JGeRvsRcd7DVoXMl3BHRO0+Hdsqx4czH55vyXhi1JiUgsRfe8YZctGvDmU/ONuW+PPaZQFwFF2IOQRByy2QyuUjqxSOvSfniZ0c+Kud7xGfRYG56+Y91Ep742ZGrVUJlL7yjrdlR+nM97WQiJwYGpMu73JdEOxhbi8Pe6rS1trfVOfQup96FEOLJSJ6EIMUEISAIAR4+JkASIYQgbKcp0zsaKRzHEYaRIhzn4aSY4IkJnowkRT01gkaTayzf1BDQXx47K7RHuw8r9FVrj76eHJD46JD73T78upvabjxX9ls2I3ky12XxXxCE3DpZdfb1s+/OSpy2pN8dvfr4dSM0S3+R+8P28r3Pjlh9zfHeblR/WFOztyV6alDYqIAeqgPawTiNFGWlKRtNO2jawaiSpXwFD4KwHVejRik7XfRVDcKwlKVR13n8744tpbs+Of/NqkHLp/XK2LMuKNVVPHfwlYl9xqwYsATWfOIEBCGH2s7QZ0euHhY+mOuyXNu5pguvHFs/LW7CPf0X9cQyNCzDlv/WYKyw9r03RqC68XJxbgdB2I7D6RMsw5b92mCqsqb/PfY6XYxdYHFZ12W/V2WsfWn0k70zQaLL6jUNb55/n8TJF0c9LuVLuC6O34Eg5ATF0O+c/iinOe+5zEeTQhO4Ls716B2GV46tt1OOF0etCRIHuvGVWZot+KKapdmUZdE99DBwQ7AfIfcwHEuYF65MkuZ/Vs1QbruZKNaWrdjxqJQv/WDaGx6eggghGV/6xsS10fKI+3c9Xm2s47o4APQ4vcPw2L7nNTbtB9PeDJOEcF2cG1AKFK9PeGl4+OC/7VxzvO60216XRSU/1SGE+t4Xw1UKupcv/Bs4FJsVKlTzCr+sZhk3ZOGvRduePLB2RcaSx4b+3XMmA10fgREPDr7vrr5zH97zzMn6s1wXB4AeVKavvH/n4/2D014e+6yYJ+K6OJ2CIeyutHn/HPP0+lMffHD2c4qhu/+aVTuarE2O5CVRHjLHuvsgCLsHQwl3RNB2pmZvS3dexuQ0P3/4Pzsr9r0/7Y0JMaPdVbpeMyN+8r/GPvvaiXd+KtjEdVkA6BFHao6v2ffC3wYsvS9jMSczBbsjPSj1kxlvVRvrHtrzVIO5qTsv1ZpjbDlvSLsvpqNZDd7Id/4lXMFJLGVZVOMxranK1rVXuNBScO/21WGS4Pemvs7VZPnu6xeU8sH0N/dUHvrP8bdctIvr4gDgNixiv7zw44Yzn7w+Ye3EmDFcF6eL5ALZv8c/Nylm3N93PX6g6mjXXsRppMp+rU9eHMWT+tSm7hCEbsCTkvFzw4u+qaEdzE39RYZlvrzww4tHXn10yP0PDL7XLbu0cChYHLhhyqsO2vnI3uc0Nh3XxQHADeyUfe2R1082nPnftDcvX1/QG2EIm5dy6+sT1n6a++1rJ96xU/ab+/ssKvmxLnyUWhbtHc3CnQdB6B4B6XJ5rLhy2020OTRamh/e82xO88VPZrw1IiKz58rWm4Sk4KXRT4yIyLx/55pCTUnPvRFDsdYmh77Y3PaftdFxs3chwBsxFGtrdhhKLW2fu6XOTjt78HNvtDQ/sPspMU/01qRX1CJVz71Rb0pSx388Yz3Lsit2PFqkKe38X2w4rqWsdOTkoJ4rG1e8+xHEo8TdFnbm1ZLQESpJ2I0nv+8q3//+2c8Xpt1+R8ocr+tvuD4MYUv63RGn7PP0wX/+faCb50Hamh2tOUZNrsHa5BCo+O2LE7LIaXTZdS5SRMhiRPI+YmWytDOfAvAKdq1TX2g2VlhN1VaHziVQ8fiK9qXwXSbK1uIQqPgB6bLA/gqpW59UzjVd+Ocfby5Omzc3+VY3vqwnEJHCp0c8cqDq6FMH/zEvedaitLk3nApMWenqnc3pq2J9ZoDM5WAeoTvVH9Fo8039Vva5zjF6h2Fd9vs1pvrnRz6WoIrtraL1iOvP46ky1Dx3+N9DwwY9MPie7s/ntWudVTuaDcXmwIGKtufvq09Iu9ZpqrQZKy26AjNLswHp8sCBCnmM2JMW/XAn355HaG10tJzTa3KNLiutSpEq4iWyGLE4WHDlp8kic71dk2toOWsQBvD7ZIVII90Qhz8X/v5d/sYXRz0+MCS9o2N8oP5brK3/Of62nXI8O3J1pCz8OkeWb2pgKDZh3vWO6WUwob6dpwUhS7Nn3yiNmx2qSr32x3Oo+tjbpz+cFjvxnv6LeAQHazG41w2/iGan5V/H1lldtn+Meao7yx7WHWqt2dsSPiYgYlxgJ+ct2ZodrecNLecNjJMNHqIMzlQKA7xjRkrn+cCF+GouM9V8Rt98Sk/Z6MABiqAMhTRK1JlbGZZmm7J11bubAzMUsbNCu7wyvp1yvJH9brWx9l9jnwmRBF/nSN+ofxaxvxVt/+LC98vSF9yWlHXNBipbizP3nbJBTyV61BgZCMJ2nhaECCHtRVPltsaBjydc8byis+vfOvVhhaH6qeEPpQWmcFU89+rMF5Fh2S8v/LCtbM/a0U/2C7rpfzjtZEp/rLNpnKl3RwuUXbl1sNTZm07pWs4aJOHC0BGqgH7yXt48pOf4xoW4HYv0JebG41p9sUXdTxYyVKWIk3ThUZ62M8Xf1bqsdMqyKL7spq/adaaGFw7/J1Ed99jQVYIbzeX1pfqvNdW/dmIDyzJPDn/46qU8Cj6vlsWIIid6Vu8gBGE7DwxChFDuhvLwMQGBAxRtv7KI3V629+PzX90SP2VZ+p3eMlO+Mzr/RTxed/q1E+/clTZvXsqtnV+emHYwef+rEIcI4+eF42S30ouhWG2esfG41trkCBmqCh2h5mR1RPfyjQsxZaGbTukajmlJAR46Qh00SEkIuzeIj0U1e1uasnXpD8be1M3TkZrj606+v7z/otmJMzpzvG/U/yUMy24u2f557vdzk29dlDb30iB2S7394sdVmc8ldfMcdDsIwnaeGYTafFP1zuYBj8UjhCoM1f89+QHFuNYMfcDbewSvdlNfxEZL84tHXg2VBD81/GEJ78abWDEUm/9JlVDNS5gf4cZOPluzo+GYtvmMXt5HHDYqQJUs9d4eRG+/EJuqbQ1/aLR5JnU/WdjIAFmMO4e61B1qbTqhS38wjteJnddcDPXRuS8P1xz/x5inOr/fmbfX/zU1WVrePv1hranhsSH3DwhJRwgVf1srDhNGTnTnUqVuAUHYzjODELHo7Bsl4VnqjbZNeyoPLu+/6NaE6T42NLTNzX4RXbTr/XOfn6g7/eLox2+whSmLCr+qRghLXhLZE6PUGCfTcs7Q8IeWstGhI9UhQ1WduVx6Gi+9ENNOpuWsofGYlrLTYSPUIUNVXd8l9LqqtjfpiszpD8Refw2UBnPTP46+ESBSPTXi4ZvaU9BL678zjtScePfMJ2mBySvillZ9oB3yXHJ3H9N7AARhO88MQoZl9mz7o/WsoWFa1YoBSxQCOdcl6ild+yIeqTn+31P/uyNl9oLU2zq6P6g/rGk9b0h/ILan+/NM1bbGY1rNBaMqVRo6Qt21fimueN2F2FJvbzyubTlnUMRJQkeqe+FxvPi7WpyHJ8zvcKzjvqojG05/vDht/tyUrJvdU9Dr6v+mOGjntxd/Me1yxgfGTF4yQuhhW6IiCMJLPDAIj9Zmf3L+60BBwNw/7khfES+N8Lhvjxt1+YvYbG195Y//4jj+7IjVV+8OY6m35/2vMmN1nFDdS/2plI1uPq1vPK5laRQyTBU8RNmFcRa9z1suxJSdbj1naMrWOY1UyDBVyDBV18Y9dQHtYM6tK429NTQg/cr7UYvL+tapD4u0pS+MWpOoiuvCi3tL/XcZZaFP/rvo6MQDJ41n7k5fMCN+ck/sa9hlEITtPCoIzzTmfJrzrYN2rBiwZHh4Zu2BVluzI3GBp2+l1B3d+SIyLPtd/sZfCjevGnTv1Njxf/45xZ7/b1nkpMDgwV2fbtFlpipr4wmdJtcojxUHD1Gq0+SeNkDgch5+IWYZ1lBiaTqt1+WblEmSkGG98Qh4NVOVreCzqgFrEvjyP29uzjVdePX428PCB68adI+QFHTxlT27/ruv7mCrtcGRuDCiUFPy0fmvmq2ty/svmhA92kM6eiAI23lIEJ5uOP9l3o8Gh3FZvwUTYsa0fUtcZurMf0qGvJjsG/t1XVP3v4gluvJXjq2PkoU/NvTvbRMNq3Y225odKUuj3FTGrqCdjCbX2Hxab661BaTLgwYqFAkSD1xQw0MvxCwyVdtazulbzxsESl7QYGXwIGUP9QJ2UvWuZmuDPeXuaISQnXJ8nPPVoerjTw57cGj4oO68rIfWv/ucfa0kYUGEvE/70LYzjTmf535ncpqX9lswIWb0DRej6WkQhO24DUKapQ9XH/8uf6OLoRanzZsYM/aKG6XCL6qVKbLQ4T6yROHV3PJFdNGuLy58v61s76pB94xTjTq/rnTAmoReazq7PqeRarug2zWugHR5YIZcES/xnGmInnUhZpGp2tqaa9TkGHEeFjhAETRIIQrq4sOWezEUe/bVksSFEeXS8tdPbOgXlPpQ5n03NS7mmjyr/t3NWGEt/alu0FNXjqE93XD+q7wfW6yaBalzpsdN6vLzdPdBELbjKgiNDtO2sj2/FW8Lk4Qs6HvbiIjMa3az6wpM1bubMx7x7hXrr8ONX8RibdlrJ96ZUjgtPSElbVZXOmx6lF3r1OQYW3ONtlaHKkUWkCZTJktJEcdNEZ5wIWacjKHMoskzafNNpAgP7K8I6C+XhHtc13jNqcaLOys+T//00aH3Dw93zxr3nlD/Pafk+zpxuCBi3LVnTVxsLfw+/7cLzfkz4yfPSZpx/SV4eggEYbteDkIWsbnN+VtKdh2vPzUmasTc5Kzr97GzDHv6X8VpK2LEProAtHsvBPpKS86nJe+kv31b35kLUud45hJ0TiOlzTdpLxoNZRZJmFCVLFUmSaXRIk4aTjm7ELPI0mjXF5v1RWZjhVUaJVKnyQLS5MJAT1wsgkXszvL9H5/7emXRyrRJ8ZHD3HbJ9uEgpO3MqZeLBj9zgzXVGsxNvxZv21m+Ly0w+daE6cMjBvfmaBoIwna9FoS1pvq9lYd2lR8QEPxbEqZOi5vQyXaVqp3NtJ2OmxPW0yXkhHsvBBferwgeomT7ut4980mFvvqhzPvcdefeExiKNZZZdEVmfbHZoXPJY8WKeIk8ViyNEvVa22lvXohZhrU2OozlVkO5xVBmIQS4MkmqSpIqk6QeOMPskiJN6VunP0QIrR6yMtwUUfhldeZzSe76gHw4CBuPa/XF5pRl0Z052EE791cd2Va6u9bUMKXPuMmx45LVCT1dQgRBeElPB2G1se5o7YkDVUc1Nu346NHT4ibc7Ads1zpz3i4f+lKyBw616D43fhFNVbair2sGP5vYVlEnG85uOP1xqCRk1aDlscoYt7xFz3GZKUO51VhmMZRb7K1OSbhQFiOWRYukUSKhmt9z4yR7+kLsNLjMtXZTtdVUbTNX2Xhysi3vFQkSD+nEvY4Wa+vHOd+cbji3YsDS6XET2zovLn5UGThAETLUPd32PhyEuRvKoyYHdbR5QEdqTfW7Kw7srTyMIWx8zKixkSOSAuJvdnZm50EQtuuJIHTSzgstBdn1Z47XnbZR9tGRw8ZFj8wITuvyEKmct8pibglRJkrdWEgP4cYvYsEX1coESdjogEt/QjH05pLtX+f9NCpy2PL0hYHigOv8dc9BOxhTtc3cFh41NtrBSMKFknChOEwoDhGIQwVu7Fl074WYdjC2Zoe1yWFtsJvr7ZZ6O2KRNFIojRbLokWyGLG3LL5jcVm/z/91c8mOWYnT7+o7T8z7c+U2fYml/Nf6QU8muuXi7KtB6DS4zr1ZOnRtSpcfnYu0pYeqjx2pOWGn7MPCBw8LHzw4NOPyD8ItIAjbuSsIHbSzQFOc23zxfFNegaY4ThnT9uElqd1wO1N3oNWmcXrUPl7u4q4voq3FkbuhYsjzSfhVS2GZnZbv8jduKd11S/yUO/vephQouv92vclloS31Nku93drosDbabU1OjECiIIEwkC8M4AsD+AIVT6Di8eW8LkxY7Fr9swzrMlF2rcuhczl0TrvGaWt12lqclJUWhwhEIQJJqKAtuT3/se8KdsqxqWT7D/m/jogYsrz/ouCr1mpACJ1fXxY9NVid5obvra8GYf1hjaXB7pY50DXGuuP1p7Prz+S3FsUqYgaGpqcHpaYH9e3MasM3BEHYrstBaKccFYaqUl1FsbasUFNSbayNU/ZJD0odEJKeEZzmlg/pz/fy3dZRd30RS3+u58vJ6GkdjmJotWm/yft5X9XhWxOmzU+ZrRJ6WRxezmmibM0Ou8Zp1zgdOpdd63LonE4jRYoIvozktf0nJkgJQYoIUkQQApwQ4ISQwHDU9jSJEVjb5FSz2SyVShFCjItlXAxCiLYzLMtSVppxsbSdpuwMZaNdZoqy0i4L7TK6nCbKZaF5EkKg5guUPKGaLwzgCQP4oiCBQMnzouXlrmCn7L+X7Py+4LeM4LTl6QtjFB3OQ209b6g/oun/kBtGJvtqEOZuKI+aEqxKcWcjlpN2XmwtOtd0Ibf5YpG2NFgcmBKQmKSOj1fFJihjpXxJF17TjfXvBetIdZPBYWyxapqtLQ3mplpTQ52podpYq7XpYhRRCarYRFX8jLjJCao+Pbc7klDNF6p4hjKrMrErH7bPc5mp1hzD4GeutwZ3oEi9esjKRWlzv734y9Itq6bEjl+QOidE4lm7o3USX0byZaQi/sovg8tEOc2Uy0Q5TRRloV1W2tbipO007WBoB0PZacQgykYjhFiapR0MQohlWQzDEEI4D8N5OEKIEOIYhpFiAudhhIAghTgpJgQqnjRCxJMSPDmPLyN5UsKX7smMTtNvRdt/K946ICT9vxP/ecMe5YD+8sptTaYqm3s3u/AZToPL1uxw+8WKT/AHhqQPDElHCNEsXWmoKdSUFGlK91cdKddXiUlRjCIqQhYWIQsLl4YGiwODxYFKobLXlrDx7ifC07XnC/UlOI5bXTaaoe203U45rC6byWk2OEwGh1Fv14t4oiBRQIgkKEQSHCELi5SFxygiQyXBvbksQu2BVrsvto665Y6s7lCrtcGReGdn22G0Nt1PhZu3le0ZEjZwQeqc3hmf5pl89Ymkk+rNjT8Xbt5bcXh01LCFfedevZ1sR2oPtNpbHAl3dLfpzyfr343top3XbG2tNtTWmurrTA0NlqZmS2uztdXkNCkEcoVALhfI5XypmCcSEAIJT4xh2Mz4yZGycHgibOegHGanBcMwMU9E8IhAIkBICiQ8sZQnUQjlSoFcKVRe2l6SQ4EZ8py3yuNvD/OlO3F3acrWJcy/ibNOLVLdP/DuJf3u2Fq6+4XDrwaK1HOTs8ZGj/SEDxr0AoZlTzec+7V4a0FrSVbClC+y3g0Q3dwo0JBM5ZlXS2Jnh/nw8odd1ppjiJrS27Pj2x4BM8MGXP6HNEvr7Ua9w2B0mAwOo42yOyiHlbL1RAG8+9oxKKT/8KhMztcavSGhmi9Q8YwV1qsbxPycqcrKMujSYoadJ+GJF6TOmZ8y64/ak78Vb9tw5pPpcROzEqZGynztsRtcorHpdpTv21a6W8qX3JZ0yz/GPC3oUo8GT0YqEiSt5w0hw3x2+cOucZooa5P720W7hsCIAJHqZu9yusa7g9CLqPvKdAUmCMIrNGXrQoapujxGA8fwMVHDx0QNrzXVby3d/dCeZyKkodPjJo2PHtW17nfggZy081jdqZ3l+/JaCsdFj1w7+snkgO62h4cMU9Xua4EgvIIu36RMknrOarq9BoKwl6hSZCU/1vXJ4rocnoR2Mq25xkFPXrmqbxdEysLvH3j3igFLsuvP7Czf//7ZzwaFZkzuM3Z4eCaHiwKD7qAY+mxTzv7KI0drs5PU8dPjJq4d/ZS7Pk1VirT0pzprk0McAl+PP+kKTeq+PruR+HVAEPYSWbTIZaYcOpdA5WVzs3pOa45BESe5fJe4biIwYmTE0JERQy0u6+HqY9vK9rx+YsOQsIFjo0YMj8h076wY0EOctPN04/kjNSf+qD0ZJQ8fHz16xYClbm8fw3AseIiq+aSuz62h7n1l78XSrL7YEn+7P3YuQBD2FgypUqS6AlPoSDXXRfEULWcNPbRHlYQnnhE/eUb8ZKPD9Edt9t7Kw+tOvp8SkDg8InN4eGbnhxeCXtNibW1bzulc04UkdfyYqOEdzYh3l+BBiosfV/XJCvXe2ZPuZaywCoP4PJk/hoI//pu5okqRtZzTQxC2oWy0ucqmWt6pVX27TC6QtSWinbKfacw9Xnfql8LfMYQNCRs4ODRjUGh/hcAfG4I8hI2yn2/KO9uYc6rxvNamGxI2cELM6KdGPNz9nQI7QxwmJAS4qRomFLbTFZrUN7m4qM+AIOw9qlRp2S/1DMV2YTEt36PJNSqTJcRVa6r1ECEpHBU5dFTkUIRQlaHmVMP53RUH38h+N0QSNCCkX/+gtPSgVG9ZztSrGR2mvNaC3Ob8nOaLFfqq1MCkwSEZTw9/OEmd0Guzpy8JyFC05hogCNto8829PH3Qc0AQ9h5SRIjDBMYyizLZBxfgvlmtucaQIUpO3jpGERWjiJqXcivDMsXaspzmi3sqD7516kMBKegXmJwamJwakJSgioVRNm5BMXS5vrJQU5LfWpSvKW61alIDk/oH9b1/4N2pAYk9t6JTZwT2lxd8Vh0LraMIOfQul5mSRfvpPQEEYa9Spcq0BSYIQspGmyqsKUs7XBCyd+AYnhKQmBKQuCB1DkKoxlhXoCnOby3eV3m4wlAdJg1JVMUmquLiVbFxyj5evcBpb7K4rOX6yjJdZYmuvERbXmWsjZCGJgckpAWlzEuZFaeM6c1Fna5PEi7ECMxca5NG+WkAXKIrMKmSpX57QwBB2KuUidLSn+u4LgX3NHlGRaLE09b1iJJHRMkjpsZOQAhRDF1hqCrRlpfoyv+oO1Wuq8QxPE4VEy2PjJFHRcsjIuXhweKg3m/N8zStNm2tqb7GWFdtrKvUV1caqs0uSx9FdLyyT5I6fmb85HilRz9bB2bIW3OMEIT6EovfdhAiCMJeJo0SOnQul5niSf265jW5xqABHv2AReJEoiouUfXnHgUam67SUF1pqKk0VB+tPVFjrNM7jGHSkAhpaLgsNFQSEioJDpEEBUsCvW6vqM6wuKzNlpYGS3OTpbnR3Fxvbqw3N9WZ6oWkIFIWEaOIjJJHDArpH6uMDpEE9dxerG4XkKEo/KK6T1YI1wXhFIsMpZZYP55J4teX496H4Zg8VmwoswRm+OC1spMYF2MosyTdFcl1QW5O22pPg0MzLv2Jg3bWmxrqzI315sYGc9O5ptwmS2uztcVOOYLFgQEiVZA4UC1UqkUqtVCpFCqUAoVSqFAIZEJSyOE/5JqctNPoMOkdRp1dr3cYdHaD1qbT2HQam7bVpm22tOAYHiwODJEEh0iCwqQhqYFJ4dLQCFmYt8/OlEYIEYv8fGa9tdFOCHB/nuIMQdjbFIkSQ6lfB6Gh1CKNFJFCT18h9oYEBD9WGXP1vj8O2tlsadHYdS3WVp3doLFqKwzVOrteZ9MbHEaj08SwrJwvlfGlUr5UwhNLeCIxTyzlS0SkkE/wpXwJH+cJSAGf4AsIPoHhov/f2lvG/0vvssVqMWGWy//E6rLSLIMQslMOiqEohrJRdhfjslMOm8vmoJ1Wl83isv7//1pNTrPRaTY5TTTDKAQyhUCuFCrUQpVKqFCLVHHKGLVQFSgOCBYHun17cc+hTJHqCkz+HIT6UouHrC/KFQjC3qZMkBYer+G6FFzSFphVKb7cGyEg+G3djR0d4KCdJofJ5DSbnGaLy2Z1Wa2UzeQ02ym73m6oMzU4aKeTybXENwAAIABJREFUdjpoh5N20Sxjc9kQQixizc6/xB7DMDj+l35WMU9MYDhCSEAKeDhJ4ISYFPFwnpAUiEghn+TLBNJQabCIFEl4YglfLOVJZHypXCATed5Daq9Rp8rqD2sixvfg5H0P5+e35giCsPdJwoWUhXIaXHyFnzZE6ApNqT08j97DCQi+QBzQ/WmLPrkfXu9TJEqKvq2hHYynjd7qJSwyllvifW631Jvilx88tzAkj5cYyiw3PtIX2ZodLMVKQv33+QN4GoKPy6LF+hIz1wXhhrnOxpOSfL9cWe0SCEIOKBMkhlI/DUJdoVmVKvOeQYXAL6hSpboCPw1CQ4lFkejvM5shCDmgSJTqS/w0CLUFJlWKv591wNOoU2W6QhPXpeCGocyiSPDrkTIIgpAT4mAB7WAcehfXBelttJMxVVmVSRCEwLOIggUYjlkbHVwXpNexyFhpVcR59xyY7oMg5AKGZDEiU5WV63L0NmOZRRop8tMhCcCzqZKluiK/eyi0Ntp5EsLP1/dAEIRckcdKjBV+F4SGUovS7xthgGdqm+DLdSl6m7HCKo+FUxKCkCPyWLEfBqG+FLrlgYdSxEuM5VaWYbkuSK8yVlhlffy9XRRBEHJFGiWyNTtoB8N1QXoPZadtzQ5Y3Rh4Jp6U5Ct4ljo71wXpVcZKqzwWghCCkCM4iUnChOYaG9cF6T3GMqssRgybEgOPpUyU+NVwbpeJoqy0ONh/15a7xG1BaLVaq6urGabDRxyapquqqux2/7rhug6Zn7WO6kvMfr6eIfBwigSJodSPZhMaKqzyWDFM6kXuCsJ33nknIiJi4sSJKSkpBQUFVx9w+vTpuLi4KVOmhIeHf/XVV255U28n7yM2VvjR7aehFKYrAY+miJcYK60s7S/dhKYKC7SLtnFDEFZVVT377LPHjh0rLS1dtGjRI488cvUxK1euXLNmTXFx8e7du1etWqXRaLr/vt5OHicxVdr8pHOestIOrUsaCR2EwHORYkIYwPefDgsYKXOJG4Lwhx9+mDBhQmpqKkJo1apV+/fvb2pquvyAgoKC/Pz8++67DyGUmZmZkZHx22+/df99vR1PQvBkhLXJL+bwGkotslgxRkArDPBoygSJ3j8mUTAuxtrokMHgNYSQW3afqKioSEpKavs5ODhYJpNVVVWFhIRcfkBERIRY3H7rkZiYWFlZ2dGr0TSdm5vb2NjY9mtYWFh4eIfLorMOq6u5msa9dciPJIjSnaviUd5aftpqdYo7dUepPU9JAzFnTUlPF8mvdL7+QSdJlEzjOTokWd+Zg726/k31rFDNUk1lXBfkpmEYTobFYIQ7FwFww2uZTKbAwD+38pJIJAaD4YoDRCLRdQ64nMFgeOKJJ3i89i2KBg8evH79+o4ONv+xA1081vWic87aT1sRSOQd5LocXcQwjK1zdyF6zR1BssOaksaeLpJf6Xz9g05iWIG5eUnr9+9h6MZ9Fl5d/3pLf4JWaL4/wnVBukI49S4yvr/Z3KmRTUKh8FKgdMQNQRgcHKzT6S79qtPpLn8cbDtAr9dffkDfvn07ejW1Wr1r166IiA43Nb0cOXEuf/pCgvDWvc5l1baSn+vC1tzBdUG6qJP74dFOpvLFwuinX4C5E+4F+xH2hKbXShR3rZNE3HinMK+uf+M3NaEpsuDMhVwXpFvcVf9uuJ3JyMg4efJk2895eXk4jsfHx19+QN++fZubm+vq6tp+PXnyZEZGRvff1wdIIoT2Fift9PFp9aYqqyRCCCkIvIKsj1/MazJV2WTR0EHYzg1BOH/+/Orq6jfffDMnJ2f16tXLly+XSCQIoaeffvrVV19FCIWEhMybN+/BBx/Mzc196aWXWJadMWNG99/XB2AEJg4V+PxiFqYKWL0CeA15rNhY6eNBSFlpykqLgmAqfTs3BKFEItm7d++xY8dWrlyZmZn5+uuvt/15ZGRkWFhY288ffPBBdHT0ihUriouLd+3aRZL+vtj5JdJosanax886YyWM0gZeQ95HbPL1IDRV26SRIphKf4l7Aik9Pf3XX3+94g8ffPDBSz/L5fK3337bLe/lY2TRIh/fEZRFpmpb0iIIQuAdREEC2sk49C6B8gYjLLyXudoqhXbRy3jrkCefIY0Smap9eQKvpdHOk5Cw4RnwGn6wXaipxgYzCC8HQcgxcbCAstCUhea6ID3FVGmVQ7so8Co+v0uaucYGT4SXgyDkGoakkUKT767qZKywymCkDPAqvt1N6NC5EIt8uOG3CyAIuSeNEvvw8oZGeCIE3kYaJbI2OhgfnddkgsfBq0AQck8aLfLVgaOUhabMtDgERmkDb4LzcHGowOyj85rM0EF4FQhC7kkjhb46ldBUY5NGwSht4H2k0SJfbacxt52V4DIQhNwTqvi0k3GZKK4L4n4mGKUNvJMsSuyrw7ktdTZJBJyVfwFB6AEwJI0Q+mQ7jLkalnECXslXOywcOhfCMb4cpjP9BQShR5BEiCx1Pnj7aa6FRhjglXx1XpO51gb7Y18NgtAjSCOE5lpfeyKEUdrAi2FIEiE01/ra7amlzi6NvPHGGv4GgtAjSCJFZp97IjRV22QxMHECeCtZtA92E5rrbFLoILwKBKFHEAcLXCaKsvtUO4y5BkbKAC/mk92E5lq7BJ4IrwJB6BkwJAnztUkUJhgpA7yZLFpk9q0nQpeZYpyMUMXnuiAeB4LQU0h8bDYhiyy1dmiEAd5LoOQhDDn0Lq4L4jbmWrs0Ugjzeq8GQegppBEiX+qZtzY7eDKClBBcFwSArpNG+dRDoaXOJoEho9cCQegpJJEiXxo4aq6BUdrA60mjRL60IL65zi4Nhw7Ca4Ag9BSSUIFd62RcPrLOr7kW7j2B15NG+tQEX3gi7AgEoafACEwUyLc2OrguiHvAdCXgA6QRQp9ZcZR2Mk4DJQqCkTLXAEHoQSThQku9T7SOsshSByNlgNfjK3gYjvnGeBlrvV0UKsBwGCpzDRCEHsRngtDW6iDFBCmGkTLA60kiRb4xnNtSDx2EHYIg9CCScKGlwSdOuTqYtAt8hNRXFlqz1NvFYXBWXhsEoQfxmSdCM8wgBL5C4ivzmiz1dgk8EXYAgtCD8KQkTvhCh4S51gYjZYBv8JF9s1lkaYQg7BAEoWeRRPjCQyHs/Al8hlDtC/tm27VOUkSQIui2vzYIQs8iCfP6IISdP4FP8Yl9sy31dgl0EHYMgtCz+EA3Iez8CXyMD+ybDR2E1wdB6FnE3h+EMJUe+Bgf2DcbgvD6IAg9izhY4NC7vHqhNTN0EALf4gP7ZlsaIAivB4LQs2AEJgry7oXWLPV2aQSccsB3iIL4LiNF27319pR2MC4jJQyExdU6BEHocbx6vAxloykrLVTDKQd8B4ZjolCBtdFbz0prg10UAourXQ8EoceRhAm995RrH5wGZxzwLV59ewrtojcEQehxxGFCS4O3No1CnzzwSf/X3p0HNHGmDQB/cyeTcB9JuBVBVq1yKB5YD6wWarFKPSrFtatbrdpLV3frtp/aY2vV1qNaaz3Q7Umrrd2q2NauZ+uNVUEUKKIECEdIArnP+f4YN6WgEMgxk+T5/TUZ3sw8CZk8mXmfeV9+BFftsYlQ22DARByyo6A0SISUg4k4Wo8dcVRTr+dDByHwOh59X5NGCjcRdgMSIeVwAllWC25Se+RIFpo6OCMEXogfwdVK9QgnO45e0TYYYLjtrkEipCJMxPXEwlHciuuaDJgIDjngbZg8BhNj6FuMZAfSY0aVGeE42w9GeuoKJEIq4os4njgfk67ZyA5gMjjwoQJeyEOvjmqlMPtS9+A7i4owsUeeEUKlDPBi/EieRybCBgN0EHYLEiEV8cUeWS+jqdfzI2BMGeCdPPSMUCPVQ8lotyARUhEm5mo8sGcezgiBFxN4ZiKES6P2gERIRUweg8FleNwMvZo6HSRC4K24IWyTxmzWWcgOpCdwpG00YEI4I+wGJEKK4os9rF7GrLFYTTgnkEV2IAC4Bg1hYq5n9Vno5UYmBvPxdg8SIUV53CGnIa7AwOBqwHvxxVyNR1Wxwa30doJESFF8kScecnAFBngzvpij9ahuQriV3k6QCCkK87TCUW0D3EoPvBwm4npWh4W2AUpG7QKJkKKwcI6u2YhbPaZyVFMPlTLAy/EjuFqpwYPKubVSAx9+ntoBEiFF0dl0tj/TY4Z0IorT4Lcn8GpMjEFn0wytnlHOjVtxXYuRFw6Tg3YPEiF1YSKOp4wvo5cbmTwoTgPe797o255ALzNyAph0FnzJdw/eI+ryoKG3oTgN+AjMc2bo1Uih295ekAipiy/ieMpU9VopFKcBn8D3nHmzoVLGfpAIqcvTzgjhkAPez4Nu8IVCbvtBIqQunpCja/GMwlEYzxD4CIw4Ki2ecFTCGaHdIBFSF51J4wSwdM1ULxy1mnG9woSFwyEHvB+dSeMEsnRNVL9Ug1twvdzEC4Oj0i6QCCnNIwpHdY0GbjCbxoDR1YBP8IhuQl2zkRPEojPhqLQLJEJKwzyhXkbTAB2EwIdgYg84KrUNej5cF7UbJEJK84h6GSgZBT6F7wkDrWlglNGegERIaR5xRgh98sCnYGIP6LCAo7JHIBFSGhbO0ctNFC9Rg7vpgU/hhrBNKrPFYCU7kK7AvRM9AomQ0mgMGjeY0iVqFoPVrLVwg2E8Q+AraHQaL5yjbaTuUWk14waFiRcKR6W9IBFSHSak9CGnlep5Qg7Mxwt8Ciai9CxpuiYDN4QFhdz2g0RIdTwRl8qJUNMA87wAn0PxKjZtowETwlHZA5AIqQ4TUrpnHvrkgQ/iiznaRuqeEcKcaD0FiZDqKF44Cn3ywAdhIkrfUw8/T3sKEiHVUbxwVAvDbQPfwwlkWY1Wi46ihaPaBrg02jOQCKmOxqBxglg6GRVHHDVrrVYzzg5gkR0IAO5FQ5iQo2+i4lGJW3CDwsQLg5LRHoBE6AEoe3VU32iE0SuAb8LEXF2Tiewo7kPbZOAGQ8loz0Ai9ACUrZcxNJtgPEPgmzARx9BMxUSog0qZnoNE6AEwIUXvoNA1maBSBvgmTMTVN1IxEUL9Wi9AIvQAlL40Cr89gU/iizk6SvYRahv1mBCOyp5xWiLU6XSXL1+uqal5UAOZTFZcXFxbW+usPfoOXjhH30LFwlF9swkSIfBNLAGTRqMZVWayA+lI22CARNhTzkmExcXF8fHxL7/88rBhw1auXNm5wVNPPZWYmLh48eKUlJRp06YZjVT8JUVZ9ybFpljhqLHVRGfQWAIm2YEAQA5uOItqA60RE9NzYWL6HnJOIly+fPlLL730888/X7ly5cMPPywrK+vQYN68eQ0NDRcuXLh9+3ZpaenevXudsl/fQcGp6jVSA1cIJdrAd3HCWFQ7KnVNBm4wTEzfY05IhI2NjadOnZo/fz5CKDIyMisr66uvvurQZtKkSWw2GyHk5+c3cODAxsZGx/frUyjYTaht1HNC4XQQ+C6ekEW1GXq1jXBdtDec8EUmkUj8/PxCQ0OJh3369Omip/C33347derUmjVrHtTAZDKdOXPGtrWIiIgBAwY4HqSnw4TcltI2sqP4A63UwBXBGSHwXdxwdmsJxY5KKBntFXsT4dq1a0tKSjqsHDp06LJly7RaLXG2R+ByuRqN5r4bkcvlTz755PLly5OTkx+0I41Gs3XrVg7n3o+aIUOGvPnmmw9qrNPp2Gw2g8Gw81V4MH+Lul6nVqvJjuN3qjpNcAKPUiH5Go1GQ6PBRTDSWAUmTb1OrVJTZxqytjp14EC+jxyVdn7+uVwuk9lNprM3EWZkZMTHx3dYGRkZiRASCoVKpdJqtdLpdIRQS0uLSCTqvIXW1tasrKxHHnnk1Vdf7WJHgYGBX331FbHlbjEYDB9JhFgcXqlo4PP4VBkwAkcGmTkw1l8gEJAdiu/CcRzefxLhOM7ktbHMHE4QVUYZNMqkwbEBmMAnTgqd+Pm3NxGOGTPmQX+Ki4sLDAw8f/78qFGjEEK//PLL0qVLO7TRaDRTpkxJTk5+9913ex2rL7MVjlKkA0AvNzIxBoMD96ECn8YXc7QNeookQigZ7TUnfJFxOJzFixc///zz33///cqVK5ubm6dPn44QOnv2bExMDNFm+vTpNTU1aWlpu3bt2rlz5+nTpx3fr6+hVOEodEUAgBDChBSajwlKRnvNOVV/q1atCgkJ2b59e2Rk5KlTp7hcLkJIKBTOmjWLaDB27NghQ4ZUV1cTD221MMB+mJCYC9Sf7EAQQkjbALMvAYAwEae16v4lEe4HJaO9RsNxao1XEhUVdeHCBTv7CH2oWAah5iutLaVtSX+OJjsQhBAq/6w2qL+A15/h5+dHdiy+S6VSwftPIpVKhRTMqgP1ycs61k+Qoub7JoRQTFY42YG4iRM//9DH4zGodWlUClNgA4AwEUfbZMCtlDid0DbCUdlLkAg9Bi+co28xUmHEUdyK65qNvHA45ICvY7DpbAHTIKfENBQwymivQSL0GHQmVaaq18uM7AAmgw0fHgAQJuZQYXwZomQUfp72DnyXeRKKzNCrkRr4MDE9AAghhDARlwpHpQ4mpncAJEJP8r/CUZJpG/QYJEIAEEII8UWUOCOEklFHQCL0JBT57amV6vnQJw8AQgghTMylwoD4cGuvIyARehKKzEGhgUMOgP+hyLzZUDLqCEiEnoQKh5zVjBsUJl4YzDsBAEII0Zk0bjBL20TypRqNFC6N9h4kQk9yr3C0mczCUW2DnhfGhj55AGwwMZfcqerv/TyFktHegkToYYhBfkkMALoiAOiA9MEudE0Gbgj8PO09SIQeBhNxNaQecjDKKAAd8MVccgtHYaQnB0Ei9DCYiEPuRRiNFM4IAfgD0su5NQ0GKOR2BCRCD4OJSK7V1kr1GJwRAtAOL5RtUpstBitZAcCtvQ6CROhheGFsQ6vZaiLnkDPrLWadhRsEJaMAtENDvDAyuwmh595BkAg9DI1O44WytY3kHHJaqQETcRB0yQPwRxh5VWwWo9WoMnNDWKTs3TtAIvQ8JHZIaBv0MMooAJ3xyeuz0DYYsHAOjQ6/T3sPEqHn4YtJq5eBShkA7gsTczX1pP08hZJRB0Ei9DyYiEvWpVFNvZ4fAYkQgI74EVxNvY6UXUMHoeMgEXoejLzR7uF2JQDui+3PRDSaUWV2/67h1l7HMckOwDmWL19eU1NDdhTuo6zUBFzD3NwrYDXhqhptQCXftuaJJ554+umn3RkDAJTFF3O09Xp2f4Gb9wsdFo7zkkT46aefvvXWWwEBAWQH4kOKiorOnj0LiRAAAibmaur1ge5NhGadxWKwcAKhZNQhXpIIEUI5OTlCoZDsKHxIfX19RUUF2VEAQBV8Mbe1SuPmnWrq9XwxF+5ochD0EQIAgBOQMuLovUQIHAOJEAAAnAATc3TNRjdPF6qV6jEo5HYYJEIAAHACOovOCXT3dKEaKZwROgEkQgAAcA6+mKupd+PVUZy4iRDunXAUJEIAAHAOvtit9/jqW4wsPpPJY7htj94KEiEAADgHJua6c/hDTT3MieYckAip5cKFC+PGjbM9fP3117/++usu2u/du3fTpk0uDwsAYAd+hFsLRzVSGPLQOSARUovRaJTJZMRybW3tv//975ycnC7az5gxY+PGjQqFwi3RAQC6wg1mm3UWs87int1poVLGSSARukRJScnJkycvXry4evXqkpIShFBxcfH69eu3bNlSV1dHtFEqlZ9//vnq1avff//9+vr6zhvZvXv3tGnT2Gw2Qujbb7+VSCR6vX737t0IofLy8h9++AEhJBAIJk6c+PHHH7vvtQEAHoSG+GKups5NJ4XqepiY3jkgEbrEsWPH/vrXv65evTowMNBsNu/cuTM/P5/D4SgUihEjRlRXVyOEjh49WlxcHBERUV9fn5qaKpVKO2zk8OHDmZmZxPLatWtLS0vVavWiRYsQQmfPnt2+fTvxp8zMzCNHjrjxxQEAHkgQxVPXuWMaCovRamoz88LYbtiX1/OeIdY6eO5nS5XKTXe2BrJphZkMxh9HOaLT6YcPH2YwGBqNZty4cVevXu3Tpw9CiMFgbN68ecuWLbNnz549ezbRuLW19Ysvvli2bFn7Ldy4cSMhIaHbvffv35846QQAkI4f4aaB1rRSAw/m43USr02EC/9Eb3FXpzWPiRidPo1paWkMBgMhVFFRodPpFixYQKxvbGwUi8UIoZs3b77wwgsSicRsNisUirlz57Z/utFo1Ov1fD6/43Y74fP5bW1tTnkhAAAH8SO59adb3LAjuJXeibw2EaaEkPxDicvl2hY4HM6OHTtoNFr7Py1atGjmzJmLFy9GCC1btsxi+UMHO5vNDggIUCgUkZGRCCEOh2Mw/D4Zr16v53DulU0rFIqwsDDXvyAAQPf4Yq6uxWg143Sma7+CNHU6fiQkQueAPkKX69evX3h4+IULF/r+D5EIGxoakpKSEEIqlerbb7/t/MShQ4farnkOGTLk6NGjxLLVaj169OiQIUOIh9euXUtPT3fHKwEAdIfGoPFC2W64m1Bdp4dE6CyQCF2OxWIVFhauXr163LhxU6dOTUpK+uyzzxBCS5YsycvLmzlz5ujRo//0pz91fuL06dO///57Yvn//u//SkpKMjIyLBZLUlKSxWJ5+eWXiT8dPXp0+vTpbns5AICuCSJ5alcXjuJIK4VE6DRee2mUXIsXL25/qXPYsGFlZWWVlZVKpTIxMTE0NBQh9MILL0yZMkUikQwePJjNZlutVoTQiBEjTp48STwrPz9/w4YNLS0tISEh4eHhZ8+evXr16tChQ0+ePBkREUG0kUgk5eXlU6dOdfcrBAA8AD+Sq6nTIRTkul1omwwsPyaTC4OrOQckQpewdRDasFisAQMGdFgZGxsbGxvboRmRJhFCAoHgnXfeOXfu3OOPP06siYqKotFotiyIEDp37tzmzZuJew0BAFTAj+TKrra6dBeaOp0gkufSXfgUSISUNmPGjPYPeTzeihUr2q+ZOXOmeyMCAHRDEMnTSPUIR66bOF5dpxdEwXVRp4E+Qk/C5/PffvttsqMAAHSFwaWzBExds6H7pr2lqdXx4YzQeSARusQHH3xA3BfhClu3bl2/fn0XDYqKip577jkX7R0A0C1BJNel9TLqeqiUcSZIhC4xZsyYWbNmuWLLarV63bp1zz77bBdtsrOzz5w5U1ZW5ooAAADdEkTz1LWuGmhNLzfSGTS2H3RsOQ0kQpewWCxmsxkhJJPJPv/888rKytdee+3tt99WKpUtLS3vvPPO6tWr7969SzSura394IMPVqxY8d577zU1Ndk2cuPGjTVr1qxdu1YikRBjbSOECgsLMzIygoKCEEI//fQTke327t2r0Wjq6+sPHDiAEKLRaHl5eTt27HDzqwYAEATRPHWNqxKhpk4viILros4EidAljh8//umnnyKEJBLJkiVLXn755djY2AsXLjz55JNz5szx9/dXKBTjx483mUwIoSNHjiiVypSUlKamprS0NGK8tNLS0nHjxnG53ODg4FmzZhFjbSOEioqKxo8fTyzv3Lnz1KlTCKEXX3xRqVSWl5e/+eabxJ/Gjx9fVFTk/hcOAEAICWJ46jodbnXJcMfqOj0fEqFTee3JtaJws1ne6J590Vjs0PmrEf3+vyrUavW+ffvCwsKmTJkiEon++9//ZmZm4jj+zTfflJaWpqSkLFy4kGiZl5dXUVFRVFT01FNPbdq0acGCBa+88gpCSCgU2spHS0tL7el9TEhIqKqq0ul0PB4cMAC4G5PLYPsxdc1GTOj8GeQ1tTrhcBfepOiDvDYRCsZOtarcNF0tjct/UBZECIlEImIsUKFQSKPRBg0ahBCi0WhCobClpQUhdPr06RdffFGv1wsEgrt3744ePRohVF5enp2dTWwhNTXVtjWNRoNhWLchEaN1q1QqSIQAkEIQjalrdK5IhOpaXfyTEd23A3bz2kTIEschcRzZUSCEEDEHReeHNBoNx3GE0DPPPPPRRx9NnDgRITRjxgxiSJrAwMDW1nv35CqVStvTw8LC5HI5sczlcvX63yvT2p//yWQyJpMZHBzsohcFAOiaXwxPJdGGDwt07mYNShPCESeI5dzN+jjoIySfTCYjppiQSCQ//vgjsfKxxx7buXOnWq3GcXzLli22xunp6deuXSOWU1NTi4qKiI5GhNC3336bkpJCLF+/fj0tLY3J9NofOgBQnIvqZVR3dX6x3V8TAj0CiZB8y5cvf+SRRyZPnvzEE08MHz6cWPnss88mJyfHx8f369cvOjradqo3Y8YMWxXMokWLgoKC+vXrp9PpMjIybt++/frrrxN/OnLkCAw6AwCJBJFcTYPBanZyvYy6RiuIgf4OZ8MpJjIysra21s7GWq3WbDbjOC4UChsaGlwZV88YDAatVovjuNlsViqVtvVyudxqtRLLra2tRqORWK6pqbl69arRaNRoNDqdztbeaDRardbvvvsuNTWVWGO1WpOTk0tKSmxtWlpaMAy7du2abU1ra2vfvn3lcrnLXh+O4/jmzZsXL17s0l2ArrW1tZEdgk/r9v2/sqFSJdE6d6fXt91WVKicu00P5cTPP1w6cwnbKNgMBiMgIMC2nrj/j+Dv729bjo6Ojo6ORgixWPcu/RuNxvnz548ePVoul2/btm3Tpk3EehqNtn379hs3bhBFNwih4OBgOp0eEhJi21pJSck777zTfl8AAPcjro468Z4/3Iqra525QUCAREhRLBZr2rRpZWVlGIYVFRXZpuFFCI0cOXLkyJHtGy9dutTPz8/2MCMjw32BAgAewC+Gp6rRiZy3QW2DgRPIYvJg9iUng0RIUTQaLTc3Nzc3157Gb7zxhqvjAQD0lCAGqz/T4sQNqmp00EHoClAsAwAALsEXcwxKk1lj6b6pfVR3tVAy6gqQCAEAwCVodJpfDNZ2R+srsBA6AAAVJklEQVSsDaprdH5wRugCcGnUVQwGw+7du8vLyxMTExcuXGirggEA+A7/vlhbtTZ4oF/3TbtjMVj1LUa+GGZfcj44I3QJHMcnT5588eLFsWPHHj16NC8vj+yIAAAk8O/Db7utccqmVDU6fiSXxnDZtPc+zGvPCN/4eUOtSuqefXGYnC2P/ItO+/1Xxc8//1xTU/Pjjz/S6fRJkyZFRUVVVFQkJia6Jx4AAEX4x/I0Ur3VjNOZjiawtipNQDzfKVGBDrw2Ef41eY7KqHbPvjgMTvssiBC6ceNGSkoKnU5HCPn5+SUmJt66dQsSIQC+hs6mY0KOukbr39fRHNZapYl+JMwpUYEOvDYRRgicePdOj2m1WrX69zSsVqvb3+cHAPAd/n34bdWOJkKrGVfX6vzioGTUJaCP0FVOnDhRXV2NEPrll1+IGXfJjggAQAL/PljrbUcLR1V3tZiIw+DAN7ZLeO0ZIelGjRo1Y8YMHo938+bNnTt3th9QDQDgO/z7YpWFdbgVp9F7303YCh2ErgSJ0FWio6O///77qqqq2NhYLhcqngHwUSwBk+XP1EoN/Mjefw+0VWkjx4V03w70CpxouxCTyezfvz9kQQB8XFB/gaK897V7uAVX1Wj9+kAHoatAInSJhx9+eMaMGWRHAQCghKAkgbJc1eunq2p0vHAOkwtjbbsKXBp1iWHDhjnydKvVarFYnDIYjcFg4HA4jm8HANBr/vH8W59ILAZr76pdWqs0AQ7ffQG6AGeErnL+/Pmqqqqu25w4caK+vr7z+gMHDjz22GOOx2CxWLhcbltbW6+3oFAoioqKHI8EAF/GYNP9orHWql4OMaMoUwX2Fzg3JNCecxKh1Wo9d+7ckSNHlEplF80MBkNxcXFLizPnJaGsTZs2HTt2rOs2q1atunz5snvi6Z3q6uolS5aQHQUAHi8wSaC41ZtuQrPGom0wBPSDM0IXcsKlUYvFkpOTI5FI+vbtO3/+/GPHjj300EP3bblq1ap33323oKBg7ty5ju+Xys6cOXP9+vXW1tY7d+6kp6fn5uZqNJrt27ffvHmzf//+S5YsEQgERUVFd+7c+fjjj8+ePTtx4sQJEybcd1MqlWrbtm2VlZUDBw5cvHgxj3dv7Pnvvvvu2LFjer0+MzNz9uzZarX6448//vXXX1ks1uTJkydPntxFeCtXrszPz9+zZ4/RaJw7d67tQu6lS5cKCwt1Ot2UKVOysrIQQtu3b1coFK+88gpCaNWqVRgG3fUA9EZQf8Gtf9cgJO7pE+U3VQGJfMdHaANdcMIZ4eHDh3/77beLFy/+5z//WbBgwerVq+/b7OLFi2fOnElJSXF8j9QnFAqDgoJiYmLS0tJiY2OtVuuECROuXbuWm5t748aNcePGWSyWqKgogUAQHx+flpYmFt//8DCbzaNHj/7tt99yc3MvXrw4adIkHMcRQqtWrfrnP/+ZkZGRk5NTWVmJEKqvr29qapo6derYsWOXLVtWWFjYRXjr16+fP3/+qFGjBg0alJ2dff36dYTQ6dOns7Oz+/fvP2bMmIULFxYUFCCEkpKS2Gx2WlpaWloakwk9ygD0El/MJaaP6OkT5TfanDJ5BeiCE77aDh48OG3aNOJMJS8v76GHHjKbzR2+NI1G43PPPVdQULBw4ULH92iPsj13tQ0G9+yLyWMMeblv+7tlExMTo6OjU1NTidrRH3/8saGh4ZdffmEwGNnZ2QkJCT/88MNjjz0WGhqakZExZcqUB2350KFDxHRONBpt0qRJcXFxp06dSk5OXrduXVlZWXx8PEKIeHpiYuKaNWt0Ol1jY+OSJUsKCwufeuqpLmJ+9dVXc3JyEEI1NTXvv//+7t27N2zYsGLFigULFiCEuFzu8uXL582bl5mZ+cEHH0AFLACOoqGg/n6KW2pxRrD9T8ItuLJCE58b4bq4AHJKIqytrR06dCixHBMTYzabGxoaoqKi2rd56623cnJykpOTu92awWD47rvvgoPvfVaioqJGjBjxoMYWi8Viuf/szwkzIy1Gq72vwTE0Bq3rMSMqKyuHDBnCYDAQQgwGIyUlpaKiwp5ymMrKytTUVBqNhhBis9mDBw8uLy/n8/mBgYFEFrSRSqUzZ85sbm6Ojo6Wy+XEvrpg+1+kpKS8//77CKGKioqXXnqJWDl06NDq6mqjsZtfr0R1a7evArhIF59/4AY9ff8DB/ClP8vDRwTY/5TWSg03jE3HaPCP7szO959OpxNfoV2wKxFeu3Zt2bJlndcXFBTExsYajUZboT+xYDD84VTs+vXrBw8evHTpkj370uv1hw4dsvWEpaamdnE11WAw4DjOYDCIC4btsfyY1JkJNyAgoH0ZkUKhCAwMtOeJgYGBHZ4YFBQUHBzc1tbW4bR73bp16enp7733HkLoiy++2LhxY9dbtpWSKpVKIpj2+1IqlXw+n81md/0BslqtHf7XwJ2MRiO8/yTq6fuP9WVp9+tVTRp2gL1nILKSVv/+PPgv35ed7z+bze62W8eu/0dcXNyqVas6rw8LC0MIiUQimUxGrGlubkYIdejx+te//hUZGblmzRqEkEQi2b9/v5+fX25u7n33FRAQsGvXrsjISHsCo9FobDabwWB0m/DdLzQ0VCKREMvjxo1btGjRlStXUlNTr169evHixb1793Zoc1+ZmZl/+9vfysrKBgwYcP78+Zs3bz788MNCoTApKem99977xz/+gRCSyWShoaEajSY0NBQhpNPpduzY0W14O3bs2Lp1q8Fg2Lt3LzFvcFZW1s6dO5944gkmk7lt2zaiWCYkJKSlpUWr1d63TIbJZEL5DIksFgu8/yTqxfsf8lCA5pYxcLx9Iw/jqK289k/zYzAMRqe6Dyd+/u0qlgkICBh7P0QQGRkZx48fJ1qeOHEiJSWlQ3Dz5s3Lzc3t27dv3759ORxOWFgY8ZXt3RYsWHDo0CGRSLR06dKoqKhdu3Y9/vjjycnJ2dnZH330UWxsLELopZde2rp1q1AoXL9+ffvn0ul04idMv379tm3blpmZmZycnJubu3fvXrFYTKfTv/jiiwMHDvTr1y85OXn69OnEpj755JOhQ4empKQMGTKE2A6NRmMymff9lcBgMAYPHhwfHx8TE/Pss88ihFasWEFccU1MTCwrK9u0aRNCKCYmJj8/v3///sHBwQqFwsXvGQBeLiw1sPnXVjsbt1ZpGBw6XwxZ0PVwh8nlcrFYvGzZst27d4tEosLCQmL9qFGjNmzY0KFxenr6vn37uthaZGRkbW2tnbvWarVmsxnHcaFQ2NDQ0PPY3cpisdTX11sslt490Wq1dlgvl8vbv2qTyXT37l3icnHX6HR6c3OzUqmUy+Ud/qRSqTqvvK/NmzcvXrzYnpbARdra2sgOwaf15v234hfW3NI06O1pe+sTSd1pWY934TOc+Pl3wu0TQUFB58+fZ7PZly9fLigomDVrFrF+0aJFY8aM6dB4yZIlPjszH51OJ87nevfEzid2QUFBQqHQ9pDJZMbExLDZbDs3GxAQEBQU1GGlQCDovBIA4Bw0FJbsL7PjpNCssShuqsKH2lVMABzknDvDYmJi1q5d22Flfn5+55Z//vOfnbJH4Ii1a9fy+TBQBQAkCEsLvLVPEj0prOtS86ZiZfBAPyYPBtp2Bxhr1Bf9/e9/t9XlAgDcSRDF4wSzmq90c1LYcF4uGgHXZtwEEiEAALhVzKPhkmNNuLXjTV82rb9pcCvy7wOXbdwEEiEAALhVQDyfHfDAk0Lcit/+jzQ2KxxR7qYwrwWJEAAA3C02W1jz4/1PCqU/y1l8ZmhyDwagAQ6CRAgAAO7m3wfjhbLvHmnssN6oMkt+ao7P7fEkFcARXjKfAIfDGT58eLcDbPos3IJbLTiDRbddbMGtuNWEM9j0Xl9+aW1tnTNnjrMiBMDX9M+Pvr71NsufGTn23gAjZo2l4lOJcHgQL5xDbmy+xksSYXFxsSPzsPuCxktK6c8tMZPC2P4sjdQg/aUlYVYkP6L3g1bQaLSAALh6A0AvMTHGwIVx17feNuusAfEYwlHll3WhQwJis8LJDs3neEkiDA0N9YVh2xzRty9SZ+gqC+twJi0sMmjg35L84xwdpk+lUjklNgB8EyeQNWhhnPQXueTHZkOrKf7JiOABMPUgCbwkEQJ7CKJ4Kcv7kR0FAOB3vHBO32nQI0gyKJYBAADg0yARAgAA8GmQCAEAAPg0z06E3333XVlZGdlR+K69e/c2NTWRHYXv2rRpk8lkIjsKH2U2mzdu3Eh2FL5LJpPt2bPHWVvz7ER45MiRS5cukR2F7/riiy8qKyvJjsJ3bd26ValUkh2Fj1KpVJs3byY7Ct9VVVX12WefOWtrnp0IAQAAAAdBIgQAAODTIBECAADwaTQcf+CcWKQQCAShoaF2jhra3NzM4/EEAoGrowL3JZVKg4ODORwYF5EcNTU1UVFRdDr8nCWB1WqVSCSxsbFkB+KjjEajTCaLiIjotmVeXt6bb77ZdRvKJcL6+nq9Xm9nY5PJxGAw4IuALAaDAbIgieD9Jxe8/+Sy8/0Xi8U8Hq/rNpRLhAAAAIA7wbkUAAAAnwaJEAAAgE+DRAgAAMCnQSIEAADg0zxmPsLm5uY9e/Y0NzdPnjw5MzOzcwOLxfLxxx+XlJQkJSX95S9/YbFY7g/SWzU3Nx86dOjmzZuBgYEzZ85MSEjo0MBqte7evdv2cNCgQaNGjXJvjN7s5MmTFRUVxDKTyZw3b17nNvX19QUFBUqlcurUqaNHj3ZvgF5u165d7YsKO3+8KyoqTp48aXs4derU8HCYZd4hTU1NxcXFEolk3LhxiYmJtvW1tbV79+5ta2t78sknR4wY0fmJRqNxz549lZWVycnJ+fn5dt5T4BlnhFqtduTIkeXl5bGxsXl5eV9++WXnNkuWLNm+fXtCQsKnn376zDPPuD1Gb/b888//8MMPYrG4qakpOTn57NmzHRqYzeaFCxfeunXr9u3bt2/fbmlpISVOb7Vv374vv/ySeG+rq6s7N1Aqlenp6XV1dVFRUU888cSRI0fcH6QXq66uvv0/S5cuLSkp6dDg3LlzGzZssLUxGAykxOlNJkyY8MYbb6xcufLcuXO2lTKZbNiwYc3NzWKxODs7+9ixY52fOHv27C+//DIhIWHLli3Lli2zd3+4JygoKBg2bJjVasVx/PPPPx88eHCHBlKplMPhSCQSHMflcjmXy62qqiIhUC+l0+lsywsXLpw3b16HBsSRr1Kp3BuXr5g7d+7GjRu7aLB58+bMzExieceOHRkZGW6Jy+dcvnyZx+MpFIoO6/ft25eTk0NKSN7KYrHgOD58+PB9+/bZVq5bty47O5tY3rJli+0zb3Pz5k0Mw1pbW3Ecr66u5vF4MpnMnt15xhnh6dOnJ06cSKPREEITJ068fv26QqFo3+DcuXP9+vWLiopCCAUFBaWmpp45c4acWL0Rl8u1Lev1+gcN5fPRRx+9//77V65ccVdcPuT8+fMbNmzYv3//feddOn369KRJk4jliRMnnjt3DqZncoU9e/bMmDEjMDCw859qa2s3bNhQUFDQ3Nzs/sC8z30vaXb4nJ85c8ZqtXZokJ6e7u/vjxCKi4uLjo6+cOGCXbtzOGB3kEqlYWFhxHJISAiTyZRKpQ9qgBASCoX19fVuDdE3nDt37uDBgy+88ELnP02YMKGlpeXmzZuZmZkbNmxwf2xeLCYmJiwsTC6Xr127Nj09XavVdmjQ/vMfHh5utVobGhrcHqaX0+l0hYWF9+2gDQgIGDRoUGtr68GDB5OSkjpfOwVO0eFzbjKZZDJZ+wYNDQ3tE0F4eLidicAzimWYTKbZbCaWLRaLxWJhs9kdGlgsFttDk8nUoQFwXHl5+YwZM3bt2tWvX78Of2Kz2T/99BOxPHv27EceeWTx4sV8Pt/tMXqnN954w7aQmppaUFDw/PPPt2/Q/gAhFuDz73Rff/11UFDQmDFjOv9p6tSpU6dOJZYXLVr0+uuvHzhwwL3R+YRuP+e9TgSecUYYGRlpS+zEglgs7tCgrq7O9rCurs6ewViB/SoqKiZMmPDOO+/MnDmz65YjR440m821tbXuCcynsFis9PT0zvUy7Q+Quro6FosVGhrq9ui8XEFBwfz584kOmi6MGjXq9u3b7gnJ13T4nGMY1uEyda8TgWckwpycnEOHDhGDcR84cCAzM5M427h69erdu3cRQmPHjpXJZMXFxQihysrKW7du2S4lA8fdvXv30UcfXbVqVX5+fvv1ly5dIj52Op3OtvLw4cMYhsXFxbk5SC9me3tVKtXJkycHDhyIEDIajcePHyfKlHJycg4ePEj0C+7fv3/y5Ml2zt8C7FRdXX3mzJk5c+bY1mi12uPHjxPnJbZ/EI7jR44cGTRoEDlRerucnJxvvvmGOOfbv39/Tk4Osf7SpUtEgnz00UevXbtG/FI8f/68VqvNyMiwa9NOqfBxNbPZPGnSpLS0tDlz5oSEhJw9e5ZYn5mZ+eabbxLLmzdvFovF8+bNi46Otq0ETpGVlcXn89P+57nnniPWp6SkbNu2DcfxXbt2DRw48Omnn87KyvL39//kk09IjdfbBAcH5+Tk5OXlicXixx9/3GQy4ThOHPl37tzBcVyv148ePXrkyJF5eXmhoaG//vor2SF7m1dffXXy5Mnt19y6dQsh1NLSguN4VlbWhAkT8vPzH3rooYSEhLt375IUpvdYunRpWloan8+Pi4tLS0sjvvO1Wu3w4cNHjx791FNPhYeHl5aWEo0HDx68Y8cOYvm1116LjY2dN2+eSCT68MMP7dydx8w+YbFYTp48KZPJxo4dKxKJiJXl5eV+fn62k9/S0lLihvqUlBTyIvVCFRUVKpXK9tDPz4+4xfXGjRthYWFEr3VxcXF1dXVAQMDQoUPhbmLnunPnztWrV41GY2JiYnJyMrHSZDL9+uuvycnJRC+IyWQ6ceKEUqkcP358+3oB4BQVFRX+/v62bx6EkF6vv379elpaGoPBkMvlFy9eVCgUkZGRI0eOhNE8HFdVVaVUKm0PExISiFpQ4kJIW1vbhAkTQkJCiL+WlpYKhULbx764uLiiomLIkCEDBgywc3cekwgBAAAAV/CMPkIAAADARSARAgAA8GmQCAEAAPg0SIQAAAB8GiRCAAAAPg0SIQAAAJ8GiRAAAIBPg0QIAADAp0EiBAAA4NMgEQIAAPBpkAgBAAD4tP8HgQR0KLa5vj0AAAAASUVORK5CYII=", "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" ], "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" ] }, "metadata": {}, "execution_count": 11 } ], "cell_type": "code", "source": [ "ρ = scfres.ρ[:, 1, 1, 1] # converged density, first spin component\n", "ψ_fourier = scfres.ψ[1][:, 1] # first k-point, all G components, first eigenvector\n", "ψ = ifft(basis, basis.kpoints[1], ψ_fourier)[:, 1, 1]\n", "ψ /= (ψ[div(end, 2)] / abs(ψ[div(end, 2)]));\n", "\n", "using Plots\n", "x = a * vec(first.(DFTK.r_vectors(basis)))\n", "p = plot(x, real.(ψ), label=\"real(ψ)\")\n", "plot!(p, x, imag.(ψ), label=\"imag(ψ)\")\n", "plot!(p, x, ρ, label=\"ρ\")\n", "plot!(p, x, tot_local_pot, label=\"tot local pot\")" ], "metadata": {}, "execution_count": 11 } ], "nbformat_minor": 3, "metadata": { "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.9.0" }, "kernelspec": { "name": "julia-1.9", "display_name": "Julia 1.9.0", "language": "julia" } }, "nbformat": 4 }