{ "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 one dimension and we\n", "show how to define local potentials attached to atoms, which allows for\n", "instance to compute forces." ], "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 a custom\n", "potential" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "struct ElementCustomPotential <: DFTK.Element\n", " pot_real::Function # Real potential\n", " pot_fourier::Function # Fourier potential\n", "end" ], "metadata": {}, "execution_count": 2 }, { "cell_type": "markdown", "source": [ "We need to extend two methods to access the real and Fourier forms of\n", "the potential during the computations performed by DFTK" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "function DFTK.local_potential_fourier(el::ElementCustomPotential, q::Real)\n", " return el.pot_fourier(q)\n", "end\n", "function DFTK.local_potential_real(el::ElementCustomPotential, r::Real)\n", " return el.pot_real(r)\n", "end" ], "metadata": {}, "execution_count": 3 }, { "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": 4 }, { "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;" ], "metadata": {}, "execution_count": 5 }, { "cell_type": "markdown", "source": [ "We define the width of the Gaussian potential generated by one nucleus" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "L = 0.5;" ], "metadata": {}, "execution_count": 6 }, { "cell_type": "markdown", "source": [ "We set the potential in its real and Fourier forms" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "pot_real(x) = exp(-(x/L)^2)\n", "pot_fourier(q::T) where {T <: Real} = exp(- (q*L)^2 / 4);" ], "metadata": {}, "execution_count": 7 }, { "cell_type": "markdown", "source": [ "And finally we build the elements and set their positions in the `atoms`\n", "array. Note that in this example `pot_real` is not required as all applications\n", "of local potentials are done in the Fourier space." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "nucleus = ElementCustomPotential(pot_real, pot_fourier)\n", "atoms = [nucleus => [x1*[1,0,0], x2*[1,0,0]]];" ], "metadata": {}, "execution_count": 8 }, { "cell_type": "markdown", "source": [ "Setup the 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", " PowerNonlinearity(C, α),\n", "]\n", "model = Model(lattice; atoms=atoms, n_electrons=n_electrons, terms=terms,\n", " spin_polarization=:spinless); # use \"spinless electrons\"" ], "metadata": {}, "execution_count": 9 }, { "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 `nucleus` 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 Eₙ-Eₙ₋₁ ρout-ρin α Diag\n", "--- --------------- --------- -------- ---- ----\n", " 1 +0.293919908373 NaN 4.29e-01 0.80 7.0\n", " 2 +0.441598831460 1.48e-01 5.35e-01 0.80 2.0\n", " 3 +0.260922834830 -1.81e-01 1.72e-01 0.80 2.0\n", " 4 +0.243052382673 -1.79e-02 2.68e-02 0.80 2.0\n", " 5 +0.242662237918 -3.90e-04 6.12e-03 0.80 2.0\n", " 6 +0.242635935708 -2.63e-05 1.86e-03 0.80 1.0\n", " 7 +0.242633436315 -2.50e-06 3.72e-04 0.80 2.0\n", " 8 +0.242633393265 -4.31e-08 1.76e-04 0.80 2.0\n", " 9 +0.242633377217 -1.60e-08 4.51e-05 0.80 2.0\n", " 10 +0.242633376619 -5.99e-10 2.01e-05 0.80 2.0\n" ] }, { "output_type": "execute_result", "data": { "text/plain": "Energy breakdown:\n Kinetic 0.0304528 \n AtomicLocal 0.0972437 \n PowerNonlinearity 0.1149369 \n\n total 0.242633376619 \n" }, "metadata": {}, "execution_count": 10 } ], "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-8, ρ=ρ)\n", "scfres.energies" ], "metadata": {}, "execution_count": 10 }, { "cell_type": "markdown", "source": [ "Computing the forces can then be done as usual:" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "2×1 Matrix{StaticArrays.SVector{3, Float64}}:\n [-0.03871253930283908, 0.0, 0.0]\n [0.03871254222307639, 0.0, 0.0]" }, "metadata": {}, "execution_count": 11 } ], "cell_type": "code", "source": [ "hcat(compute_forces(scfres)...)" ], "metadata": {}, "execution_count": 11 }, { "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": 12 }, { "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+gvaeTAAAgAElEQVR4nOzdd2BUZbow8Of0qZlJbySEJIReg9KbBZBV6oKuZa27ei246nqV3cunrrvq6orXtmtZRNcrYltdFIFVxBKKCEiXEkogIT2ZPnP698dhxxgCTJIp50ye31+TkzMzbybznue87XkJVVUBIYQQ6qnIRBcAIYQQSiQMhAghhHo0DIQIIYR6NAyECCGEejQMhAghhHo0DIQIIYR6NAyECCGEejQMhAghhHo0DIQIIYR6NAyECCGEejQdBcJdu3ZJkhThyYqiYHK4xFIUJdFF6NFkWU50EXo0vAQlXBQvQToKhBdddJHH44nwZJ7n8UKQWIFAINFF6NHw808svAQlXBSrgI4CIUIIIRR/GAgRQgj1aBgIEUII9WgYCBFCCPVoGAgRQgj1aBgIEUII9WgYCBFCCPVoGAgRQgj1aBgIe4pgA9+8N9J8BQiheHIf8XurMENCwtCJLgCKB1VRD7x5UvLL9VtaS+bncalMokuEEAIAkPzysVV1rsM+VYWRD5TSJirRJeqJMBD2CDVfNjM2evhvSmq+bNr5zJFR/1NGsdgZgND5SZJ00003hUKhdscVRSEIgiCIbr6+93iAtlKmdDbYKMAXYMnmuvmCSexnP/vZ9ddfH4tXxkCY/EItQs2XTcPuLiYootfFme6jgZY9nsxyZ6LLhZAB+P3+Dz744PXXX090QXq6r7/+ev369RgIURcd+6iu19QMUzqr/ZhV7mzY1oqBEKEIMQyzYMGCRJeipwuFQp999lmMXhz7x5KcIiiuw76ccWnhI+mD7d6qoOiNdMcrhBBKbhgIk5z7iN/Wy0xxP/6jSZZMG2Rv3OlOYKkQQkg/MBAmOdchn7Ofrd3BzJHOxu2uhJQHIYT0BgNhkms96E8tax8InWVW3iUGG/iEFAkhhHQFA2EyEzyS4BGtvUztjhMkkTYopfWgLyGlQgghXcFAmMxcB33OvjaC7GCpk73Q7DsZjH+REEJIbzAQJrPWjgYINfZCs/cEBkKEepCdO3eOHj06/OPDDz+8bNmyc5y/bNmyRx55JPblSjwMhMlLBfdhn7PM2uEvzdmc4BaloBznQiGEEkUUxZaWFu1xXV3d8uXLr7nmmnOcf80117z22msNDQ1xKV0iYSBMWoEGnmRIUxrb4W8JkrDmm33V7RNHIYR0bv/+/evXr9+6deuSJUt27NgBANu3b3/iiSeeeuqp48ePa+e0tra++eabS5Ys+ctf/lJVVXXmiyxfvnzWrFkmkwkAVq1adfz4cVEUX3nlFVVVKysrP/30UwAwmUyXX3758uXL4/e3JQgGwqTlrwnaCsznOMFeYPadxIT3CBnMl19++atf/eqhhx7KzMxUVXXZsmXXXHONzWaTJGnChAkHDhwAgA0bNuzcubO4uNjr9V544YVnxsLVq1dPnTpVe/yXv/zl+++/53n+1ltvVRRl69atzz77rParqVOnrl69Op5/XUJgirWk5a/lrbnt54u2ZSs0N+/CZfUIddorB5T3jinxea+5vcnbB7Zvsaiq+vHHH9M0zfP8JZdc8u2335aVlQEATdNLly595ZVX5s2bN2/ePO1kr9e7YsWKxYsXt32FPXv2aE85t7Kyst27d0fpT9EvDIRJK1Abyhmbdo4T7IXm45/Uxa08CCWNqblEsT1O+yUVdDTdbcSIETRNA8CRI0e8Xu8dd9yhHW9sbHQ4HABw8ODBO+6448SJExzHtba2zpo1q+3TVVX1+XwWi+W8726z2Xw+n6qq3d9nQ88wECYt/6mQJfdcW7qY0lhZUASvxNrxa4BQJ/R1EH0diSwAx52u2mazmWXZt99+m6JOB2YtQN51112zZs1atGgRACxevLi1tbXt0wmCyMjICB80mUxt95kKhULa2CEANDc3Z2ZmJncUBBwjTFZSUJaCsim145kypxHaMCEuokDIqHr37l1UVPSvf/0r9T94ngeAhoaG0tJSAPD5fB988MGZT7zgggvCfZ7Dhg1bs2aNqqoAoCjK6tWrhw8frv1q165dF154YZz+mMTBQJicArUhS64JzncbZ8NAiJCRkST59ttvP/300+PHj589e3ZZWdk//vEPALjzzjtvuOGGuXPnjhs3buDAgWc+ccGCBWvWrNEeL168uLq6euTIkQAwcODAYDB43333ab9au3btz3/+83j9NQlDaHcBepCenn748OG0tHMNa4UFg0GGYbROAHSm2o0t/lOh0gV55z6tea+nfkvrwFt6d+EtfD6fzdbxan0UB16v1263J7oUyc/tdvfu3dvl0lGSep7nJUmyWn9cIizL8tGjR5ubm8vKysKX0Orq6pMnTw4ePJhhGEVRLBaLJEkej0c7IRgMDhkypKKiIicnRzv/8OHDZWVl1dXV+fn52pHa2tpJkybt2bMn3FOaQG+++eZnn32mhXlNFC9BGEiSk/9U6NxTRjWWbC5Qj6m3ETISjuPCY4QaiqL69u3bt2/ftgd79erVq1evtkdomg6HSbPZ/PTTT2/atCk8uTQ3NxcAwnERADZv3vz000/rIQrGGgbC5BSoDWWOPP9ovimdFdyiIqkkneSD4QihdmbPnt32R5ZlH3jggbbzYsIxMunhGGEyUsFfF1GLkCAJLo0NNQtxKBRCSM9Yln3iiSdIsicGhZ74Nye9UKtAmyjaEtE6J3MmixsTImQgr7766o033hijF3/55Zf/9Kc/neOEDRs2XH/99TF690TBQJiEArW8JYLmoMacxWEgRMhAJk6cGKNQFAwG//jHP956663nOGfq1Kl79uzRcpwmDQyESchfG1G/qMaSxQUbsGsUIcMQBEFb/97U1PTCCy9UVlb+7ne/e+ihh5qamlpaWh5//PEHH3ywsrJSO7m6uvr555+/5557Hn/88ZqamvCL7N+/f8mSJY888khVVdWf//xnRVEA4P333x81alRGRgYAfPLJJ9u2bQOA5557zuVy1dTUvPrqq9pzr7vuupdeeinOf3VMYSBMQsEG3px9zqX0bZgz2UAjtggRMoyKigptFUFDQ8PixYvvvvvu0tLSAwcOzJkz55e//GVqaqqqqpMnT9aC5eeff+7z+caMGRMMBkeNGtXc3AwABw4cmDx5stVqLSgouPbaax988EEtEH766adTpkzR3uXdd9+tqKgAgEcffbS5ubmqquovf/mL9qvJkyd/8sknifjTYwVnjSahULOQMyY1wpOxaxShzvJ99VFg2/r4vJd5xGT7RWdd0h4IBJYtW5aTkzN//nyn07lmzZoZM2YAwAcffLBz584xY8bccMMNACDL8rRp0/bu3fvxxx/fcMMNzz333A033PDggw8CQF5e3mWXXaa92t69e6+99trzFqlv3761tbUul8vpdEblb0w4DIRJKNgkmDLOlWW0LcZGA4DolxlrnJIII2R05mET2OJB8Xkvyn6uYJOenq6t/HM4HCaTafDgwdrxrKwsrfFXUVFx5513CoJgsVhOnjxZXl4OAIcOHbrpppu0M7WEMppAIGA2n2vvNo2Wrdvr9WIgRDol84rMK53Ko61NHGX6nD8VPUIIAChnBuXMSHQpAADCubY14cUPBHE6a9jNN9/8zDPPzJw5EwCuvfZaWZYBwOl0hjNuh/esB4DMzMzwj2azORj8Mf9iMBgMx8jm5maSJDMzM2P0R8UfjhEmm1CzYEpjz5tltC1zJvaOIpScWlpasrOzAaCmpkbbdx4ALr/88ldffdXtdququnTp0vDJo0eP3rVrl/Z45MiRn376qZbCGwA++OCDcNtx165dw4YNS6aMMxgIk02oWTBnRDpTRmPO4oKNOHEUoST0wAMPzJw5c8aMGVdcccWYMWO0g9ddd92UKVPKysqKi4sLCws5jtNalgsXLgwHy5tvvrm4uLikpMTlck2ePPnQoUOPPfaY9qvVq1cvXLgwIX9OrKi6kZaW1tzcHOHJgUBAFMWYlsegTn7RePSj2k49pXGXe/9rVZ19I6/X29mnoCjyeDyJLkKP4HK5HA5HokvxE6IohkIhVVUlSXK5XOHjra2tsixrjz0ejyAI2uPa2trdu3eLouj3+4PBYPh87eS1a9cOGjQofHDMmDFbt24N/+jxeJxO57Zt28JHfD5fSUlJY2NjTP62s/vHP/5x3XXXtT0SxUsQjhEmm1CTYM3rXJcFJpdByEBomtY23qEoStuPXtN26krbnUlycnK0CTVtt+u58sorJ0yY4PV6//rXvz755JPh4y+88ML+/fsvuOCC8OvQNN32lffs2fPHP/5RW2uYNDAQJptQs5A+JKVTTzFncqFmQVVUgsTU2wj1CNddd93+/futVuvq1auHDRsWPl5eXq7NLA1btGhR20A4ZsyYcBdr0sBAmGy6MEZI0gSbwvCtoim9c09ECBnU5Zdffvnll0dy5pIlS2JdmITDyTJJRZVVwS1yqUxnn8ilMaEWMRZFQgghncNAmFRCLSLrYAiq0z2cXCrDt+LEUYRQT4Rdo0kl1Cx0rXvTlMbyrdgiRMgYeJ5//vnnDx8+3K9fv9tvvz2ZlvQlBLYIk0qoSTB1coBQw6UyoRZsESJkDHPmzNm9e/fMmTMrKirmz5+f6OIYHrYIk0qomTd3tUXY8J0r6uVBKCm9f+Djfx/bEJ/3uqho4lUD5rY9snXr1gMHDnzyyScURU2bNq2goEBL9RKf8iQlDIRJJdgkpBRbu/BELg1bhAhFanLh2CFZA+LzXqmm9omt9+/fP3ToUC0XjNlsHjhw4MGDBzEQdgcGwqQSau5q16iTEb2SKqtdmGiDUE+TacnItCRsRXkgEPB6veEfPR5PSkrnlg6jdnCMMKnwLaIprSuBkCAJxk7zbpwvg5ABbNy48cCBAwCwdevWqqqq0aNHJ7pExoYtwuQh+mWCISiuizc3pjS2y3EUIRRP48ePv/7662maPnjw4CuvvJKaGulG3KhDGAiTB98qcM5OL6UP41KZUKvggK4MMSKE4ik7O/utt946fvx4fn4+x0W6Czc6GwyEyYNvFbsTCLUWYRTLgxCKHZIki4uLE12KJIGBMHnwLpFL7XrHJpfKeI4GolgehFAsjB07tlevXokuRVLBQJg8eFe3WoRcGhvajksJEdK7ESNGjBgxostPl2VZkqSodKgGg0Gz2dz910k4nDWaPPjWrqTbDjOlMdg1ipAhbN269fDhw+c+56uvvjp58uSZx1evXn3JJZdEpRgpKSmNjY1dfrrH4/n444+jUpJuwkCYPHhXtwIh52QEj6gqahSLhBCKhRdffHH16tXnPufRRx/dvHlzfMrTNdXV1b/+9a8TXQoA7BpNJt2cLENQBGOjBbfUnWiKEIq1LVu27Nixo6ampq6ubuTIkQsXLgwEAtrO8qWlpXfddZfD4Vi3bl1lZeWKFSt27NgxZcqUGTNmdPhSXq/3+eefP3ToUL9+/e666y6bzaYdX7Vq1dq1awOBwJQpU2644Qa/3//6669v376dpunp06efO7vpkiVLFi5cuGzZsmAw+Mtf/nL8+PHa8W3btr311luBQGDmzJmzZ88GgL/+9a9er/fBBx8EgN///vd2uz2aH1NnYIswSaiyKvklNqVbdzZcGouJ1hDSuczMzLS0tIKCgvLy8j59+qiqOm3atF27di1cuPDEiRPjx48XRTE/P99utxcXF5eXl+fn53f4OrIsT5kypbKycuHChQcOHJg6daosywDwhz/84cEHH5w8efLChQu1ztWGhoampqYFCxZMnz794YcfXr58+TmKt3Tp0ptuumny5MmjR4+eM2fO1q1bAWDz5s0zZswYMGDApZdeet9997344osA0L9/f4ZhysvLy8vLGSaR99/YIkwSvFtk7DRBditB2ulhwpJoFQqh5FS9vrFuS2t83iv7AmfBtKy2R0pKSoqKioYNG7ZgwQIA+PLLL48dO7ZhwwaGYS677LKBAweuWrVq/vz52dnZY8aM0c7p0Lp161wu16uvvkpR1PTp00tKSj777LOJEyc+9thjO3fu7N+/PwDMnDkTAPr06fPQQw/5/f76+vq77rrrnXfeufHGG89R5gcffHDu3LkAUF1d/eyzz7711ltPP/303XffrXWE2u32m2+++Y477rjooov+9Kc/naOEcRNpIKyvr9+2bdsPP/wwZsyYCRMmdHjOvn37Xn/9dUmSrrnmmlGjRmkHZVl+4403vvvuu759+952220WiyU6BUc/JXRvyqhGW1MflfIglMRyxqZlDHfE571oE3XuEyorKwcPHqy1qAiCGDFixKFDhyJ55cOHDw8bNkxL3k1R1LBhww4dOpSfn282m7UoGFZfX//zn/+8tbU1NzfX5XIJwnmuEkOHDtUeDB8+XJsOc+jQoRtuuEE7WF5eXlNT4/f7IylkfEQaCG+99dbGxsZTp075fL4OA+GhQ4fGjRv329/+1mKxXHzxxevXr9di4b333rtx48Y77rjjvffeW7du3bp166JZfPQfoVaR63Z2NM7J+E4Go1IehJIYbaFoy3niU9w4nc7W1h+bp62trRFmXEtNTT3ziWlpaT6fj+f5tusrli5dOmLEiOeeew4A/vnPfz788MPnfmW32609cLlcWmFSU1NdLlf4jUwmk9lsJgi9pPiPdIzwo48+2rhx49ixY892wvPPP3/VVVctWbLkvvvuu+uuu5YuXQoALS0tr7766sqVK2+88cYPPvjg22+/3bFjR3QKjn6qmzNlNJyT4d1SVMqDEIqdzMzMEydOaI8nTpx48ODBTZs2AcCePXsqKiq01RFtz+nQ1KlTt23bpl2Tv/vuu++//37q1Kl5eXnl5eWPPfaYdk5DQwMABINBLWiFQiFteO/c/va3v6mqyvP8smXLtHk6l1122csvvxwKhQDg2WefnT59OkmSGRkZLpcrHDUTKGqTZSoqKi6++GLt8cUXX/zNN98AwPbt27OyskpLSwHAbDaPGzeuoqIiWu+I2urmanoN62AEFy4lREjvbr755m+++aaoqOj+++/Pzs5+4403FixYMGTIkEsvvfTFF1/ULrmLFi16/fXXi4uLn3rqqbbPpShKa+0VFBT8/e9/v/zyy4cMGTJ79uxly5b16tWLIIgVK1asX7++qKhoyJAh8+bNA4A77rjjX//617Bhw4YNGzZq1Kjwrk9na9U5HI5BgwaVlJRkZWXdfvvtAHD33XcXFRWVlJSUlJTs2rXr+eefB4DMzMw77rhj2LBhJSUlWsRNFEJVO7Fu7Oqrry4rK+uwXZyTk7NixYqLLroIAPbt2zdixIhQKPTWW28999xz3333nXbO9ddfn5eX9/jjj3f44ikpKZMnT2bZn/TvvfTSS1ZrB2mgg8EgwzA0jZN9Tqv8R13m6BRHv24NwcpBZe/TJ4b9T1EkJ/t8vvBkaxR/+PnHh8fj6d+/f7hbT8/q6+uzsrI629+oqmpDQ8OZT/R4PIIgZGSc3nZRUZS6urqMjIx2l+gzWa3WgwcPOp1OSZKczp/sKhwKhUKhULuDEXrzzTfXrl376quvho9EWAVMJhNJnqfJF7VAwrKsJJ3uVRMEgWEYkiTbHtSOnyOvD0mSc+bMafeH2Wy2Dj93RVEwELYleWRbpqW7aZM4AAVoYCLZy0kURUx7n0Dnrk0oWs573deP7OzsLjyLIIgOn9hus1+SJPPy8iJ/2Q5DlMlkMplMnS1h2zK0/c5HeAmK5M4gaoEkPz+/urpae1xdXa2tXMnPz6+pqVFVVStKdXW11mTsEEVRc+fOTUtLi+TtqP+IRtmTgeCSzOlc9z8Q1sHIPoW1nL+XFT//xMLPPz7wQ+6s559/Pkb7IxIE0fbfEcUq0K0xQpfL9emnn2qP58yZs3LlSq2jdeXKlXPmzAGA0aNHUxT1+eefA8Dx48e3b9/+s5/9rNtlRu3JIUVVVdoche8E58R96hFCXXfTTTd1OJ6lZ5EGwqeffnrUqFHr1q175ZVXRo0a9f777wPA/v37w4Ht1ltvraurmzJlyvTp07ds2XLPPfcAAMMwTz755NVXX3311VdPnDjx/vvv71TjGkUo1CpEKy8a62AEF04cRQj1IJF2jS5cuHDKlCnhHwsLCwFg+PDhe/bs0Y44nc7vvvvum2++kSRp8uTJ4b05rrvuuokTJ37//fe///3vBw0aFM2yo/+Iymp6DedksEWIEOpRIg2EBQUFBQUF7Q5aLJbBgweHf+Q4rsPdPYqKioqKirpaQnR+vFuKViBknYz/VCgqL4UQQoaASbeTgeAWWUeUWoS4lBAh1MPg8oNkILhFW2F0kriy2DWKUBssy2orwdsdD0+GPzdFVAiK6GY2/NPvqKiqrJJMD229+Hy+WbNmxejFMRAmA94tpTui86/kHDS2CBEKM5vN69evPzMNWCgUYhjmvNP39y8/0Xt6ljWv64vnwoJNwpH3Tw2+raj7L2VQsZtriYEwGUSxa5Sx0rKgKqLSY288EWonPT09PT293cEIk1s1UUK/oaXd3ChUI+cpng+V4uLi7r8UagcvdsmAj14gBALYFFrA1NsIdZsqq1JAZmzRWfRNmUgAkHklKq+G2sJAaHiKpCq8wkRvUxjOyfDYO4pQtwkeqfvbZbfFOhgBh/BjAAOh4QlukU1hIHobe+F8GYSigneJXLS6agAAgHNg4qeYwEBoeIJbZKM0U0aDKygQigrBLbJRWuCrYR0MDlvEAgZCw+PdUtQGCAEA040iFCV8a5RbhNg1GiMYCA1PcMegsmG6UYS6jXeLrDO6vTV4kxoTGAgNT/BIUe4axckyCEWD4JZi0CLEm9Tow0BoeHz0Mm5rsPsFoajgYzBGiC3CWMBAaHhRXE2vYe20GJBVRY3iayLUAwluMSpL6cMw8VOMYCA0PMEd5a5RIICx0YIHe2AQ6gYVRK/EpkTzJpWx06Ifb1KjDwOhwakgRLuygXbjiYEQoW4Q/RJlokg6eit8AQiSYKyU6MW6GWUYCI1NDMgUS0S3sgEOEyLUbYJbim6/qIZ1MjzOl4k2DITGFvUBQg2mG0WomwRPTOomhzepMYCB0NhiFQgdjODByoZQ1/ExahE6aAyEUYeB0Nh4t8RFd6YMAGCLEKFuEzxRznSh4RzYNRp9GAiNLWZdowyPLUKEuiFWY4TYNRoDGAiNLWZdo9giRKhbYjRGiF2jsYCB0NgET0zuOjkcI0Soe2JUN9kU7BqNPgyExhbNvenboC2UIqqKiHthI9RFgitmLUK8SY02DITGJngkLgZ3naDNl8E19Qh1iaqoYkBmbFTUX5k2UaCAzONNajRhIDQwVVGlgExbo1/Z4PTEUbzxRKgrBI/E2GiCjHKmCw2TQmNymejCQGhggkdi7bGqbKyDwRYhQl0T9c3R2mJTcA+KKMNAaGCCR2Ji0y8KuAUoQt0Q9e2y28Jhi6jDQGhgglvkop1uO4xJwRYhQl0Uo0WEGg7ny0QbBkIDE7yxrGyYXAahrhI8YtT3hAljUhgRb1KjCgOhgcVooZIG040i1GUxHiOkeQyEUYWB0MAEt8jE7K4TE1gg1GW8O4YtQjaFEfEmNaowEBqY4IlJxm0Ni7l9EeoqwR3jFiHWzajCQGhgMR2HoFiSIEEKyTF6fYSSWIy2ntBgcpmow0BoYDGdmQantwDFG0+EOkeRVEVQaHNMMl0AJpeJAQyERqUqqhSUGVsMAyGbgjeeCHXa6a6amCS6OI3BpYRRhYHQqLS0MjGtbGwKtggR6rSYDhBq2BSc1B1NGAiNKqZrJzSsA+86Eeq0mA7eazC5THRhIDQqIZbzszV414lQF8R68B4wuUy0YSA0qni0CPGuE6HOi0PdxOQy0YWB0KgEj4iBECEdEjwx2ZK3LUwuE10YCI1K8Eixr2wMJpdBqLPi0luDyWWiCQOhUcWpaxT3/0Sok+IzWQaTy0QRBkKjikNlI1mSJAlMLoNQp8Rhsgwml4kuDIRGFYe1SqCt28UbT4QipkiqIsYwrYwGk8tEFwZCQ1JlVQrKjDXmgZDDpYQIdUYc0spomBRaxJGLKMFAaEiCN+ZpZTQMLiVEqDPi0C+qYVMYHueyRQkGQkMSPBITl8rG4QoKhDojplvytsViizB6MBAakuARuRjPlNEwKTTO0kYocnGYxabBZb5RhIHQkOLWImRTcHtehDohDuuaNJgBMYowEBpSHCsbdr8g1AlxbRHiTWqUYCA0JDF+lQ0H5BHqhPhNlsEZ3dGDgdCQ4tYixOUTCHUKdo0aEQZCQ4pbZSNZkiABk8sgFCGcLGNEGAgNKQ6bEYaxuOELQpFRJFURFNoS27QyGtpMaW8Xh/dKehgIjUdVVDEgM7Z4VDbA9L4IRSxuaWU0rB3T4kcHBkLjEX0yY6UIMk61jXXghi8IRSRuYxYanC8TLRgIjSdugxAa3AIUoQjFbcqohk3BPSiiAwOh8cT7rhOTyyAUGdErxifThYZNYbBFGBUYCI1HcIvxDYRY2RCKCO+WOEdce2uwbkZFJ66ntbW1H374IUmS8+bNy8rKavfbffv27d+/v+2RuXPn0jS9ZcuWkydPakdYlp09e3Y3S4wEjxTnrlGsbAhFQvCI5ixr3N6OTWECDb64vV0SizQQHj169MILL5w3b54syw8//PC2bdt69erV9oR9+/a9//772uPKyspTp07NmzcPAJ599tlDhw6VlJQAgNVqxUDYfYJHsuWb4vZ2mFwGoQgJbik+2fA1mGUtWiINhM8888y8efNeeeUVALjuuuteeOGFJ554ou0JCxcuXLhwofZ4/vz5l1xyCUWdnt9/880333777dErc08neER2oD1ub4fJZRCKkOCJ67AFg701URLpGOG6deuuuOIK7fEVV1yxdu3as53Z1NS0evXqG2+8MXxk3759//d//7d169buFBSFCR6JtcevspEsSZIEJpdB6LwEj8TGcYyQS2EEL/bWREGk19NTp07l5uZqj/Py8k6dOnW2M994441Ro2emursAACAASURBVEYNGDBA+9HhcNTW1v773/9+4IEHysvLP/zww3BLsR1JkpYuXWo2m9sevPPOO02mDroBeZ5XFEWWe+LVmXeLqknheT5u78ikUL7GgDmL/UkxeJ5h4lfnUTs8z7Mse/7zUGyceQnS8rzIpCTz8WqlUaDwStAfIul4reHXkwgvQSzLEsR5Pp9IAyFBEKqqao9VVT3H6y5fvvzee+8N//jSSy9pD9xu99ChQ99+++1rr732bM91uVztru/hN0WnqSDHMa2MhrHToldqFwgRQm2JHomx03FLKwMAQABtoySfzDrj10WUlCL9+HJycurr67XHdXV14dZhO1u2bKmqqlqwYMGZv3I4HBMmTGg3s/QnRaHpP/zhD2lpaZGUR1EUhmFousf9+0WvRJkokyV+k2UAwORk1QDBcdxPSiKK7Y6geBIEAT//BDrzEsSHZM7JxPmfwjlYCBJcdk/8JkTxEhTpGOH06dM/+eQT7fEnn3wybdo07fHhw4eDwWD4tNdee+3KK6+023+cyhFu0vE8/91332nTR1GXxXk1vYZ1YAILhM4jzimfNLi6KSoivaT+5je/GTNmDEEQsiyvXbt227Zt2vFBgwZ99tlnkydPBgC/3//OO++sWbMm/CxBEIYOHXrJJZeYTKa1a9dmZGSco18URULwiPEcjdewKUyoRYjzmyJkLIJbYh0JuUnFQNhdkbYIS0tLd+/ePWjQoOHDh+/evbugoEA7/u677w4aNEh7HAgEli9fPm7cuPCzGIZ59dVX+/fvn5+f/+STT1ZUVGBnTjfFecqoBpcrIXReCaqbuD1vFHTi35aXl3fnnXe2Ozhnzpzw48zMTG0RfRhBEBMnTpw4cWJ3iojaEtxiIu46sbIhdB6CR7Tk2OL8pqyd9hwLxPlNkw/mGjWYOOdX0+A4BELnFeetJzR4kxoVGAgNJkGTZRjBLQKuZEHo7BI0fo/DFlGAgdBg+ERUNpImSIaUgj0xfQFCEUpUi5DHFmG3YSA0mIRUNtAmp2HqbYTOQhEURVZpc1wzXQAAY6EUXlEk7K7pFgyEhqKC6EtQIMRdCRE6u4SMWQAAEDhxNAowEBqJ4JMYC0WQCcgryDpo7IFB6Gx4dwLGLDQ4l637MBAaiZDIysaIWNkQOouEtQhxvkw0YCA0kgRXNgyECJ2F4BG5uC/w1eAKiu7DQGgkiWwRaisoEEIdET0SE/cFvhq8Se0+DIRGktgWIY+VDaGz4D0Sl7C6iTep3YWB0EgSPQ6BlQ2hjgluMWEtQsy73W0YCI1EcCdgnxcNm8KIXgmTyyDUoYTepGKLsLswEBqJ4EnAPi8akiYojhT9eOOJUAcEt8g5sUVoVBgIjSSBk2Xg9HwZrG8ItSeHFCCA4hJzOaXNlCKrsqAk5N2TAwZCw1AVVQrKjDXeOZzCWAfDYw8MQmfg3SKXuDtUAGDttOjFm9Suw0BoGIJHoq10QtLKaDhMN4pQRxLbVQO4uqnbMBAaRgIHCDXYNYpQhwR3ousmLiXsHgyEhiG4RS5BU0Y12DWKUIcS3jXKOWism92BgdAwEjg/W4Ndowh1KOFdowymAu4eDISGIXjExAZC7BpFqEMJ7xrlMPFT92AgNAzBIyX2rpPDrlGEOpLwrlEG19R3DwZCwxDcCW4R0lZKEVVFxOVKCP2E4Epw1yiHa+q7BwOhYSS8RQg4OQ2hM6iKKgZkxpawBb6AWda6DQOhYfCuBHe/gDZx1IX1DaEfCR6JtSdygS8AUCYStAQ3qEswEBqDIiqKqNKWRN51wumJo9giROhHgltM7EwZDefEIfyuw0BoDIJbYlNoSORNJwAmsEDoDLw78WMWgHWzezAQGgPvSlhu+7ZYXLeL0E8JiZ4yquFw2KIbMBAaQ8JX7Go4XEqI0E/ppGuUxWGLbsBAaAy8W+J0Udmw+wWhn9BL1yiOEXYDBkJj0EmLELtGEWpHP12jeJPaZRgIjYF3i6wOxgg5ByN6JVATXQ6EdENHXaM4RthVGAiNQXCLeugaJSiCMlGiH4ciEDpNcEtsQreF0XAOhscxwq7CQGgMvEsX4xBwuncU6xtCAABSUCZIguISfyFlbLQckhUJu2u6IvH/P3ReqqKKvgTvwRSGQxEIhSV834kfEcDYacGDdbMrMBAagOiTGSuV2BxOYayDwaEIhDQ6mcWm4Zy4uqmLMBAaAO/SxUwZDZeK63YROi3UKnKpeqmbuLqpyzAQGoCgg3TbYZwTAyFCpwn6SPmkweQyXYaB0AB4nXW/YGVDSKOT3Ica1klji7BrMBAagKCPtDIazsnwrVjZEALQAqFuukZxBUWXYSA0AF0NyLNORnCLuKYeIdBbixDHCLsKA6EB8G4dVTaSJigTJfrwxhMh4F16ukl10Dhs0TUYCA1AJzmcwnCYECEAkPwySetiNb0GMyB2mV7+hegcBH2ktw/DFRQIgTZ4r5uuGghnQMTems7DQKh3UkAmKIJidfSfwhYhQqCzAUINh5sxdYmOLq+oQ3rrFwVt5zMMhKjHEzySfjJdaHB73q7BQKh3oVaRS2UTXYqfwBUUCAEA3yqadLN2QsOlMnyrkOhSGA8GQr3TZ/cLphtFSHDrrkWIwxZdg4FQ73g9JTPUYGVDCHSWX03DOVnsrekCDIR6p8MWIeugRZ+kKjhNG/VofKvu6iaXisMWXYGBUO902CIkSIK20qJXTnRBEEocFUSvvtY1AfbWdBUGQr3jWwW9BUIA4JyM5MFAiHou0S9TZoqkdbFLaNjp3hoZe2s6BwOhrqmKKnol/ezBFIZbgKIeTnDpZm/6NgiSYOw0LiXsLAyEuiZ4JNpKE5S+7jrh9MRRDISo59JbWpkwnNTdBRgIdU2HM2U0XCojYtco6sFEt8Q6ddciBAAuFSeOdhoGQl3T4fxsDYtdo6hnE9x6HLMATAXcJRgIdS3UInJpeqxsJicjurFFiHouwSXqtEWIiZ86DwOhrum3azQdu19Qj8a3SKY0feU+1GCLsAswEOqabrtGGSulyqoUwkYh6qF4l8imYoswSWAg1DV9LiLUsE6Gb8H6hnoimVdAAsamy0CYyoQw73YndeIf2djYuHnz5ry8vPLycoJoP6G/vr6+uro6/OPgwYM5jtMeV1VV7dy5s2/fvgMHDux+iXuUkEt3W0+EsU6KbxWteaZEFwSheAu1CIyTSnQpOkabKQCQQjJt0mkJdSjSFmFFRcXAgQPfeOONa6655rrrrjvzhLfeemv69Om3/kdDQ4N2fMWKFaNGjVqxYsW0adMeeeSRqBW8B1AERRFUxqLTbzOXxoSa8cYT9UR8i+4SH7aFvaOdFWmLcPHixUuWLFm0aJHb7e7fv/+mTZvGjRvX7pwrrrhi+fLlbY+Ionj//fe/+eabM2bMOHr06ODBg3/961/n5uZGp+zJLqRlGdXdYvrTGAeFO5+hninUIrB6DoSpLN8qWnOxtyZSEbUIGxsbKyoqrr76agBwOBwzZ8788MMPzzzN5/Nt2bKlqqoqfOTbb78VBGHatGkAUFxcXF5e/sknn0Sp5MlPh7nt2+JSmRCOEaIeiW8ROF3OlNFgcpnOiuh/WV1dbTKZMjIytB8LCwsPHjx45mk7duy4//779+/fP3r06Pfee89qtdbU1PTq1YskT4fbgoKCmpqas72LLMsffvihzWZre3DWrFks28EgmSzLJEmeOVSZTIItPOukZVmnMzOpFCLYzOu2eElPlmX88BMl2Cyk5JplWdbnJYhxUKFWIem/HhFWgUgiRUSBUBAEhvmxacIwTCgUanfObbfddu+99wKA1+u9+OKLn3jiiUcffVQQBJqmz/3EMEVRPvroo3Zh76KLLrJarWeezPO8oijJ/Z8ONocoG8HzfKIL0jHVIvPNom6Ll/QEQcAPP1FCzbzDyvG8Tm8EKRvhORxM+q8Hz/NtA9PZmEym6ATCnJwcn88XCoVMJhMANDY25uXltTvHYrFoD+x2+y9+8Ys1a9ZoT2xqagqf09jYOHLkyLO9C8Mwb7zxRlpaWiRFIgiCYZi2UTb5yN4WR19r+IPVGyVNIYhWluC0WWoozmRZ1u13I+kJLtmRazdZOH1egsRstWW7L+m/HoqiROtvjGiMsKCgoLCwcMOGDQCgquqGDRvGjx+v/UpVO9j46ocfftBmxIwcObKuru7o0aMAEAqFNm3aNGHChKiUuycINQv6TF0RxqXhMCHqcaSQrCoqZdbvImxTOoszujslotsZkiR/+9vf3n777UuWLNm0aVMgEJg3bx4AbNq0afz48VosvOWWW0pLSzMzM7ds2bJy5cqKigoASE9P/9WvfnXVVVfdcccd77///gUXXFBeXh7TvyeZhFpEU7quA6EpjeVbBFs+Tk5DPQjfrPeKyTkYKSArokIy+o3WuhJpu/7OO+/Mzc1dv359fn5+RUWFtli+qKjoz3/+s3bCvHnzvvzyy+PHj/fu3Xvfvn2FhYXa8Weeeea1117btGnTpEmTbr/99lj8DUlJkVTJL7Epeux4CePSmFAL3niiniXUot98T6cRp1dQmLO4RBfFGIgO+zYTIj09/fDhwxGOEQaDweQeIww28PuXVZUvLkt0Qc7K5/N5dvChZqF4Li4MTQCv12u32xNdip6o5qsmvlXMne7U8yVo3yvH8yampw5I5m+Iz+drt8qgy7DhrFOhZkHn3S9weikhtghRz8K36jqtjMaUxoaacfw+UhgIdSrULHL6nikDAKZ0FvNuo55G/7PYAIBLZ/EmNXIYCHUq1GKAymZKw8qGehy+VcS6mWQwEOpUqMUAXaOUiSQoQvTrcU0xQjESahG4NN13jaZjTvxOwECoU6FmwZSu98oG2oqlJqxvqKcQfRJBEfpPIoFLCTsFA6FO8bpfRKgxZ7LBpiTP5IRQWLBRMGcYYE0CbaYIgpAC2FsTEQyEeiQFZFBB/3edAGDO4EKNeOOJeopgE2/ONMAdKgCYcMfQiGEg1COj9IsCgCkDW4SoBwk1CqYMgwRCnDgaMQyEemSImTIacyYbxDFC1GMEm4zRNQraCgpcShgZDIR6FGoWOKMEwgwu1ICBEPUUoUYDdY1iizBSGAj1KNRsgIVKGtpKAQm4ggL1EMEm43SN4hhhxDAQ6pGBukYBwJzBhhpxmBAlP8ErkQxpiFlscDrxEwbCiGAg1KNQs2DS/YrdMFMmh8OEqCcINfJmgzQHAYBLY/lWUVX0sq2CnmEg1B1VUQWXARKNhpkz2CC2CFEPEDTOlFEAIGmCsdG8C+fLnB8GQt0JNQusgyFpItEFiZQ5kw3iUkLUAwSbBGPt8GfO5II4ly0CGAh1J9hosMpmyuAwyxrqCYzVNQoA5izsrYkIBkLdCTYYZn62xpyJlQ31CAaaMqoxZ3JYNyOBgVB3go2COdNILULaTBE0IXqlRBcEoVhSIWS8QMhi12gkMBDqTrCRN2cZqbIBgDkDJ46iJCd4JZIjaZMx1k5ozFnYIowIBkLdCTbwxmoRAu5BgXqAoNEGCAHAlMaKPlkRlEQXRO8wEOqLzCtySOEchllEqDFn4eQ0lOQMN4sNAIAAUzqDvTXnhYFQX4KNvCmTBcMsnTjNks0F6kKJLgVCMRSoDVlyjBYIcb5MZOhEFwD9hBH7RQHAkmMK1GFl67SgBC286hLAJYBHAI+ougXwiRCUwSuqIRmCEkgKeEUAAEEBv/RjlhBZZinqxwlKZorQRq8sNHAUcBRYaLDShJmGFAbsDNgZIoUFJwtOFtJNhBWrficF6vi0gfZEl6LTzFkcLvM9L6wN+mLI7hcAcwYreiVZUCgW+xh+JChQG1BP+qAmoNYGoDag1gWhMaTWBaAxBM28qqqQxhFOFpwcOFiwM4STBSsNZhpSWUILZjQJdgYAgCHBRv/48QYCvMXy45BVUFZDMgCAXwRBAV6GgAR+SXUH4LAIHhG8guIRwSVAKw8tvCqrkM4RGSbINkO2mcg0QS8rkW2GAhuRZ4FeVsJQk0LiIVBn0BYh6670J7oUeoeBUF+CjXzqAOPddQIB5kw2WM/bCsyJLkpi1PjVSg9UetSjXvWYF6p86nEvNPNqtpkosEK+lci1QJ6FGJQKmSYyywyZJsgwEZZu1D+vV7Hb2/ahd64/PSRDc0ht4qE+CA1BtSEINQF1RzOc8Ck1fqgJqE4WCm1EkY3oY4c+dqI0hShNgUIbQRqt3z4qJL+siCprtMF7ADBncnWbWxJdCr3DQKgvwQYhb5Lx7joBwJJr8teFekIgVFQ46lV3t6gHXLCvVT3gVg+5VSsNfR1EaQpRYid+VgBFdrLIBjkWgtJr2DBRkG8l8q3aTx2UsjYAJ3zqcZ96zAs7mtV3jyqVHmgMqX0dRH8HMcAJA1OJwalEmYNgekAvgL8uZMk1ZMXEiWyRwECoL8Emg6WVCbPkcMk6TBiSYXeLuqNJ/b5Z3dms7mtVs83E4DRioBNmFBC/GUyWOQiHIf9p55JrgVwLMTrrJzHSL8Eht3rIre5vVd8/pj60XTnhV/s5iGFpxIh0ojyDGJ5O2IzXajq/QB1vyTEluhRdwWg7hvokxoZX+7PCj0ZHBI+Rdjtrx5Jjqt2YPD0wh93q5gb120Z1a4O636WWOYiR6cTIDOKXfcmhaYQ9Ga/1kbDSMCKdGJH+Y3QMybC3Vf2+Sf2+WX37qLK3Re1jJy7MJMZkEWOziYHOJOlKNegAocacyQUbBQyE54AfjY4YdMqoxpLDBWoNvIJCBdjTon5Zq35Tp35Tp7AkMT6HGJNJXFtKjkjHmSNnZaJgVAYxKuN0uJMU2NOqftugbqxX/7JHaQiqE3KIiTnk5ByiPIOgDduJGqjj04c6El2KLjJnssEGPqWPJdEF0S8MhDoSqAtZDXvXaUplpZAsBWVjtWiPedXPatTPa9QNtUo6R0zJJeb0JpaOoQusSdGQiTuaPN1kvG0AAEBDECrqla/r1FsrlOM+dWIOcUkeeWkvYqDTYB+voeumJceEy3zPDQOhjvhreWu+IcchAAAIsGRxwXreXqT3G09Bga9q1U9PKmtOqm5BvSSfvLyQ+N+xdJ7FYFdn/csyw7wicl4RAEAzD1+cUj6vUZ/dp8gqzOhFzCwgLs0nuzN1Nj60hPKMXfcFPQtrrqnmoC/RpdA1o/5rk5L/VChrlDPRpeg6S47JX6ffQOgRYc1J5cPj6rpqpb+TuKKQfPsicng6gdEvPtI5WNCHXNAHAOCgW11zUn1+n/LLL+VJucSc3uSs3mSmXm8C/XUhS65eCxcBa57JfwpbhOeCgVA3VAjUG3hAHk5PHNVdffOIsKpKee+Y+lWtMjGHmNObfG4sk5X8qzx0rZ+D6OcgfjOYdAuw5qTyYZX622/FkRnEwmJyXpHuImKgjjd0xWRTaAAQvZJxG7Wxhp+LXoRaBNpMGWuArR1Ljqn1gF56YAQF1pxU3qpU/12jTM4lF/Yh3pzCpPTU2Z665WDhqhLyqhIIydTaauW9o+qDW8Wx2cTVJeTcIlIneeACdUYeswAAAEsO568NOe22RBdEp/TxRUMA/lMhq5G7X0A3PTC7WtRlB5WVR5SBqcS1peTLE5hUA9/N9xQmCub0Juf0hoBErapS3jqiLNosz+5N3lxGTshJcO+10ccsAMCaa/LXhpxlGAg7hoFQL/y1IWuesQMhm0ITJPCtIpeagJaXT4QVR5RXDyiNIbipH7l1Nl1kx+E/47HQWhuRrA/CiiPKbRtlSYFb+pM39CUzElE/VFkN1IUM3yLMNXmPBRJdCv3CQKgXgdqQcRcqhdkKzN6TwTgHwgMu9fn9ysojyuRc8k8XUJfkJcki7h4u2wz3DCbvGUxublBf/kHp+654eSF51yDywsy4/nf9dTyXyho9m7w1z1S3KXnyXUQdBkK98J/iC2cY+64TAOwFZt/JYMbQlDi8lwqwrlp9Zo+8u0X9dX9yz3xc/5CcxmYRY7OoVp567ZDyiy/kTDPcM5icX0TGZ3m+70TAXmj4uVXWHC7QwKuKijeJHcJAqAuKqPBu0aBZRtuyFVpqvmyK9bsICqyoVP6yR2FIuGcwuaqY5Aw8xwhFJJWD+4aQ9wwmPz6hPLNXeWCr8pvB5C39yFinNvWdDCZBKnmSJTkHHWwULNk4YN4BowZCjwhP71ZpUnawBAA4WLDQkMJCKktkmCDdRGSbjbTNe6CWt2SySXCzZisw+04GQe3spkCRCkjwygFl6R5lgBOeHUtdnGf4Twx1CknA7N7k7N7k9ib1yd3KYzvF2weSiwZRaTG7vPtOBnPGpMXq1ePIkmsK1IaMFQibQtAYUptD0MKrHhECErTykMbBr/pHuTfAqIGQBEhlCSCJVkEFgGNeCEjgEaGFV5pD0MyrLTxkm4neNuhtI/rYoTSF6O8kBjh1ukuAvzZkMfhMGQ1jpWgzFWwSot669Uvwtx+Up3fL43PIjy6lRmZgCOzRyjOIdy6iKj3kE7uUsnfFWweQ9w6h0qN9kVckNdgoJEfd1CaOZgzX40QEvwQ/uNRDbvWwG4561eNe9aQfTgVUGw2ZZiKdgzSOSGHBSoOTBSIGOTCMGghtDNw/BOizjxIICtQF1BM+OO5Tj3rh81PqX39QDrhUJ0sMSQNty5jRWYRORpX8tYZfOxGmNQqjGAhDMrz8g/Ln3fKEbPLfl9FD0nTxL0N6UJpC/H0i9T/Dycd3Kf3eFe8cRN4zmIriza6/JmjO4kg6Gb5y1lxT/XetiS7Fac08fNugbm9Sv29WdzWrdUG1n4Po6yD6psDUXKKojCy0QZ4lfsnujRoIz4slodBGFNpgQptOOhWgyqvualG/b1ZfO6TcWqGaKGJCDjExh5iUQwxKTdjX3V8TShtowI3pO2IrMPtOBjJHRuHGU1bhzcPKQzuU4enEmhn0MAyBqCNFduLlCdSDw8hHdihl74kPDKXuGBidYWPfyeTZa9qab/J/lMhlvke96te16td16sZ6tT6oXpBJjMogriom/nwBWZyS4C2skzYQdogAKLITRXZidu/TR4541Ip69eta9andSkhWL84jZxYQ03qRUe9jOQdVUX3VyTAgr7EVmE+s83b/ddZVq7/9VnZysPIiamwWhkB0Hn3sxOuTqf0u8nffKc/tk/40iry6tLuj7t4TAUeJNTrlSzRTOqtIquAWWUf8Vjf5RPj8lLLmpPp5jRqU1am55MQc4r4h5ACdbVTZswLhmUpSiJIU4vq+AABHvernNeo7R9XbKsTh6cTs3uTcIqJP7BdlB2p5LpUxdHK1tmwFJn9NsDsTtfe71Pu2yEe98OSF5Ozexl6/heJsoJP46FKqok6991v5uX3KM2Oocdldr8K+k8H8KRlRLF5i2QvN3hPB9CExD4R1QfhXlfLRcWVTvTo6i7isgFw0iExgl9t59fRA2Faxnfh1f+LX/YGXqS9OqR9VKWNXyYU2YmEx+YtiIj9mG9R5TwTshTrdsaELaBPFOphgPd+FhP0uAR7eIa+oVH4/nLp9IMn0jCAo84rgFgWvJPok0SdLfkkMyHJQkUKyzCtySFYkVeYVVVZlXgk/S1XVtrMGSJYkaYJkSJIhaDNFsSRlImkzRVso2kIxVppNoRkbzTropLnlOocJOcS3s+kVlcpVX8iTc4k/X0h2YTaAzCu8SzR0uu127L0t3qpA+pBYLfNt4eHdo8o7R5XdLeplBeTN/cj3Lo75+paowEDYAY6CywqIywqov46nvqpVVx5Vhv1THp5O3FBGzi8izdH+zLxVQXvvJOkX1aQUWdzHAp0KhCrAm4eVB7+Trygk9/2c0dv+A1EhBeVQkxBsEkLNAt8ihFpFwSXyrSIQwDoYxkYzNoq10bSVMqWxlImkTRRlIimOIhmC4kiCJCjTj7cGPp/PZvsxdaQWKRVBUSRVCmoRVJFCsuSXg42851hA9EqiTxI8kiKpnJPhnDTnZLl0xpTGmjNYUybHWJMqQBIA15SSc4rIx3fKw/4pLR5GLRrUuTX4nmMBWy+zvrrwusfe23zys8aov6ykwKcnldcOqV/VKjMLyHsGk9N7GWxpLwbCc6EIuCiPuCiPen4s9ckJ5bVDym82y1eVkLcNIAdHr5nvrQrkTUqP1qvpgaPU2rLPkzsu0tVXP7jU2yrkoAyrptGjkmVdhBSS/TWhQD0fOBUKNPDBel4WFHMGa8rgTOmsrcCcPtTBpTKck6G4rrR8KYls27aLvJ2niArfKvJukW8RQy1C6wFfbRMfbBIAwJLNWbI5S47JksNZ881JEBqtNPxxFHV9X/KuzfLrh5WXxneip9Rd6XP0TZIBQo2t0Oyr7tawRTvHveorB5TXDyulKcRNZeSbUxi7Edp/Z8JAGBGOgvl9yPl9yGq/uuygOmOtXGKHuwaRc4vIbk52kkJyknW/AICj1HpsVV0ky+p5GR7fJf91v/LwSOq2AaShb75lXvGeCPpOBnwnQ77qoOiTLLkma67JksOlD0uxZHHxnKRwDiRDmrM4c1b7r5zolwN1oWAD768NNe/1+GtCJEvaeplsBWZ7gdleaKENGxf7Ooi1M+j3jikLv5Bn9yYeGxXREgtXpb94Vk7sSxc/tIninEygju9mfn8V4PMa9bl98pYG9Zd9yQ0/o/s5jFx1MRB2Vi8r8dBI4vfDyY+qlGf3KfdvVe4e1K08T74TwSTrfgEAzsnQZuq8+2lsaVBv+loe4CR2zjNqmlDeJXqO+D3HAp5jgVCLYM0z2Xtb0oek9J6ZZc7gjJTcCICxUo4Sa9tJkqEWwV8d8lUHa75q9p44ydqZlCJLSrElpdhqxHSAC/qQl+aTi7+TB38g/W08dXnhuf49UlAO1vP23skzeK/Rhgm7HAi17IZP71EIAn4zmHz3ougPFSVEUvwRcUeT8PM+5M/7kNua1Kd2K4/tFG8bQN49uCuJLZJvgFDjKLW6K/1nq29BCX6/TX7nqPrsWPLnfQw2JUbyy67DnUaE0wAAIABJREFUPtchn6vSL4cUR7ElpdiafWGqNd+U4MVQ0WZKY01pbLqWQl2FQF3IcyzgqvSfWNegquAstTrLbM5+Nm0DdENwsvC38dRVxeot38grjxLPjT1rbjbP0YC9tyXJ/qEAYO9t9lYFc8Z2+on+/2Q3HOiEZ8ZQl+Qn1SdjmG+wPo3KIN65iDriIZ/crfR7V/xVf/K+IVSndk3zVgWyR6fGrIAJ4+xrbdju6nDsc0uDesNXcnkGsXs+Hc/1mt2igq862LLf23rAG6wXUootzn62vEkZlmyDNfu6jgBLrsmSa8oZlwYAoWbBddjfst97dFUt52BS+9tSB9hT+lgM0bcxOZfYNY/+n23ykA+klyaQVxR2cCvmOuxzJtcAocbe23Lq6+ZOPcUnwov7lf/dK0/MIVdNo0akG+Bf3FkYCKOgJIV4eQL1PyPIJ3Yp/d8T/2sged8QyhlJ15EK3hPB0oX5MS9i3DlKrZXvnWo3LC8o8PB2+fXDygvjqHlFBmgIqrLqOuxv3uNp2eehzVTaIHvRz3JS+iRhQ6GzTOlsTjqbMyZV+w63HvAeW1XHtwip/e3pQ1JS+9tIfW/gZ6Fh6RhqbpF649fyR8fV/x1LtZvl4a70ly5IwoppyeF4lyiFZDqC9GVBCV46oDy5S56aR67/GT3QmbRfewyEUVNgJV4cR/33UPIPO5Syd8X7h1J3DSLP/WULNPAUSxqoZylyjI1mHYy/5scMVftd6rUb5EIbsXMuk6XvzmBVUV2H/E073c17PZYsLn1oytA7i00ZxhsViwcC7L3N9t7mwulZgkdq2eup3dRyeGWNs58tY4QjbYBNz6tBJ+YQO+fS930rD/+n9MZkakLO6Qu96Jf5VtFWkISLeAiSsBeaPUcD587pKKvwxmHl4e1KeQbRExL8JuElOLF624hlk6gDLvJ325QX9kuPnTPPk+ugz1mWhN0vGmdfq+uw31ZgVgH+ul95ZIf8+AXUzf30e1kEAG9VsHGHq/F7tzmDzRju6D0zOylvU2KETaFzxqXljEuTAnLzHk/dppbKd2rSB6dkjnI4S2367EO2MfDyBOqTE+qVX8g39SMeGkHRJLgr/Ubp5u0CZ5nNddB3jkC4tlr97bdyhgnevZga0zOyG2Ilj4n+TuKfl1Cb6s+T58l1yJd1gTP+xYuP1P626i+a2LEZN30tNYZg0yy6NEWnlUrwSg3bXA1bW1UFskY5h/2m2JSG7b+uoy1U9ujU7NGpoldq/N59/ON60V+TdUFq9oVOfX6wlxcS38+lb/pamvCJ9NYUSvnBm9rfdv6nGZOzn+3g/1V3+Cstu+FxHzx14Xlm1SYZDIQxNC6b2DyLfvuIcqWW5+kCsm2eNlVW3Uf9ZVf3SmAJY8rR17b3zeqr3gvNGcg8PJLSYw+ZCq7DvrrNLa5D/vShKaUL81P6JNt0+cRi7HTepPS8Sen+U6H6ra27/veoNd+UMzYtfbBdb+2tLDN8PJ1+YZ8y8V/iyh88Yy7LSnSJYsWWb5YCMu8SOeeP46IuAR7ZIa84ovzPcOq2AXruz44JDISxRQBcXULO7k0+tlMe/qH04DBq0aDTXzLPsYAli6MtRl2kfG6SAv/zvZJitb2Y4Z8ySnd5c+SQUv9da+3GFpIhcsem9b2yV9vsZSjqrHmm4jm5RZfnNO92137TfPTD2txxaTljUxmbji5BBMBdg8hxAX5PJbtsJ/niODBEnsxOI8DZ1+o66NPmq6sA/zisLP5OnlVI7pvPdGrSe9LQ0bcwiVlp+NMo6sYy8q5N8uuHlL+OpybmEK0Hfc6y5Ox+qfKpv/hCTuPghZlO7+YmAB0FwlCLUPtNc/02V2qZre+V2ASMK5ImMkc6M0c6/adCtRtbtj9xOH1ISt6kdF3tSp1y1DNlqrOChFEfSe9cTCXlLpjOfrbWA77s0al7W9X/2ijzyZXdsAs6FwiPHTtGkmTv3r07/K0oilVVVWazOT//x2nHfr9fEATtMUEQTmfSDomdV2kKsWYG/cEx5eoN8rR84r8O+vrOTqoETpqPqpTbKuT7h1L3DiFVybb1g2rRKzH2xN9y+WpCNRsaXQd92WPSRtxX2rZfCMWZNc9UuiCv6GfZdZtb9r1SZc0z9Zqa4ShN/MQxVVZb9nqGT8/6u5NacUS59FPpkXLqvwYkW2+Bs8x27OO6322V/35I+UM59ev+xs5u2H2RXp78fv+sWbOOHTsmy/KgQYP++c9/mkw/uYlbtmzZAw884HQ6g8Fgbm7uBx98oMXLW265ZdWqVRzHAYDT6Tx69GjU/wZjmd+HvLQX+ceNQkutsFo2XZ3o8kQRL8N/b5U/PqGumkZfmEkAAEETaQPsTXs6kYA7FjzHAic/bwzUhvImp5cuyO9akmsUdbSF6nVxZv6UjIZtrsr3T9EWquCSzLQB9gTOL3Ud9pszOe0m6eoS8oIM4qoN8oZT6qsTI0pPahRf+2mXQgVqQrvnW3P0vZYpPiK9Irz44ouyLB86dKiysrK1tXXZsmXtTujbt+/3339fWVl54sSJkpKSBx54IPyrp556qqWlpaWlBaOgJoWB/3YGbX0sT+5VL1srHfeqiS5RFFR61PEfSzUB2DH3dBTUZAx3NO10J6pUnmOBvX87fmhFdfqQlFG/L8ufnIFRUG8IisgenVr+QN/8yRlVa+p3PnOkZZ8XElQnmna5M4Y7wj/2dRCbrqCzzTDyQ+m7xmSop808/PJL+dYKOWOA/T57AKOgJtKLwsqVK2+55RaaphmGuemmm1auXNnuhEmTJhUUFAAARVFTpkw5efJk+FeyLLe2tkarxMmhaZe7pDxl2xx6ci55wb+k5/YpipFr2cojyviPpRvLyPcvbp9Sx9nPFqgNhVqEOBfJWxXY+9LxQyuqM0c5yhf3zRmTiulgdI2AjGEpI+4tLZiWeWJt/c5njrT+4I1zEWReadnrOZ1b9T84Cp4fR/1lNHn5v6Vn9hq6msLbR5QhH4gZJtgznx45JqV5tyfRJdKLSLtGjx8/XlJSoj0uKSmpqqo625mSJL355ptXXHFF+Mjvfve7P/zhD4qiPProo7fffvvZnqgoyq5du1JSfvItHD58OEUl27xKOaS4DvpLF+bTJDw4jJxXRPy6Ql55RPn7JMpwSYwCEty9Wf6mTl03gx7eURJCkiayx6Sd+rq5eE5unIpUGzr+ab3/VKjg0qzsC516m6aPzoWA9MEp6YNSmvZ4jq2qO/l5Y9HlOXGb0FS/tdXR19bh4PHcInJEOvGLDfIXp5Tlk2jDTa2s9qu3b1SO+9SPLv1Ph02JVfCIgXrekm2UhL8xRKhqRLc4HMdt2bJlxIgRALBly5bLLrvsbI28RYsWbd++fcOGDSzLAsDx48cLCwtJkvz6669nzpz56aefTpo0qcMn2u32Pn360PRPYvPq1avt9g4yIASDQYZh2p1sFC07fe4fAn1+8eNCJRVg+RHq0d3UrWXKvQMkfadpPM3v9x8TbTduYkamq0vLJSt91i+S6JEOvHhq4D0xX6IguKW69S5PZTB7oiP9AjtJJ3MIbLdDfRJSoWWXr+4LlymbyZuWZsqM7eQmVVF/eLamaEGWpddZBwNFBR7dQ71XRb88RrzA7jfEJajtteW+AVLbBYKn1rUSDJF7kVEnMPr9fqv1/BOsLBYLSZ7nyhPpfzE7Ozsc+Zqbm7Ozszs8bfHixRs3bly/fr0WBQGgqKhIezBp0qSZM2euX7/+bIGQZdkvv/wyLS2iWRUURRniW9ih43ubcsemt7uKLRoG80rV2zcqkz+j/z6RGq3vzEYqwEuHqCf300vHUNeWkgDnvKm0Qfogn3eP0GtqRozKI4eU6i8aaze35I5N67ewsCcsClRVNckDIYBtgq3XmOzajc1HltenD0kpnJHFxmz6cdNuj8nJZvU/z/Vn6Xi4rEi94Svi2mLqkZGUidX1JeiQW/3VN7KowFdXdNDblDeaOvDGyb5X9NJn9rtIRKsKRHq9GDFixObNm7XHmzdvHjly5Jnn/L//9/8+/fTTf//732dbI1FfX9+u57MHEjyS72Sww0R/vazEqmnU74eTcz+XfrNF9onxL11E6oJw+Tpp5XFy8yz62tKIvkL5kzNqv2lW5eiPsKiKWre5ZfsThwSPNPK3pb1nZveEKNhzkDSRPzmjfHFfykR+/+Th6vWNihSTcbpTXzblT47oRu3SfOL7ufSuFpj0qXrIrdNBQ1GBP+1Uxn8sze9DVlzR8cYRtl5mgiK8JwLxL57eRHo7s2jRoiuvvHLIkCHy/2/vvqPjqu59gZ/pvXeV0Yy6ZFmyZNmy5F6wIdiG5QWEGGd5hefYECA3tBvyHO4iD15IbkIeSagmF4xfQvK4sROwwQVi3OWiXqw2aiNpuqb3mXP2+2NydY1sZEkezYw0v89fo6Mtzba8Z39P2QXH33rrrWPHjsWOl5SUfPjhh8uXL3/zzTdfffXVl1566fDhwxiGcbncnTt3RiKR559/ftOmTRwO58iRI21tbX/605/m6p8yT9haXOIy3hRLGD2cR96cRX7uCr7ocPSNultvlpZEfx8mfnAR31NEfqY2IuRP9+kCJ5PJktEtDc74br7o7PMO/t1E5VAWfV/DyZxvz23AtFFZFO02papOPHTM1PSLPs1WhbRCEMfrGJfOF/ZExWVTbchwIzkLO7KeeL+fsupo9OVqyt7i1HoQXW9B+y7gag7WeD9VzZ2qarIqgbXJxctJ92UlphuEGzdufPfddw8cOEAikQ4dOlRXVxc7vm7dutj1H5/Pf/TRR/V6vV6vxzBMIpFgGEahUNhs9jvvvBOJRAoLCxsaGrKyFuzSmtNkuebQbLvNPHoxA3t/DeUrI3r8Av5BL/pdLTmLk/wPmiOE/bAev2pFhzdRa+Ukr3dmP67Zruo8MCRexIvLqlpBW3jwU5PPGNRuV0oWp/tthjTBlNCLd6vdA76Bv5sM58dz71dNbPJ1J4go0v3VoL1POaNxVSQMe6wY25RF3X0WPzJIvLeaMnXkJIYjhP3kGn5sBP2mhvxQ7u3PoeVVwtbfD2i2KRf2A/Xbmu5gmQSQSCR9fX3TfEY4TwfLOHu9A38zVv1rwTRPZkM49ss24ved+PPllB+VJXOv08ODxL9cJh7Ukn5eTWFRMWxWgzWGjplCzkjRruw7qUk0gI98YbU0ODPXSTPWSNL2A+zxeG45jiwtIMx81TF8wiIq5Obce6dbZQ0fNwcsoeLd6hn91EQXhCPsV23E/+nAX6qi7EveEi0Ewg72Efuv4Q9oya9Uz2D6//U/DIsX8ZS1yVzyYnbiOF5svgah1+Eb+WScTCVTGGQylUSmk2lsCpVDoXGpDCGNLqDROKk46aL9rUHFcpG8embDtAY86Ef1RK8L/bomCXuj6L3oqXqi340OrPraZlKzaIVEhGj6d13eDpWoZDY9OMKRqd6u/8IqKePn3CNPqfWaEy+tgxDDMAzDQ8ToP6ymertqlSRzvZQyq/NEvzHY/vZQ5XP5M03TSefiXU609wJOIOzdVZQyUaI/pJfM6OnLOJWMvVFHqbzVLKYpuAf9vR+NLv1JQQpONIoG8LArGnJGIp5oxB+N+nA8TFCYlJy75Vhcg3C+diUkCiYu55IxMh4iiCgiwkTYE/WbQhFfNOSMhJwRhCOmlM6SMThKBlvJ5GQxk74RmmfYH7JHZFWC2xf9ulwe6dPNlOMj6Lkr+Osd2K9qZtzWZycQxX7dTvyuE/9RGeU/N86un/kaMo1c8O3Mnj+OlD2mZStnMnsJYbZW1/BxC1NCK3tMk1JrNINkoTDIOd9SKGvFQ5+ZG1/tU2+WKZbPbNmEsDva9eGIZmsctl8uEZLObaW+101s/Dz6cC75pSqKKCHT8/rd6H82EPVm9L+XkXd98x7gU+Br2Qwhzdbski1N8jyKsDvqHQ34jUG/OeS3hIK2MCIQQ0hjCGk0LpXKptA4VBqXOhdnwPP1ivC2t0ajATw4Hg5YQn5TyGcM+saCeIjgadi8HJYgl8PLYSV+x63r/zEsKuHdyaqbUQJ7r4d4pZlYqSD9bCm5ZM5m3+MI+799xEtNxHIZ6dc15Fs+/Jj16Zi1yTV41LT4cQ1LPo2uAmHjnW79SQuJQtJ8S7FQ9+uYBbgivJF3JDD0uTloC6s3y2VVgunEYcQTbX9rUF4tzNoom8U7flMXNB7CXmzAjwwR/1pO+UEpmTlnd6ZGfOiVZuLIEPEviyjPLCaz7yAdHN3ewU9NVc/nJ3geBSKQVx9wDfg9w37PkB8hjJPJ5GYwWQoGW8FkSmhTZx7cGp3NM8KIJ+rRB9yDfveAz2cMcrNZoiKusJjHzWAm4L/fow90vT9c/dOiO3+m5Y9ib1wnXmvH1yjJP64gx3fzlCiB/WWA+HkLIWdhry6j1H7zdMY7aYWWBufw5+b8BzOmuEeKhwlbs2vsjI1MJ2ffJZMs4s/f2U5zAYLwZq5+38gX1oAllLFWolgmmmKzT8+wv+8vY7JKQfbmWW7AO3UXdN2JftpANFjRjyvIjxaSWXG9hul1oV+1EUeGiL3F5OfLKeJ4XHq2/KY/a4P0xnVW507AEnJ0ex09XvegjyWh8/M4fA2bl8NmiGa2ZgIE4Z0OlsHDhLvf5+jxOrq9RJgQl/LEZXxhAWeOlqMkoqjlNzr1Znkc25kviv1HD/FaO6HmYo8Vkx/Qkhl3du5pC2KH+ojfXye0XOwnSyh3Zd7mT3GHrdDZ6+0/YmSKaZnrZbwc1sSN17A76h702zvd9k4PP5eduTYldudJQRCE38Q7GjCcHbdf9wjyOZJyPl/LnngsQkSRdzRgvGh39/s0WxWyqtnfDJxOF3TNil5tJS6ZicdKyP+jiJx9Z2O/cYR9PkK800U02NATpZQnSsmS+N199QwHuj4Yrnwuf64evSPMPeQf73DbO9xEFImKeaIirqCAQ2XNvtuCIIznqNGANWzvcNva3QFrSLKIL6sUCAo48X1uPHTMHBwPF+++o9GStxQlsKN64p0uotGG7teQH8olr1WSZpSI3gh2fJT46yA6NUpszyE/UUq+ce+IqX7wjlshwpHxkt3W4vKNBekCGiIQESYIAvE1bGEhV1YpSPPhMFODIJwaHiRsrS5Ht9c96CNwRGVSSBRS2BlhyRniMl7mulkOrpkw/S6o24nevE78uZ+okZMe1JK355BndA2HI+yyBX08QPznIJHDJT1WQn5IG+dLzJg56aYQ5h7225pdtlYXjUuVlPPFi/jcOE35hSCck+kTIWdkvM1tbXYG7RHZEoG8WhiXWUqeIX/XQf0cnmphGIZhYz708SD66yDRbkcrFaRVSnKlhLRYjGWwJ1/l4ggb8qAeF3bFQpw3oUYbqlOQdmjID+aShTMZThTHVohwFBwPkygkMo1M51HhFuh0QBBOX8QbxUMEwhFdSIvDoC8Mw2beBQWi2N+GiSND6ItRokxMWqUgrZCTSkWkXB5p0nAFhGHmANZhRy12dMmMzhiJHC5ph4b8cC6pQDCHn43Yjavsu+SyyjjcuApYQpYGp7XJRaaRZFVC6RI+Sxbn4UMQhHM7jzA4HrY0OK1NThKZJK8WypYKZ72buc8Q7DwwlP9Q5i3XVJsLzjD2lYG4bEHN46jTgdmCSMYiMcgYl4aFcMwZRs4wpmSRCgVYtZS0SklepSTxZ/WPW/iLPqc2CMLkmnUXFIhil63onBFdsxI9LmzEh3g0TMQgcaiYO4xFCMwcQGIGViQkLZGQlklJmzLJikTtGugdDXQeGC7alTXrUWlRP25tclkaHCFnRLZUKF8q5GTM1RhvCMIETah3D/kt15y2Nhc3kyVfJpSU8We0rat3LHj9wFDujgxpRdKWPokSmDmAQgTmjWAsCsajk0R07A6fJsZAECYXBGFyxasLihCYM4w5QyiAY1waRiVhSjYpiUtnuAf9XR/oC3dmiYpn8OlGOLJf91ganC6dT1TCk1cLhYVxfsB0MwjChK4sQ0SRvdNtueZ0D/iFxVzpEoGomDv1DRZEIOMF+8iX1vwHMxbqAmAQhMkFQZhc83Rxq+nwDPuvv69X1oiyNsmm7uiIKHLpfLZW13iHm6NiyquF0nJBwla9hyBMTiuM+nFbm9vW4vLo/XwNW1DA5WvYnEzmjW0lOB52dHuNF8cZQlru/appTZWbnyAIkwuCMLkWcBBiGBZ2RQaPmt0DvozVEmEx98YlLIgo8puC7kG/q9/n6vWxVQxJOV+2REAXzO2GkTeDIExyK8SDhKPX6+73uYf8flOIQifReFQigqIBnEwji4r+edWYlLolDARhckEQJtfCDsIYz5Df0uh09HgjniiNSyXTSFEfHvHjLBmDr2Xzc9miIm4Sh3bDEmtJRmGSpeV8afk/73lGvNGoHydRSBQGGUb8AwAWBp6GzdOwMQyLBvCoH0c4ojApC3JcN/TacTBHy98BAEAqoLIodzLzPfWl1qavAAAAQIJBEAIAAEhrEIQAAADSGgQhAACAtAZBCAAAIK3BWEcAkolAyBfx+SOBEB4ORAOBaBAncG/YhzDkjwRwhBOI8EX8E+WD0VAEj2AYFg6H6XQ6hmEUMoVN++/FKNlUFoVMIWEkLp2DYRiXxiGTyBw6m0lhMKgMDo3NprEopIU8AhCAmYIgBGBO+CL+8YDdGXQ5gi5nyOUMut1htzvk9YQ9nrDXE/b5wj5vxB+Khrh0DovKZFAZbCqLRWVSyBQunUPCSCwai0qikEgkLu2/d2RkUhlMBhfDsBAWYjAYGIZFiagn5J0oYPJaCEQgDHnDPgzDPGEvgQh/JBDEQ8FoyB/x+yJ+GpnKoXO4NDaXzuHSuXw6j8/gChh8Pp0vZPLFTKGAKZAwRXwGTNgHaQGCEIDZcwRdZp/F7LNa/DaLz2oL2K3+8fGAfTxgp5KpEpZYwOALmQIRUyBkCJQcRaEoj8fg8ehcLp3Do3E4dA6LOsu1+e9kZZlYInojPk/Y6wl73SGvJ+x1hdzD7pFWi8sRdLpCbnvAGYgGxSyhjC2RsMQytkTGlsrZUjlbpuTKJSwRaeFNqwbpCoIQgNsjEGHyWUY9hjGPccxjHPOYjF6T0WdmUBhKjlzOkSk4UjlbViwtlLMlEpZYwhIzKDPZ3TGxmFQGk8oQs0RTF4vgEXvQYfXbbYFxm99u8Vu7bL0Wv83ks3hCXiVXruIqMrmqDJ4yi5eRxcvI4CqpZLjpCuYfCEIAJosS+IhnbMipH3Tph10jevfoqMcoYgqzeKosXkYWT1UhL8vkqVRcxayv5+YFGoWm4MgVHPnN3wrjYaPXbPCaDV7TqMfQYGwZcY/ZAnYFR6YRZKv5WRqBWitUa/jZNEqi12IGYKYgCAHAvGFfr71f5xjocwwOOIdH3KMKjkwrzNEI1GvVdWp+VjY/M5Wv8BKPTqHnCLJzBNk3HowQUYPHOOQa0btHL49d+/P1w2Meo5KryBXmFIhyC0S5BeJcEVOYrDoD8E0gCEE6CkaDvfb+rvG+7vG+HrvOGXTlibT5Im2lYvEDxds0AjXE3izQyNRJ6Rgl8GH3yIBzSGcf/EvX33rt/Uwqs0icVyIpLJYUlEgLOTR2EisMQAxswwRmad5tw2Tx29osnR3W7g5r14jHkCfUFEsKiiUFxZL8LF4mmTTPhn7M022YDF5Tr72/e7yvy9bb6+iXs6WLpMWL5aVl0uJsfmayazcD0AUlHWzDBMC0GLymFnNHi6Wj1dwRwsNlspJyWclm7boCUS48u0qKDK4yg6tcp16JYRiO8EHncLu1q8HY8kHbRxEiWiFfVCFfVKko13z9pisAcwquCMEspewVoTPkajC2NJpaG01tOBGtVJQvUSwql5ep59UFx23N0yvCKZh91lZLR4u5o9ncHogGq5TlS5UVy1SVcrY02VW7BeiCkg52qIdWmHwpFYQ4wjut3VcMjVeNzUaveYli8VJl+VJlhZqfleyqzZWFF4Q3Mvssjaa2RlNro6mFz+AvV1XWZCytkC+ip8yzW+iCkg6CEFph8qVCELrDnstjjfVj1xqMLSquYnlGVU3G0lJpYTosIbawg3ACwlCffeCqsfmKobHfMbhEUVabuaw2c5mUNa2OYu5AF5R0EITQCpMviUFo9lnOj1y5MHq5195fpayozaxekVEtud308AUmTYLwRp6w96qh6dLYtavGpkyualVWzersFTlJepoIXVDSQRBiRqf56TM/JZFIEysxMqgMGpnKobFZNBaHxuLRuVw6V8jgCxh8CUskZoqETOG8GxmYyhIfhGMe4xn9xbP6S2afdWXW8lXZK5YqK9J2nkMaBuEEHOGt5s7zo5cvjFxmUpnr1CvXqmvzRbmJrAMEYdw5gi5n0DkecDhCTlcwtiSv1x8NBCIBT9iLIyIQCWAYxmNwX9vwvzAIQgzDvH6vI+KiUCgIIW/Eh2FYMBqKElFv2BeIBn0Rvzf8z+UTnUH3eMDuCLk8IY+QKVRwpAq2TMlVxFaHyuSpFBzZHP/LFqaEBaHJZ/lq+MLp4fPjAfvq7Np16pUV8kVkUrrvIJbOQTgBYah7vO+s/tJZ/SUKmbIhZ9X6nNVagToBbw1BOGuOoHPEPTbmMRp9ZqPXbPJZrX6bzT/OoXFELKGIKRAxhUKGgEfn8uhcDo3ForF4dC6FRGbRWBiG8ejcDK4SgyDEZtUKcYSPBxxmn9Xss5i8FoPXNOYxjnoM3ohfzc/M4WdrhWqtUJ0n1EI0TsdcB6Er5D49fP7LobOjbuNadd2GnFXl8jK4pp8AQThJz7ju9PD508PneXTuJu3aTZq1czrcFIJwmhxBp84xOOAcHnTph116vXuMSqJm8lTZ/AwVV6HiKJRchYwtkbEkM53RBEEYz1bojwRG3GODLv2QS9/vHOp3DEbwaIE4t1CcVywpKJEUQi7e0hwFYRgPXxq7dnLgq3br9drM6o2aNdXKSljK+WYS2fH/AAAT10lEQVQQhLdEINRuvf6PoXNn9Be1wpwt2vVr1XVzsX4NBOE3sQccsTWbeu39vY7+KB4tEOdqhWqtQK0R5KgFmXx6fNotBOHctkJH0NVr1/1z/YvxXoRhi6RFZbKSxbKSInEBdMoxcQ/C67ae4wP/OKO/WCjO26LdsCZ7BXNBL2l9hyAIpxYholcMjScGTjeb22ozlt2Tt7FSUR7HOwoQhBMIhAacg62W65227g5rdyASKJEWFkvyC8X5heK8ubsuhyBMaCu0+G3t1q5Oa3e79fqox1AiKVyiKKtUlKfJMP1vEq9W6Ai6Tg1+9Xn/Fzgi7s7duEW7TpaSE6hTDQThNLlDni+Hzh0f+NIV8tyTu/FbeRtvuZ/GTKV5EBII6RwDLeb2ZnN7m/W6lCVeLCstk5UskhYlbKk8CMKktUJfxN9mud5sbms2t495jBXyRdWqyuWqyvm1TGJc3GErJBC6Zmw6pjvVbG5flb3i3rxNi2WlcazeggdBOFP9zqHPdF/8Y+hcvli7NW/zquwVNPLsO5D0DEKL33bV0HTN2NxkbhMxhVWK8iWKsiWKMiFDkPjKQBCmRCt0hdxNprYGU8s1YzOZRK7JWLoio7pKWZ4mA/pn3Qqtftvn/V9+1v+FmCnamn/Xhpw1bBor7tVb8CAIZyeCR86NXv5Md2rAObRZu35r/pbZrb2XCl1QYkQJvN16/fJYw2VDgzPoXqaqrFYtqVZWSNmS5FYMgjDlWuGgS395rOGyobHXrlsiL1uZtbwuc/ltdwCf12baCglEXDE0fdp3osPatVGzZmv+5nyRdu6qt+BBEN4hg9d0THfqxMA/snmZ2wvuXpNdO6NRi6nWBcWdN+y7bGi4MHqlwdiSyVPVZlbXZiwrEOelzshtCMLUbYWesPeqsfni6JVrhuZsfubq7BVrsmszeapk1yv+pt8KbQH7Z7ovPus/JWGJtuXfvSFnNZPKmOvqLXgQhHERJfALo5c/7Tsx4By6O3fj1vzNWbyM6fxgynZBd8jmH78weuXcSH33eF+lYnFd1vLajOrUPKeHIJwHrTBK4C3m9nMj9RdGLwsY/LXqujXZdbnCnGTXK25u2wpjTwGP6k62mjvX56zaXrAlwWt/LGwQhPE15jEe1Z08MXA6V5izLX/LbZ8gpn4XNCMGr+ncSP1Z/aVRj6E2c9nqrBXLVJUpfsIKQTifWiGB0HVbz9mRS+f0l+gU2lr1ynXqlQvgruAUrdAWsH/e/8Vnui8ETP72/C0bNWtZMBEi3iAI50KEiJ4fqT+qOzno1N+du2Fb/pZvup0zj7qgKYx5jF/pL57VX7T6bauza9dm1y1RLJ4vM8QgCOdlK0QY6hnXndVfOqO/SCaR1+WsXK9eOX8vkm5uhQQiLhsaP9OdarNcX5ezclv+lkJxXrKqt+BBEM6pUY/hmO7UiYHTGkH21vzNa7JrJ+0ANR+7oAmjHsMZ/aUzwxfsQcea7Lp16rryebhsIQTh/G6FGIb12vvP6C+e0V/EMGxtdt069coiSX6yKzUzN7ZCvXvseP+XJwdPKzmKrfmbN+Ssgrnwcw2CMAEiRPTi6JXPdF902/s25qy5J29jkfifn9P52AXp3aOxZeudQdcadd169coyWWnqDH6ZKQjCedkKb6nX3n9Wf/GM/lKUiK5R163NriuVFs2Lpun1egkaOj18/tTgGZPPfJdm3bfyNiVrT5w0BEGYSGaf9cTAP04MnGZQ6FtyN2zSrOGSOPOlC9I5Bs6NXD6nv+SN+Ndk165T183r/JsAQbhwgnDCgHP4rP7S+ZF6R8i1KqtmddaKSsXima5Cmxj+SODS2NWTutOd9p4VmdVbtOuXqSrn3X2V+Q6CMPEQhjqsXScGTp8bqdfw1Zs0a9blrBQw+Mmu1y0QiGi3dl0YvXJ+pJ6EkdZk165V15VIC0nYvM+/CRCECzAIJ8SGb10avapzDFYpy1dkLK3JWJoKq465w55Lo9cujF5uMrVVyMtqFdV3FayDUTDJAkGYRBE8cm6o/oLh8lVTc7GkYE127crM5UmfYI5hmCvkvmpsujzWeM3YrOTKV2YuX5W9Ik+oSXa95gQE4UIOwgnukOeKsTHWpsVMYbVqSZWyvFy2iEvnJKwOBCL67APXjM2XDY0DzqFq1ZJVWSvqMpdx6Zwk7lAPMAjCZIt1QTiJuGJovDB6pX7smpIjr8lYulxVWSotTuTAy2A01GnrbjK1XjO1jHmMlYryFRlVKzKqUyGY5xQEYVoE4QQCoT5H/zVjS4u5/bqtR8VVlstLSqVFxeKCLH5G3O91hPFwn2Ogw9rdZulss1yXsETVqiU1GUsr5ItuHDgHQZhcEITJNakLIhDRaeu+PNZwzdQy6jYskhUvlpWUy0oLxflzsYKg2Wftses6rF0d1u4B51C+KLdKWV6trFgkK06fnQAgCNMrCG8UJfAee1+ntbvT1tM93ucJe/NEWq1ArRGos3iqTJ5KzpbO6MliGA8bvGaD1zjsGh1y6XWOwRH3WI4gu0xWUiYrqZSXfdOiEhCEyQVBmFxTdEGesLfV0tluud5u7ep3DsnZ0gKRNleoUQuyMnmqDK5yRg8UcITb/OOjHuOYxzjkGol9SGlkaqE4v1RauFheWiwuSPGZ73MEgjB9g3ASd9jT7xgadOqHXPoxr3HMY7T5x3kMnogh4DN4AgafTWMxKIzYB49AhC/iJxDhjfi8YZ896HQGnd6wX8WVq7hKjSA7R5CdL9TmCnOmE6UQhMkFQZhc0+yCcIQPuUb6HUNDLv2wa2TMazJ6TTQyTcISiZhCLp3LpXOoJAqLxqKSKBiGBfFQGA97wz5fxO8IupxBpyvkFjGFGTxVFk+Vw8/SCnPyhJrUXPMsweLYBaV7kMx3fDqvUrG4UrF44gjCkD3gdASdrpDbE/b6Iv4wHg5EgySMRCFTsqgZJIzEpXN4dK6IKRQxhSJmEvZPASBNUEiUPKFm0nAVd9gT+5B6wl5v2IcjPBANEoggECGliBkUOpfOZdNYAgZfwhSJWSIYkj3XIAgXGhJGkrBEEjhhBCBV8ek8Pp2ngUm3KQNONAAAAKQ1CEIAAABpDYIQAABAWoMgBAAAkNZmFoSRSGTqAgRB4Dg+ix8EAAAAkmK6QWgymTZs2CCRSGQy2YcffnhzAYTQM888IxQKRSLRnj17otFo7HhnZ2dFRYVEIsnJyfnyyy/jVnEAAAAgHqYbhM8++6xGo3E4HCdPnnzyySeHh4cnFfj444+PHj06MDBgMBiam5sPHDgQO7579+6HH37Y5XL99re/ffjhhwOBQDyrDwAAANyZaQWh1+s9fPjwCy+8QKFQqqqq7rrrrj/+8Y+Tyhw8eHDv3r1SqZTL5T711FMHDx7EMKyjo6Orq+tHP/oRiUS6//77FQrFsWPH4v5vAAAAAGZtWhPq9Xo9juMFBQWxL0tLS3U63aQyOp3uySefnFRAp9NptVoWizVxvK+vb4o3cjqdpK9vFykUCkm33EASIcLvJdJ+ibUkQgEfAWOtkgcFfARl4ewtN++gYJCIUKELSgIyicyM8w480/pfdDqdHA5nIpB4PJ7dbr+5zMSybzwez+l0EgQR+8GJMjwez+FwfNO7hEKhqqqqSbHX0dEhENxiDTC/1RQ58OPpVB7MEYSQe/5vcj1/IYQ88PcH6YcsUvB+8O8Yhvl8vumUZ7PZZPJtztmnFYRSqdTj8RAEEft1DodDoVDcXMblcsVeO51OmUxGJpOlUqnb7Z4o43A4SkpKvuldGAzG9BfdplAyaS//P1h0O4lg0e3kgkW3kwvW/U8F8eqCpnVvKzs7m81mt7W1xb5saWkpLi6eVKa0tLS5uXlSgeLi4oGBgYmAbGlpmSIIAQAAgMSbVhCyWKzvfve7L774oslkOnLkSH19/a5duzAM6+rquvvuu2Nl9u7d+84777S0tPT19b322mt79+7FMCw/P3/16tX79++32Wyvv/46QRAT5QEAAIBUMN3r+l/+8pdPP/308uXLlUrl4cOH5XI5hmEEQYRCoViBLVu27N+//zvf+U40Gn300Ud37twZO37o0KEnnniiqqoqPz//2LFjcCcBAABASoGNecEswTPC5IJnhMkFXVDSxbELmq/j3994443PP/882bVIX01NTS+++GKya5HWdu7c6fV6k12L9AVdUHI1Nzf/9Kc/jddvm6+nM93d3XA6nETj4+MTg6dAUly4cAGW8E2inp4euCOSRHa7vbW1NV6/bb5eEQIAAABxAUEIAAAgraXQYJm77767t7f31guq3cTj8dBoNCaTOde1ArcUiUT8fv8tF/0BiWGz2SQSyTQ/LyDuoAtKrkgk4vP5hELhbUseO3bstvPXUygIXS7X+Ph4smsBAABg4cjKyqLT6VOXSaEgBAAAABIPnhECAABIaxCEAAAA0hoEIQAAgLQGQQgAACCtzYMg/PDDDzMyMng83vbt22/eEBjDsM7OzmXLlrFYrLKysvr6+sTXcAGzWCwPPfSQXC5nMBg1NTWXLl26ucxrr72Wd4OJXbdAXDzzzDMTf9uKiopblnn99dflcjmfz9+5c6ff709wDRe2Sc27sLAwHA5PKlNbWztR4PHHH09KPReYl19+efPmzXl5eSdPnrzx+CuvvCKVSgUCwZ49e27+j8Aw7PTp04WFhWw2u66uTqfTTfPtUj0I+/v7n3zyyb///e82m43FYr3wwgs3l3nkkUd27Njh8/meffbZBx98ENadiiOfz7dy5crW1lav17tjx45t27YFg8FJZRwOx/r167/4L7D0XXxZrdZdu3bF/raffPLJzQWuXr36yiuvnDt3zmQyWa3WV199NfGVXMC+973vTbTtTZs2FRYW3jwWf3h4+K233oqVeemll5JRzYWGRCLt2bMHIXTjNvSnTp16++23GxoaRkZGOjs7f/e73036Kb/f/9BDD/385z/3eDwbNmx49NFHp/t+KLX927/92wMPPBB73dLSwuFwgsHgjQUaGxt5PN7EQY1Gc+zYsUTXMj34/X4SidTT0zPp+P79+5999tmkVCkd7Nq167e//e0UBfbt2/fEE0/EXp86dSojIyMh9Uo7OI6r1eojR47c/C2VStXR0ZH4Ki14paWlhw8fnvjy29/+9v79+2OvP/7449LS0knlP/roo7Kysthrt9vNYDB6e3un80apfkXY19dXVlYWe71o0SK/3z82NjapQGFhIYPBmCjT29ub6FqmhxMnTigUCq1We/O3Dh06pFKpampq/vznPye+YgveL37xC5VKtWbNmlOnTt383Rs/I2VlZQaDAXalmAsnT570+/333nvvLb9777335uTkPPDAA9O/HQdmalJT1+l06Ovz4Pv6+hYvXhx7zePx1Gp1X1/fdH5zqgehw+GYWOKdSqWy2exJjwkdDgeHw5n4ks/n3/I5IrhDAwMDjz/++DvvvEOj0SZ96/777//qq69aW1ufe+6573//+8ePH09KDReqffv2nT17trm5+ZFHHrnvvvtuXnH/xs9I7L60w+FIdC3TwPvvv7979+5brlHy3nvvXbx48cyZM2KxeMuWLYFAIPHVSweTmno4HJ50zjfrOEj1bZgkEonb7Y69jq1vKZPJvqkAhmFOp3NSAXDn9Hr9xo0bX3rppfvuu+/m71ZXV8dePPjggxcuXPjrX/96zz33JLaCC9mqVatiL/bt23f8+PFPPvlk0pAZiUQyMUDJ6XRiGCaVShNcyQVvfHz86NGjTU1Nt/zuxGXi22+/LZVKm5qaVq5cmcDapYtJTZ3JZE7aCUsikRgMhokvpx8HqX5FWFxcPHEK3NbWxuPxVCrVjQWKiop6e3tjp2AIofb29qKioiRUdOEaGxvbsGHDD37wg8cee+y2hWEN6Dl1yz9vcXHxxN6QbW1tOTk5LBYrsfVa+A4dOlRZWVlaWjqdwgjWrZwbk5p6UVHRpE/EjQWcTufIyEhhYeG0fvUdP86cW3q9nsvlnjx50m63b9269amnnood/9nPfvbRRx/FXtfU1Lzwwgsej+c3v/mNVquNRqPJq+9CY7VaCwsLH3nkkYb/4na7EUJXr17dtWtXrMwHH3zQ19dnsVgOHz7M5XI/++yzpFZ5QSEI4t133x0aGjKZTH/4wx+YTGZzczNCyGq1bt++3WazIYSam5uFQmF9fb3FYqmrq3v55ZeTXesFqLy8/MCBAzceefXVVw8ePIgQ6unp+fTTTw0Gw/Dw8GOPPabVav1+f5KquXD09PQ0NDRotdpf/epXDQ0NHo8HIXTmzBmZTNbS0jI2NlZRUfH73/8+Vnjv3r1nz55FCAWDQYVCceDAAbfb/cMf/vCuu+6a5tulehAihP72t7+VlZUpFIrdu3fH/hwIoaeffvq9996Lve7v79+8ebNEIlm9enVra2vyaroANTU1Lf26xsZGhND58+e3bNkSK7Nnz57c3FyZTFZTU/OnP/0pqfVdaAiC2LZtm1qtVigUa9euPXHiROy42WyuqamxWCyxLw8ePFhYWKhSqZ544olQKJS8+i5MOp2upqbG5XLdePDHP/7xG2+8gRDq7OxcvXq1UqmMDZa5eVg1mIV9+/bd2O20tLTEjr/55pt5eXmZmZnPP//8xDXPjh07Jj4aV65cWbFihVQqvffee0dHR6f5drD7BAAAgLSW6s8IAQAAgDkFQQgAACCtQRACAABIaxCEAAAA0hoEIQAAgLQGQQgAACCtQRACAABIaxCEAAAA0hoEIQAAgLQGQQgAACCtQRACAABIaxCEAAAA0tr/B12yMsmzJmFJAAAAAElFTkSuQmCC", "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": 13 } ], "cell_type": "code", "source": [ "ρ = scfres.ρ[:, 1, 1, 1] # converged density, first spin component\n", "ψ_fourier = scfres.ψ[1][:, 1]; # first kpoint, all G components, first eigenvector\n", "ψ = G_to_r(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": 13 } ], "nbformat_minor": 3, "metadata": { "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.6.2" }, "kernelspec": { "name": "julia-1.6", "display_name": "Julia 1.6.2", "language": "julia" } }, "nbformat": 4 }