{ "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.293850389478 NaN 4.29e-01 7.0 \n", " 2 +0.441885799934 1.48e-01 5.35e-01 2.0 \n", " 3 +0.260884086606 -1.81e-01 1.72e-01 2.0 \n", " 4 +0.243049495235 -1.78e-02 2.67e-02 2.0 \n", " 5 +0.242667188386 -3.82e-04 6.38e-03 1.0 \n", " 6 +0.242636270413 -3.09e-05 2.32e-03 1.0 \n", " 7 +0.242633841186 -2.43e-06 8.43e-04 2.0 \n", " 8 +0.242634592319 7.51e-07 1.52e-03 2.0 \n", " 9 +0.242633378403 -1.21e-06 5.26e-05 2.0 \n", " 10 +0.242633376348 -2.06e-09 3.51e-06 1.0 \n" ] }, { "output_type": "execute_result", "data": { "text/plain": "Energy breakdown:\n Kinetic 0.0304511 \n AtomicLocal 0.0972479 \n PowerNonlinearity 0.1149344 \n\n total 0.242633376348 \n" }, "metadata": {}, "execution_count": 10 } ], "cell_type": "code", "source": [ "Ecut = 500\n", "basis = PlaneWaveBasis(model, Ecut, kgrid=(1, 1, 1))\n", "ρ = zeros(complex(eltype(basis)), basis.fft_size)\n", "scfres = self_consistent_field(basis, tol=1e-8, ρ=from_fourier(basis, ρ))\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 Array{StaticArrays.SArray{Tuple{3},Float64,1,3},2}:\n [-0.038684800610639045, 0.0, 0.0]\n [0.038683434756066504, 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+gvaeTAAAgAElEQVR4nOzdeXxU5dk38Ossc87sM8lkD2QhCfsOAoKCC4iKBUWpS10eQYtL3Vpr1T4+Fm2tqFW01qW1oFX0VXFv3XBBAQUUEAXZCQQSIOvsc/b7/ePQGCBAljlzzkyu74c/yGSWO5nc8zv3ThFCACGEEOqpaLMLgBBCCJkJgxAhhFCPhkGIEEKoR8MgRAgh1KNhECKEEOrRMAgRQgj1aBiECCGEejQMQoQQQj0aBiFCCKEeDYMQIYRQj9bdINyxY0csFuv4/TVN6+YrouQihOA2e1aD1cRqVFU1uwjoMMl9R7obhHPmzPn22287fv94PI4fu5Yiy7Isy2aXAv2EEBKPx80uBToMviNWk9x3BLtGEUII9WgYhAghhHo0DEKEEEI9GgYhQgihHg2DECGEUI+GQYgQQqhHwyBECCHUo2EQIoQQ6tEwCNNMtCZBNNyRAKH0Fq0VNAUrslWwZhcAdUKiQfrhqWp7Dlc5q8hT6jS7OAihTlMFbc8HB/d/1Vw8KafsvHyzi4MAMAjTCYEdr9WWTsu3udjNi2qKJub0OiPH7DIh1INs2LDh66+/7s4zyDGlblmTs4D393W/+0xj0Q8BPsuWrOJltpycnIsuusigJ8cgTBv7v24mKik6JQAUeEodGx7fVTwpQDGU2eVCqKdYunTpSy+9NG7cuC4/A1EJFEALQ9XuBFIELTv30yxW4RNramrasWMHBmFPJwblmg/rh9xYDhQAgD3AOfL4li3R7EEes4uGUA8yefLkRx55xOxS9Djr16+fPXu2cc+Pk2XSQ8PaYM5wnzOfb70lb7S//tugiUVCCKHMgEGYHlq2RrMHuNvekjPcG9wWVRJ4TBpCCHULBmEaUCUtujfhrXC1vZG1M/6+7sYNYbNKhRBCmQGDMA2Etsc8pU6GO/LNwt5RhBDqPgzCNNCyNerv5z769qz+7kS9KIXwfHmEEOo6DMI0ENwayWovCCmG8pY5w7sTqS8SQghlDAxCqxOaJVXQXIX2dr/rLnFEa+IpLhJCCGUSDEKrC26J+vu74RiLbj0ljkgNtggRQqjrMAitrmVrtN1+UZ27xBHdh9twI4Q6JxqNtu4MsG7dui+//LLduzU1Nb300kspLJc5MAitLrI77u3jOtZ3WTvD+23xA2Iqi4QQSnfhcPiee+4BAFVV58yZk5eX1+7dAoHAk08+uX79+tSWLtUwCC1NjiqaSnj/8bbldWPvKEIIgBDS1NSkqoc22YhGo63/b70lkTjys+L999/Pzc3t37//sZ72uuuue/jhh5NeWkvBILS0WK3gKmp/mkwrT4kT58sg1JPddttt119//dChQ0eOHPn9999/8803w4cPP+mkk3r37v2nP/0JALZt29a/f/8xY8b07dt3ypQpTU1NrY9dvHjxjBkzAODnP//5a6+9BgDZ2dnhcPiDDz44++yzAWD69OlvvfVWPJ7JHzK46balxeo6EoSOA183p6Y8CKEjfNdE6lKbEWcWUTxz2C3xePyDDz5YtWpVUVFRJBIZNGjQokWLzjzzzFAoNH78+AkTJgwfPnzlypWBQEDTtBtvvPGhhx665ZZb9MeuXLny1ltvBYBoNCqKIgC0tLQQQmRZjkQiAJCdnV1UVPTtt99OnDgxpT9nCmEQWlpsv+CrPOYAoc5VZBeaJFXUGB7b9wil2kf7yBcHtFS+4tg89oggBICf//znRUVFALBy5UqbzUYI+eSTTwBg2LBhn3766WmnnbZ69eqnn346Go02NDTs2LFDfxQhZP/+/ccaIGyVn59fW1ub/J/EMjAILS1WJxRNDBz/PhRDOQvs0X0JX8UJIhMhlHS/G0b/bpj516A5OYeO6a6vr1cURU9BACgpKRk5cuRrr732v//7v3fccUevXr1cLtf777+vf5eiKJ7nJUnSvyTkp/nnhBCKOrRsSxAEh8ORop/EDBiE1kVUkmiQ2h69dCyeUke0BoMQIQQDBgwQBOHee+9tG12XX375Lbfccs011wDAt99+2/b+VVVVu3fv7t+/f3Fx8Z49e1pvr66uLi4uBgBCSE1NTVVVVap+AhNgEFpXvF60Z9to24kvNt29HC1boikoEkLI4k466aRTTjll5syZt9xyi81mW7Nmzbhx4wYMGLB48eIhQ4Zs3br1+eefz83Nbb3/lClTVqxYcfbZZ1999dUzZswIBAIAsGjRovnz5z///PMAsHnzZrvdPnDgQLN+ohTAILSujsyU0Tnz+Lovm058P4RQJpo6dWp2dnbrl6+99trzzz//2muvybI8aNCggQMHjh8/nmGYv/3tb1VVVf/v//2/NWvWeDye22+/HQCuueaan/3sZ/fdd9/48eM//fTTF198EQB27tz5/vvvjxgxAgBeeeWVa6+9trWbNCNhEFpXvMNB6MjnEw0iEDjWTmwIoQw2c+bMtl8yDDNnzpw5c+a0vfHOO+9s/f/o0aMB4P777weAfv36TZgw4Y033pg1a9bQoUMffvjhRx555I9//KPP5wOAaDT65ptvfvXVV6n4McyDQWhd0Tqh6NQTzJTRMTzN8LQYko+/9B4hhI727LPPKorS7rfsdvvq1avd7mPu8pgZMAitq+MtQgBw5vOJehGDECHUWRzHcRzX+mXbuaMsy2Z8CgLuLGNZHdlcrS1HHp+oxx1HEUKo0zAILSq2XzjWGYTtcuTxcQxChBDqPAxCi0rUS468E68gbOXI5RP1knHlQQhZ1po1azZv3pzEJ1y+fHnb/Ujb2rhx486dO5P4WlaAQWhRQqNkz+FOfL//0scIjSsPQsiy/vnPf7733nvJerbq6urrr7/e6/Ue6w6XX35523HEDICTZSwq0Sh6+zg7fn/eb1PiKu44ilAP9Oyzzybx2R555JE5c+bYbO1PUBg8eDDP80uXLj3rrLOS+KLmwg9Ni+psixAosOdyiQZsFCLU48ybN2/hwoX/+te/7rrrrosvvjg7O3v06NG7d+/+3e9+l5eX169fv6+//hoA7rzzzvLy8sLCwnHjxq1atQoA9u3bd9ZZZxUUFJxyyikPPvjgvffeqyjKyy+/fOGFFyYSidGjR4ui+Omnn954440AcO+99y5evBgALrroIn3TmYyBLUJLIiC0yPZAZ4IQwJnHJ+pFd69M3hsXIatRg42aEEvlK9ryegF92PETtbW1drudoqgnn3zyo48+euWVV6677roJEybcf//9Bw8efOqpp2699dbVq1efd9558+bN43n+vffeu/jii3ft2jV79uzhw4d/+OGH1dXVEydOnDx58oYNG+x2e0lJSSwWW7t2raZpLS0t27dvB4Ddu3frq+zHjBnz4IMPpvJHNhoGoRWJQZl1MAzXufa6Iw/nyyCUatEv3xZ+XJPKV8z91UO029/ut6ZMmTJ+/HgAmDFjxtKlS2fPng0AF1xwwW9/+1sAGDNmzMcff1xdXR2Px5ubm7du3frJJ58sWbKEpumKiopLL720oaFh7969BQUFxy9AUVFRXV2dLMvH6j5NOxiEVpRolByd6hcFAABHLt+0MWxEeRBCx+Kbfo1v+jVml+IQfctsAOB5vnX3UbvdnkgkBEE4+eSTBw0aNGbMmKysLJvNVlNTw3Fc66SYvLy8hoYGhmE0rf3jFVsPZlIUhaZpms6ckbXM+UkyidDUyQFCANB3HMWJowih9qxbt06SpJdeeunmm2++5JJLIpFIbm4uTdN6tycAfP/99wBQXl6+b98+QojT6XS73QcOHGh9htYjfPft21daWsowR50OnLYwCK1IaJQ6O0AIAM5cTmiUIKNmNSOEkqOgoGDfvn3Lli3buXPn3LlzWZZlGOaWW2658sor33333QcffPCLL74AgEGDBnEct23bNoqiZs2adfPNN2/ZsiUSiTz55JPr168/88wzAWDVqlWnn3662T9QMmHXqBUJjWLOcF9nH0VzNOtkxKDMZ2VIxz1CqCPGjBmTn59PCGnt58zPz588ebL+f47jZs2a1adPn7/97W/33nsvy7K/+tWvsrOz/X7/Aw88sHDhwnfffXfIkCFXX311PB6nKOrqq69+9dVX/+///u+pp5566qmnFi9eXF1d/eOPP3799df68OGrr776xBNPmPbTGoF0z8SJE5ctW9bx+0ciEU3TuvmiGW/dw9sje+NdeOCGv+4K7oh26iGiKIqi2IXXQgbRNC0SiZhdCnSYcDhMCHn44Yd/85vfmF2WZGpqalIUhRDS0tIyaNCgN998kxBSX1/fr1+/WCym3+f111+fMmVK60OWL19+9tlnp7ic69atGz58eNtb9HckWbBr1Iq6NkYIAPZsm9CME0cRQh2yfPnysrKy0aNH9+vXb8qUKeeffz4A5Obm/uc//1FVtd2HlJWV6Yf3ZhLsGrUcKazQNpq1d2Ug2p7Nic1y0ouEEMpIM2bMOPvss0OhUCAQaDv5paKiovX/55xzzoQJE1q/7NWrV0qLmBIYhJYjNIpdWDuh47NtoR0pXduLEEprPM/rc0GPxeVyuVyulJXHFNg1ajmJzm6u1oY9mxNbsEWIEEKdgEFoOZ3eZbQNe4ATmnCMECGEOgGD0HK6tq2MjvOxclQhKq4lRAihjsIgtByhWbIHOnEkb1sUTXE+G/aOIoRQx2EQWo7YLNuzu74iHldQINQDrV69+owzzujfv/+1114bDuOew52Ds0atRZM1VVBt7q6/L3w2J+AKCoRSRVZlhbS/5M4gdpangGp7y969e2fMmLFo0aIRI0bccccds2fPXrJkSSqLlO4wCK1FbJH5LNvhf+SdY8+yidgiRChVHlnzt2U1X6XyFV+e/mzAkdX2lldeeeWcc84555xzAOCxxx7Ly8sLh8Ot262hE8IgtBYxKPNZXZwpo+MDXMvmaLLKgxA6vrtOvvWuk281twx79uwpKyvT/x8IBLxeb21tLQZhx+EYobUIzd3dMtuezWGLEKEeRVGUffv26f8Ph8PhcLiwsNDcIqUXDEJrEYMy7+9mEOJkGYR6nDfeeGPNmjWiKP7+978/++yz/f72j7BH7cIgtBaxReK7MWUUADivTUlomoJLCRHqQa666qp58+ZVVFTU1tYuWrTI7OKkGRwjtBaxReb93RojBAp4Pys2S468Li5GRAilnfz8/Mcee8zsUqQrbBFai9jSrUWEOlxBgRBCHYctQishIIVkztfdIMT5Mgj1KLfeeqvD4TC7FGkMg9BCpLDMOhma7cYqQgDA+TII9TADBgwwuwjpjQaA++67r1evXr169Zo3bx4hBAAikcjll1+em5vbv3//N9980+xC9hRiS3cXEer4LDyMCSGEOop9/fXXFy1atHz5cpqmzzzzzH79+l1yySV33313KBTauXPnN998c8EFF5x00km9e/c2u6iZT2jp7iJCHe9nxSAGIULoxBobG51Op9Pp7NrDa2trCwoK2p5un47o55577qabbiovLy8tLb3pppuee+45SZL+9a9/3XvvvV6v98wzz5w8efILL7xgdjl7BLFFticjCDm/DYMQoR5lyZIldXV1x7nDq6++Wl9ff/TtV111VXe6/fr163f81z2W1atXr1mzpsuvm1z0li1bhg0bpn8xbNiwrVu31tXVhcPhI240r4Q9iBiUktMi9NnkiEI0XEqIUE9x//33b9u27Th3uOeee3bt2pWy8pzQG2+88fbbb5tdikPY5uZmj8ejf+H1ehsbG5uamux2u8126BPZ5/M1NDQc6/HBYPCMM86gqMPmd9x00033339/u/ePxWKEkCPuj3Sx+oSj1BaNJmGnUMZBBw+GbZ4T91dIkgQAHJeEsUmUFISQeDxudinQYWKxGEVRoiiaXZD2vfTSS/v27fvLX/7yyiuvzJo1a/LkyR9++OGLL74oy/L06dMvv/zyRYsW1dfXz58/Py8v77LLLps0adLRT7JkyZLXX3+dpulLLrlkxowZACBJ0l//+tevvvoKAGbNmjVz5swFCxasXr1aVdWJEyfedNNNrTHRShTFm2+++eqrr16wYAFFUddff/3EiRMBYPfu3Q8//HBNTc2oUaN++9vfbt269eOPP6YoqqmpadCgQTfffPMJf0ZN09p+NurvSEd+OU6nk6ZPsFCQzc7Obj28KhwOBwKBQCAgCIIkSfqHYygUysnJOdbj/X7/0qVL9R+1FU3Tx3lhl8uFQdguNXrAV+hxue3dfyp7FsdKnNt94hnVGIRWo18put1uswuCfkIIcbvdPG/RTSpGjRrl9/tPO+20YcOG9e3b97333vvlL3/5zDPPuN3uG2+8sampafLkyW63+4wzzhgwYEBFRcXRz7Bo0aL777//6aefVlV17ty58Xj80ksvnTlzJsMwd9xxh6qq1dXV+mfF7bffDgB/+MMfGhsbH3jggSOeR1GUv//97zt37rzvvvv27NlzwQUXLFu2rFevXmPHjr3tttuuvPLKxx9/fPr06a+88krfvn0Zhpk1a1YgEOjIz0jTdNtKob8jXf+VHY6tqqratGnT5MmTAWDjxo1VVVWFhYUul+vHH38cPnw4AGzatGnQoEHHeQqGYVgWl2EkgdjtHbdbcX6bGJI9gEuLEDLWjtdqG79P6UG4I39XxXkO+8gdMGCA2+0eNWrUaaedBgCzZ8/+/e9/r7fqHn300RtvvPGWW25xOp0nnXTSuHHj2n3OBQsW/PnPf546dSoA3H///QsWLOjXr9/XX3+9d+9efR6N3tq544476urqDhw4cOWVVz7wwANHB6Huj3/847hx48aPH7927dpnnnlmwIABQ4cOvfPOOwHgueeeKyoq2r9/f1lZGcuyevSYjr366qvnzZt38cUXUxT15JNP3n333TzPX3bZZX/6059efPHFDRs2fPjhh/Pnzze7nJlPFTRCCOtIzuQr3m+TcL4MQsYrn1FYNq0gla/IOk/wKVFdXT148GD9/0OGDKmpqVEUpVMP2bVr1/bt2/v169d2NmlLS8u0adM0TSsvLweAgwcPHuvZWtc1Dhw48PXXX+d5vrU15XQ6+/TpY6nRSgBgL7vssh9++EH/FcyZM+fyyy8HgPnz58+ZMycvLy87O/uZZ57Rf2xkKKFFSsoiQh3vs+FSQoRSgOFpsECPadvxppycnNaJHfX19X6/n2XZ4w9IHfGQnJycvLy8I6Ju0aJF5eXlixcvBoAVK1Z8/PHHx3q2xsZGn88HAA0NDTk5OTk5OT/++GPrdxsaGnJzcwFAX7ZuBTRFUQ8++GBjY2NjY+P8+fP1sb2srKw333wzHA7v3r37sssuM7uQPYKYpEWEOj7LJoYwCBHqKfLz8zdt2qRHy/nnn//EE0/E43FFUebPn3/BBRfod9i4ceOxHn7++ef/5S9/kSRJEIRHH330ggsuGDduHCHkqaee0u+wY8cOTdOCwaCmaYIgPPTQQ8cpzF/+8hdCSFNT08KFC2fMmDF9+vT33ntPz8IXXniBoqhRo0YVFBRs27ZNli3xMYWbbluFFOruSYRt4VJChHqUe+65Z/HixZWVlQsXLrztttsGDx5cWVlZUlJCCNFDa968ef/4xz8qKir0Jl2rnJwcl8s1b948n89XVlbWp0+f3r17//73v3c4HO++++7ixYuLioqKioqeeOKJOXPmhMPh4uLiQYMGTZw4saysDACKi4uPniPSq1ev8vLyAQMGnH/++RdeeOHgwYOffPLJc889t7i4+LHHHnvjjTccDsfll1+uKMrAgQOvvvrqVP2SjonqZuN00qRJ9913X7uTcdsVjUZx1mi79nxYT9FQclZeUp5NDMrfP7HrpP/rd8J74qxRqyGExGIxnDVqKZFIxOPxPPLIIwcOHHjkkUfMLk6HEEJUVe3UTEZVVSmKOmLOvyRJNput9UNbFMXjzJ7V/3QFQbDZbJqmHfHqsiwfveKiI9avXz979uz169e33qK/I114qnZhi9AqpKDMd/vciVacl5WjuKYeoZ6LoqjOzudnGObolW8cx7VtunRwDQlN00e/etdSMAUwCK1CTMYBTK0ommJdrBQ+wVQxhBBKIp7nX3vtNcsG3rFgEFqFFEzmGCHgCgqEUMqxLDtr1qwT7uRiNWlW3AwmhRTOl8x9CXicL4MQQh2AQWgJqqgRLWmr6XV4GBNCCHUEBqElSCGZS2q/KOAKCoQQ6hgMQksQgzKf1H5RAOD9HI4RIoTQCWEQWoIUUpI4ZVTH+1ncXAYhhE4Ig9ASxKRuK6Pj/bjdKEIInRgen2QJUlB2FibhGMK2bB5WiatEJRSD+/gglASVlZWPPvrop59+2u53NVkDmqKTXd1UQWPsPb3FkkgkHA4DD5XDILQEKaz4+yX5vaBoyuZmpbCSxL28EerJzj///D59+qiq2u53d729PzDE66twJfdFN/19d9WlvY44gLAHansgVNL19F+uRYjBZB490Yrz2aSQIc+MUM80dOjQY32L+txbOa7I3SvJDRe6v7dPaaGn1MAYQD29xW0Ryd1frRXvY0XcZQ2hlDCsFtvEINZiY2EQmo+oRE2onDv5rXO9RZj0p0UIHYGoRImrNlcy98TQYS1OAQxC84kh2eZhwYAZLZyPxSqEUApIYYXzsBSd/GqMQZgCGITmk0JKEg9gaovz2cQQdqogZDjJmH5RwAXBKYFBaD4j9lfT8XgtiVBKiCEl6ZtD6TifTcLLWYNhEJpPDCV/fzUd52OxCiGUAlLYuBYhbhpsOAxC80nB5O+vpsPRBYRSQwrJyT1GrdWhWkyMeG50CAah+YzYX03HcDTQoAjtr/9FCCWLaMB2wTqapWieluNYiw2EQWg+KSRzXqN2NuBxgAEh40khmfMatXMF77PhSTKGwiA0nxRSDJosA9g7ilBKSCGZ9xt1Ocv5bFIYa7GBMAjNRkAKy8ZtJMj5cO41QoaTQoqRLUIW10EZCoPQZEpcpXmathn1RmDXKEJGUwSVoimGN6oWc15sERoLg9BkYljmDbuQBL0KYYsQISNJIcWgKaM6zovroIyFQWgyw6sQLiVEyGDGbSujwzFCo2EQmszQyWagb12PLUKEjCQatohQhy1Co2EQmkwKK8atnQAAzo9dowgZy7jtgnXYIjQaBqHJjNuZSWdzM3JMJRruS4GQUYzuGrW5GSWOtdhAGIQmM7pFSNGUzc3KEexXQcgoUsjwWsy6sBYbCIPQZAadat0Wh4uQEDKSGJIN7RoFnPVmMAxCkxl9LQl4GBNCBjPuJLVWvBdnvRkIg9BUBOSo4UGI59QjZByiETmm2tyMoa/CeVkpjC1Co2AQmkmKKqyToWjK0Ffh8Zx6hAwjRxSby/BajBNHDYVBaCbJ+KEFwP2ZEDKSaOSm+a2wRWgoDEIzpWCAEHCYHSEjpehyFkf6jYRBaCajFxHqsEWIkHGMXgGlwxahoTAIzYQtQoTSnRQ2dpdEHYd7JRoJg9BMqWkRsg5GU4kqaUa/EEI9UGouZ21ORhM1TcHNZQyBQWgmMSVVCAA4DytjvwpCBkjN5SxQOMZhIAxCM6WoCmG/CkKGSU2LEPAMCiNhEJpJCil8yqoQtggRMoCYustZFluEBsEgNA1RiSqoNncqgpDHKoSQATSFaKJmcxq7rYyO89qwRWgQDELTSGHF5mHB2P0oDrF5bdgiRCjpDk0ZTUkt5nzYr2MUDELTpGb5kY7HFRQIGUAKK7ZU1WLOi2vqjYJBaBqjD/NsC+ebIWQEKSTzvlQFIbYIDYNBaJrUrMPV4Zp6hIwghZTU1WK8nDUMBqFpUtk1ip0qCBlBCitcKluEeDlrDAxC00ih1FUhhqeBBlXAzWUQSqZU9uuwdoZoRBWxFicfBqFpxLDMeVJUhQD7VRAyQCovZwHAhguCjYFBaJoUVyHex+LxvAglVypbhADA4/G8xsAgNI0UlvkUViFsESKUdGJITuXlLO6yZhAMQnNoCtEkjU3JhhQ6zsfifBmEkkiVNKIBa09hLcbLWWNgEJpDCqVuQwodh5vLIJRUqVxEqMNNgw2CQWiOVM661mGnCkLJlcpFhDrcd9sgGITmSPEYOwBwOMyOUFJJ4ZQOEAL26xgGg9AcKZ4yCtgiRCjZzGkRYi02AAahOVK5rYzuUIuQpPI1EcpkJtRinCxjDAxCc6S+a5RmKZqjlbiayhdFKIOlvhYzPA2AW0QlHwahOaRQqq8lAYD32US8nEQoScSUD3AADvYbg3777bepw61bt+6BBx5oe8vevXvNLmemkcKpO4OpFQ4wIJREUkjmU1+LvbhFVPLR559/Pvmvl156qby8fMSIEQBw3XXXtd7eu3dvs8uZaUxpEeIZFAglkRRJ3am8rbBFaITDukYXLlx4zTXXUNShZd6KgtcdhtAkTVMJ60jdhhQ6PNgToWRREirNUAyX6tElHqd/G+Cnd7G6uvrLL7+84oor9C+ff/55r9dbUFBw//33E3LMuYaaptXU1Gw5XENDg+EFT2epX02v43zYIkQoOaRgqhcR6mxeXFOffD+9kYsWLZo6dareC3rJJZfccMMNfr9/3bp15513XlFR0Zw5c9p9fDQavfvuux0OR9sbr7jiiltvvbXd+8diMU3TWhudPVN0v8C66UgkkuLXVW1yrClxxOtKkgQAHMeluDDoWAgh8Xj8OFefKPWi0eiRtxwUaBeV+lqscUq8WUj961rN0e/IsTidToY5QffboSDUNO2FF1547LHH9C/79Omj/2fkyJFz5859//33jxWEXq93wYIFkyZN6mCZKIpyuVw9PAgFRXNk8R6PJ8WvS+WzjfHoEa+LQWg1hBCapt1ut9kFQYc5ouLEJcWV7Uh9LSb5dMvaROpf14KS+Es41DX68ccfx2KxadOmHX2PUCjkcrmS9XoIDi0/MqNr1I9dowglBw5wZJJDb+TChQuvuuoqnuf1Lx966KGRI0fm5OSsXLny2Weffeedd8wrYQaSwoottetwdTY3I8dUohGK7tEtcoS6TwrJznx76uVT4kcAACAASURBVF8XD6AwAgsAiqLwPH/ttde23qqq6gMPPBCNRnv37v3OO+9MnjzZvBJmICmkuItMqEIUTdlcjBxRUr+EEaEMI4UVf18TWoS0jaZZSkmoqZ92nsFYAGBZ9sUXX2x761133XXXXXeZVKTMZ8pqeh3ns0khDEKEuiv1+6u10jfQxyBMItxizQRSyJyJ14C7rCGUJGLQnDFCwDX1BsAgNIEUTvXpLa3wMCaEkoCAHDVhcygd52NFnC+TVBiEqaYIKlCHdpFPPZxyhlD3SRHF5mTMmnTGeW14OZtcGISpJoWU1G/U2wqnnCHUfVLItGF+AOB9uLlMkmEQppoUMm2MHXB0AaFkMGsRoU6f8mbWq2ckDMJUk8w4w6wV52PFIFYhhLrF3BYhHiOTdBiEqWbi2gkA4L3YIkSou8SwaTNl4NBkGbycTSYMwlQTzTiJsBXrYjSZaLJmVgEQygCmHMnbivOyclQhGm7LnjQYhKlmbhUCnC+DULdJIXO2C9ZRNGVzMnJUNasAmQeDMNWksGmr6XWcl8UBBoS6QwqbvD0TroNKLgzCVDN9hzPOZ8MWIULdYe5kGQDgfNivk0wYhKmlb0jhMbNFyOO2FAh1g6YQVdRsTjO3+sQWYXJhEKaUHFUYO0MxZp6ChNtSINQdUli2eVgw9Sgz3ouXs8mEQZhSYkjm/Saf/IDXkgh1hxS0Ri3GrtHkwSBMKcnUtRM63s+KQQxChLpINHuYH/ByNtkwCFPK3NX0OtyfCaHusESLEI+RSSoMwpQyd3813aHtRnExLkJdIpp3nmgrPFg0uTAIU0o0dcdtHc1SjJ2RY3g5iVBXmHuAjI51Mpqo4RZRyYJBmFJSWOHNvpYEfQUFDhMi1CViSObM7hoFCjgvzpdJGgzClDJ9Ha4OhwkR6jIpaPIuiTpcU59EGIQpZZUg9NtwERJCXaHviWH23G/AiaNJhUGYOlbYkELH+2wSdo0i1HlSRGEdJu+JocMtopIIgzB1pLDMeW3mbkih4/xYhRDqCilogQFCAAAcI0wmDMLUscJqeh2PY4QIdYkYNvMAprY4Hx4jkzQYhKkjWuda0mfDWaMIdYEUlPks69RivJxNDgzC1LHChhQ63o/D7Ah1hWiBRYQ63o8j/UmDQZg6VtiQQsfwNFCgCHjCNUKdIwUtMfEbcIuopMIgTB3JAkdPtMJhQoS6QLTGCiho3SIqirU4CTAIU0cMWqVTBfTLSexXQaiTpJDM+y3RrwMAPC4IThIMwtSRrLAz03/xfpwvg1CnSSHF9O2CW3G4V2KSYBCmCNGIFLHK8gnAudcIdZ6SUCmGYnirfGzyfpuEE0eTwSrvaMaTworNzVK0BZbTA4B+jAuOESLUGaJlZsroOB92jSYHBmGKWGpoAQA4XEGBUCdJIcVStRjXQSULBmGKWGqmDOC1JEKdZ8UWIY4RJgMGYYpYaqYMAPA+FmeNItQplloBBQC8H0f6kwODMEUkC5xN35bNxap4wjVCnSGFrLLRqI7Dkf4kwSBMEdEy+6sdQgHnx70KEeoEoUXmszizS/EThqNpllJiuEVUd2EQpogYtNZkGQCwZ9nEoGR2KRBKG2KLVXbcboVr6pMCgzBFpJBiqTFC0M+pb8EqhFBHWW2MEHC+TJJgEKYE+e+pvFbCZ3EYhAh1kBJTKdpCq+l1uIIiKaz1pmYqKaowDoZmrbKaXoe7rCHUcYJlTiJsi/OxOF+m+zAIU0EKypZaRKjDIESo46xzJG9bPO6enwwYhKlguSmjAADAZ+EYIUIdJbZIVgxCvJxNBgzCVJBCMmexKaOAVQihzrDm5SzulZgUGISpIIastb+ajuFxERJCHSW2WDEI8XI2KTAIU0Gy2BaFrXg/7kyBUIeIFltNr2N4GihQBLyc7RYMwlQQrbf8SMdn4Ug7Qh1i3Vrss0k42N89GISpIDZbcb4ZAPBZNhlbhAidCNGIHFE4n+VG+kGf9YaXs92DQWg8AlLYqteSOMCAUAcoEc1SB2u3hTtjdB8GoeHEkMw6GYqxZBXycxLuu43QiUghxZqdOgDAZ9kEDMLuwSA0nBiU+WzLjbHr+CwbBiFCJyRbOwixRdhNGISGE5st2i8K2DWKUMdIQcWytdieZRNb8BiZbsEgNJwYlOzZFq1CnI9VYipRidkFQcjS5LBq2SDEMcLuwyA0nDXX4eoommJdjBzBRUgIHY8UVCy4iFDH+Vg5quDlbHdgEBrOmutwW+Hca4ROyMotQoqmbB5WCuNgf9dhEBrOgqdat8X5WVxTj9DxSUHrTpYBAN5vE5pxmLDrMAgNZ/Ug9LE4cRSh41BFjaiEdTJmF+SY7Nk4TNgtGITGUhIqIYR1WLcKcX6biEGI0LGJzZIFT49pi8Pp392DQWgsiw8QAoA92yZipwpCxyY0y1yWpYPQno1LCbsFg9BYYlC27NoJHZ9tE5uxRYjQMQnNks3aLULej0sJuwWD0FhWXk2v4/QDKHDqNULHIDRJfLa1gxDHCLuHFQShrq6u9eu8vDy32w0AoVDom2++CQQCI0aMMK94aU8MSlaeKQMANEsxDloKW/TERIRMJzbLniLe7FIcD+/HrtFuYdeuXXv66af37t1b//rRRx+dMWPG2rVrzz333NGjR2/btm3s2LEvvfSSuaVMX2KL7CpymF2KE7Bnc0IzBiFC7ROapYDfZXYpjofhaYql5Jhqc1l3Xp6VsQBQWVn5448/tr317rvvvvXWW++6665QKDRw4MAvv/xy4sSJJpUwvQnWXjuh47JtQpPkLXeaXRCErEholiw+WQb+u+OozWX1y25rogFAVdUffvhh9+7dmqYBQDgcXrp06eWXXw4APp9v+vTpb775psnFTFsWX0So47NYXI2LULuUmEpRFGO3+nQK3HG0O1gAqK+vnz179r59+woLC5csWSKKIk3TxcXF+j1KSkrWrl17rMfLsvz555/X19e3vbFv375Dhgxp9/6apmmaRlFWPJwv6YhKlJjCumn9CsOaNE3j/Gy8RrJyIXsUQoheTcwuCAIAiDcK9gBn/XeE87NCcw+qxR1/R2j6xBcx7MiRI+vr6202myzLc+bMufHGG//85z+zLNv6YJ7nE4nEsR4vSdKnn366YcOGtjdOmzatqqqq3fsLgsAwTA8JQrFZtnlYURLNLsjxSJJEuUmiSRQEweyyIAAAQoggCCxr9b64HiJyMGbz0aIo2myW7tphPFS8MdFzanHH3xG73X7CLGQdjkN9yjabbe7cueeee25BQYEoipFIxOPxAEBDQ0NhYeGxHu9yue67775JkyZ1sPSapjmdzh4ShFJt1B7gnE5Lj72xLMsU2uqCQYuXs+cghBBC8O2wiOZo3JnncDgcFn9H3PlKw74eVItVVU3iD3tYTm7fvj0vLy8/P79Pnz7Lli3Tb1y2bNnJJ5+crNfrUYQm2R6w9LYyOs7HyhE8xgWhdgjNkj07DWqxPcAJTTjS30Xsww8/nEgkKioqduzYsWDBgkceeYSiqNtvv/3mm28Oh8PffPNNbW3tJZdcYnY505LQJKVFEFI0xXltYjA9YhuhVBKbpcAgj9mlODEHBmE3sJMmTVqyZMlHH32Un5//73//e8KECQBw/fXXZ2dnf/DBBwUFBStXrnS5LL2GxrKEJikw2Gt2KTrEHrClS2wjlEp6v44CVs8Yxk5TLCVHFZsbR5c7jR0zZsyYMWOO/sbFF1988cUXp75AmSSNooXP5oRmnHuN0OEIiEGZ99sU0epBCP/tHcUg7AKrL45Ja0KTZM9JjyDEAQaEjiaFZdZB01x6fE7as7EWd1F6vMHpSBFUopJ02fEID2NC6GhCs5wWM2V09gAnNGG/TldgEBpFaEyb5iAc2m4UgxChwwhNEp8moxugByHW4i7BIDRKGg0QAgCP15IIHUVolix+nmhbOMDRZRiERhGa0qlThXOzmqSpQk/ZnwmhjhAa0+ly1pHDCY0YhF2BQWiUNJopAwBAgT2HSzRaejc4hFIs0SA5ci19EmFbnI+Vo4qm4M4YnYZBaBShKT02pGjlyOUT9Xg5idBPEg2iIzdtajFFU5wfZ711BQahUdJrjBAAHLnYIkToJ0pMBQ3Sa1keDhN2DQahIYhGpJCcRsPsAODI5YUGrEIIHZJoEO156XQtC7iCoqswCA0htsg2D0sx6XTIhj2XSzRgixChQ9JrgFCHKyi6BoPQEEKT5EijmTIAgGOECB0u0SCmXS3GrtGuwSA0RNoNEAKAzcUADXJUMbsgCFlCojE9W4QYhJ2HQWgIoUnm02rKqM6RyyVwmBAhAEi3KaM6PIypazAIDZGOLULQe0dxmBAhACAgNEr2dAtCxk7TLCVHsF+nczAIDZGoFx15adapArgzBUL/JYVlmqNZe3psmt+WHS9nOw+D0AAEEmk4WQawCiH0X+k4ZVTnzOPjOOutkzAIk09okWxOhuHT73eLY4QI6dJxgFDnyMN1UJ2Wfh/W1peol9KxXxT0McJGCXCrQtTjJRqldA3CXD5xEIOwczAIky9NBwgBgOFp1k6LIdyZAvV0iXrRkZOWtdiRhwMcnYZBmHyJetGZbjsztXLk8tg7ilCiIW1bhDmc2CLjGRSdgkGYfIkGMU2H2QHAkcsJeDmJejaiEbE5rY5Ra4NiKD4LVxN2DgZh8sUPpmvXKAA48vk4DjCgnk1olDifjbal68ejI49L1GMt7oR0factSxU1RdB4fzqdO9GWq9AePyCYXQqEzBQ7ILoK7WaXousceTwGYaek01FbaSFRLzpyOUinYycO4yywx+owCLtC1qBRgEaBNInQIpJmEYISBCUSliAqQ0SGsEziCuj/RBUEFRIqAYCYDJJ2xJPxAIdNWfJzQFFAA/g4CgCyeGAo8HLgZimXDdws+HnwcZSfAz8HWTyVzUOOHXLtlCddL8nMFN8vOAvTtVMHABy5fKQ6bnYp0gkGYZIl6kVn2vaLAgDnZQFAjig2D/5tHEklsD9O9kRhX4zUxWFvlBxMwL4YqRfgYIJEZMjhIcdOBeyQxVHZPPh58HNUHw+4beDlwGujHSy4WHCwYGeAp8HJUgDgZIFvs4EJISQWi7nd7rYvHZSAEFAJhGVCCAQlUDSIyBCRSUyBqAwhCYISqYlCUIJmQWuRoD4BjQJRCeTaqSIX5DuoQgcUOqnebih2Ur3dUOqmXPgmtyd+QAgM9Zldiq5z5nEH17SYXYp0gvUgyeIN6bqIsJWzwB47IPg97hPfNXNFZNgeItvDZGcYdoZJdYTsjkJdnOTYqVI39HJRRU4ocVOjc6HYSec6IM9O5RjZl+b/77yNHPsRvQ0n6HxIKFAvkP1xqE+Qujjsj5PlB6Aupu2NwZ4ocTBQ5qHKPVQfD1R4qUovVeWD3q607dBIkth+seSsNK7F2DXaWRiESZY4KAaGeM0uRbc4C+3x/aK/qgcFYVCCjc1kU5BsaiGbg2RLEJpE0tdLVfmoCi+My6Muq6TL3FDiptJu/oSDhVI3VeqGdiOzQYDqCNkdITsjsKaBvLxT2xokERn6+6kBfmpgFjXID4OzqXIP1XOyUVOIGJTT+nLW5mYBQI4q+n/QCeGvKckSDaIjbRcR6lyFfHRvwuxSGKs2Rr5tJOubyPom+L6ZNAlkUBY1OIsamEWdV0L380GJu0d89OfaIddOjck97GcNy7A1SH4Mkh9byLNbyKYWaBbJkGxqeIAaEaBG5VCDs9LvgqDjEgdFR4CjmPR+/x15fKJBwiDsIPw1JRUBIQ0P8zyCs8B+cE3Q7FIkWUyBbxrI1/VkdT35poEohIzKoUYGqCsrqeEBuo+3R8ReB3ltcFIudVKbdAxKsKGJfNdElh8gj2/UdkfJkCxqTB41Npcan0+VujPqlxdL85kyOmcuHz8oesudZhckPWAQJpPQIrHpud12W85CPn5AAHLC4SeraxZh+QFt2X6y8iD5sYUMC1Dj8qhfVFKPn0xn2Ge30fwcTCqkJhUe+qVFZVjfRFY3kCXV5DerVZqiTsmnJhZQpxVRA/xpf0kRPyA6C9J47YTOkY9LCTsBgzCZ4vtFZzovP9KxdoZ1MEKLZM9Ovz7ehAJfHCCf1mqf1pGdYTI+n5pUSD86lh6dS6Xh0XIW5bbBqQXUqQUUDAEAZmeYrDxIvthPHt2oRWRyRhF9RiF1Vq90bSnG9guF47PNLkV3OQvs+1c0mV2KtIFBmEyxOsFVlPZBCP+dL5NGQbg1RP5TQz7cp62qJyNzqDOL6CfH02NyKTa9G+fpocJLVXipK6sAAGqi5NM68lkduWet6ueoqb2oaSX0pAKKT5+rkPgBwVmQ9l2jriJ7rBYXBHcUBmEyxeqEdJ8yqnMV8rH9QvYgj9kFOR6VwFcHybt7tLf3EEGFab2pXw2k35xMu3EJuXlK3NTVfamr+wIBZkMT+XAfuW+d+kMzmVxMzyilppXQAWtHjCpoSkxNo0vAY+H9Nk0jOHG0g/B3lEyxOqFkap7ZpUgCZ4G9ZUvU7FK0TyXwWR15o1p7a49W7KRmlNKvn0kPD6RlL1wGowCGB6jhAerOYXSTCP+p0d6tITd9JY/Joy4so2eW0XkOs4vYnpjeHMyIvyZXoT1WJ/j79qB1UF2GQZg0mqyJITlNj245grPQXvt5o9mlOAwB+OogeXmntqRaK3NTF5XTX09n+3gy4hMr0wV4uLKKvrIKEgrzwT7tjWpy1zfySbnUpRX0heW010ot+Ph+IQOG+XUYhB2HQZg08f2iMzftlx/pnPm80CRpCqFZ83+cXRHywjbtpR3EwcIvKjH/0piDhZll9MwySCjMf/ZqL+8kv14lT+1FX1lFT+1liaoTqxPServttlxF9tDOmNmlSA8YhEkTrRNcRZbs7uk8mqUceXysVvCUmvYTCSq8Ua09t1Xb1EIuq6SXTGZGYP9npnCwcFE5fVE5tIjMq7u0+9ar1y6Hq6qoa/rT5l7lRGoSeaP8JhYgiVxF9rrlOHG0QzAIkyZWJziLrD0ToDPcJY5ITdyUINwRJk9v1l7cro3MoX41kP5ZKc3h5M8MlcXDdQPo6wbQPwbJP7dq495RRgSo6wfSPyuhU99A1BSSOCi6ijOkRegs4BMNElGJJdra1oYfMEkTy6AWIQB4ejtSvNEaAfhwHznnQ2XCewpHw+oZ7IdnsxeWYwr2CAP91F/GMnsvtV1ZRT/8vdbnVeXBDVpzaleEx2oFRx6fvufxHoG20XyWLY7L6jsAW4RJQiC+X3Cl/85MrTwljn2pmi8jqPCv7dqCjZqdgVsG029NoXHxe8/EM/CLSvoXlfT6JvLEJq3yNfniPvSvB9NVvlS0aSI1cXdJ5lzLgr6aMINGPY2T6iBsEKjFezWKovQz2Jws5WLBY4NsHgJ2ylLzxzpFDMq0jcqkJTvOArsUlpWEyjoMDKUWEZ7arD25ST0pl35qAnNaIfbhIACAEQFq0UTmYIJ56kf1lH8rpxbQdwylj9gcPOmiexO+CpehL5FiriJ7vE6AUWaXo6tEFZpE0ixCRIaYDC0SAYCzimlfsufmp/qDO6bCuiYAIPoh3XFFiyoQkaFZhCaBCCrkO6heLih0UiVuKHNTfTxQ5aMqvJTF+8didYKrOKOuJYECd7Ejui9h0HlM9QlYsFH9+xZtein96TR2oB8jEB0p3wHzRjG/G8Ys2qZd/Jla6YX/Hc5MMuxqKVKTKD49x6AnN4WrKD02Wtsfh20hsjNCdkfInijURMnBBNTGSEKFHDtk85THBm4bZHEUAJySD2kfhGUu8swEmjrGxryiCg0C2Rs7dA54dYQsrSU7wrA3Rnq7qMFZ1JBsGB6gRuVQVjs7NCP7H9wljmhN8oOwQYCHv1cXbtUuqaDXXsCm6Y6UKGWcLNw4kP5lf3rxDu2XK9RCJ8wbmfw4VARVCsnO/MwZ3QCrbrQWV2B9E/muiXzXRDa2kM1BwjNQ5aWqfFSZmzqjCHq76EInFDopf6pWZVurK49noJeL6uWCIw4+kDXYHiabWsh3TeQfW7TrGwlFwZhcenw+dUo+NTrX/PZirDYRGOozuRDJ5untaFgfSuIThiR4+Hv1mc3apRX0hplsscWuZpCV2Wj4n770FVX0yzu1a1eoJS7402hmbF7S/oSiewVXsYOiM+pvUt9oTQornNfkj/rdEbL8IFl5gKyqJzvCZFAWpTdprqqiB2ZR2WZfflgrCI/FRsNAPzXQT80qP3RLTZSsbiArDpCbvta2h8jJedTpRfTUXtTwgDmnwIT3JMqmF5jxygZylzh2vXMgKU8lqPDXTdojP6g/K6HXz2St1qBH6YKh4IpK+tI+9AvbtVmfqqNzqQdG0/2T0a8eNWmxkNE8Jc7InrgpeyA3CLC0Vvuklny+n4gqOSWfPqWAmtOPHhYwv+lyhPQIwqOVuKkS96FcDEnwxX7t0zpy6WdaSCLn9qZ/VkqdVUw7U/XDiS0yUUkGbNR7BHs2R9TuXk5qBF7eqf3vt9qoHOqLaWxSPrNQD8fSMKcf/YtK+m8/apP+o8wso/8wksnvXopFahK5IzKtUwcAvKWOFAfhhmby7h7y7xptW4icXkRPLqLuGJqcixXjpGsQtuXjYHopPb0UAKA6Qv6zl/ztR+2qZeqZxfRF5dR5JYZvZhjZE/eWZeZJ0J4SR6QmHhjcxVr01UFy6yqVoWDx6cyEfEvXBEOpoiZHFDmuKnFVFVQloSqCpkmaJhNFUIGAKqpE/en+iqKwbEvrl4ydpmiK4WiKpVgHzXAMY6dZB8M6aNbJsi7G5mIyrE+vI+wM/GYIPbsv/afv1MFvyLcPYW4dTHf5vKfo3kT5jEzr1AEAT5mz5uMGo1+FAKypJ0uqtTd2E5qCGaXU/DHMhHwqXdZkZkIQtlXuoX41kPrVQLpFhHdrtFd2ajesVKf2on9RQZ3T26iVsuHdcU9ZBnaqAICn1Bmu7koQ1sXJb1dryw+QB8fQl1b0gA9pAlJEEZoksUUSg7IUUsQWWYoqUlCWowrFUDYXy7oY1smwdoZ1MoydZmw066T5bBtFAcPTrdt/EEJEUbTbf5p7pQoa0YiqB2dcE4OymtAUQVXiqhJX5bgqRxWbi7W5Wc7Hcl6Wz+J4H8tncfZsG5/NWWHDWONk8fDIWOa6AfRvV2uD3lAeHUtPL+10PReapYzs1AEAT4kzVpswbn+ZrSGyeIe2eAfhGJhVTr09hRmanX5/b5kWhK2yeLiqir6qim4W4fVq7eEftF+uUH9RSc/pRyd9pn5kT6LsvPzkPqdF+Kpcu96o69RDZA0WbNQe+l6d25/+x6lsyjqoU0kVtcRBMXZQSNRLiXox0SAJjSLrYPgAZ8+y8X6bPYfzVbo4D8v5bDY306lLMEJILBZzuzszWZeAHFPkiCKGFCkiiy1ydG+i8fuw2CyLLZLNzdpzOEcu78jjnPl2ZwHP+9N2xe4xVHqpt6Ywn9SSm79Wn92iPX4yU+ntRDUPbY/5jFkmZDrGTtuzuVid4O6dzIv1qAz/b5e2aJu2K0wuq6RfP5MZmZN++dcqEz+lDpfNw9z+9Nz+9PYQeX67NuV9tdILcwfQFyVp7y5NIfH9giepf2TW4SlxCE2yHFNtrg51Oa04QK5fqZa44evpbKc+iSxODMrRfYnYPiG2X4jVClJUcebxznzekc/njvQ7cjl7DseYOAGAApubtblZZ+GR3yEaEYOy0CglGqT4QbFlczR+QNAk4izkXUV2V7HD3cvuKrRnxnaUk4upDTPZBRu1k99Vbh7E3DG0oz2lwe0xf1VGLaVvy1PmjOyJJysIf2gmT2/WXt2lTSyg7xxGn9OLZtOk//M4Mj8IW1X5qD+NZv4wEt6r0Z7erN2+Wp3bn7luAN3NMfbYvoQjj8/UDTEpmvL2cYZ2xHKGnaB3NCjBb1erH+4jj42jLypP+9+GJmmRvYlwdTxak4jUxIkG7t4Od7E9d6S/7Dy7I4dLl7NbKZqyZ3P2bM7f96cblbgaqxNidUJ4V6zuyyaxWXIW2j0lDk+p01vu5LPSuL1oo+G3Q+lLKqibv9KGvan8/VRmYsGJ3ioCoe3R0nMz4UjtdnlKHcFtscJTuvUkGoF3a7QnNmpbQzB3AP19Zi1/6kFBqLPR+olo9KYW8tdN2sAl8kXl9G+G0H27uplheE/GDhDqfFWu4Pbo8YNwSbV2y9faBWXUpovY9N0nTxW18K5YaGcstDMW3y86i+zeUmfuKF+fCwrTOhuOxjoZX6XLV3moDaRKWmxfIrIn0bghtOud/TRD+Spc3gqXr8KVpgdN93ZRb01h3tmjXf65ek5vav4Y5jhLs+MHRZqjM3KAUOctc+5d2vX5MgkFXtiuPbpRC/BwyyD6wvKM2Zb8Jz0uCFsNyqKeOYX542jmqc3axH8rEwvou4fTwzt/4l1kdyIw2GNECS3CX+U+8FXNsb57MAE3rFS3BMnrZzLj03BeKNFIdG+iZUs0uC0aqxM8JQ5vhavsvAJPqTOz55i0xXC0t4/L28dVDAAAQqOkXw3UfFRP0eDv687q5/b3dbPONNsKfUYpfXoRfecadcgbytMTmPNK2n9Dg9ujGdwvCgCOXF5JqHJEsXk694EfkeHpzdqCjerYXHrRxEye+N1zg1CXY4f/G0HfPoT++xbtZx+rIwLUvSPpUZ0Z9Y3siWfqTBmdq9CuJDSxRT66VfTKTu22VersvvTLp7NdnrZuCiWhtmyONv8Yadka4X22rP7ukql53nJnBl7rdp49h7PncPljswAgUS8Gt8Xq1wZ3vFbnLOCzB3myBnrSaDdBrw2emsD8vA+5Zrn62i7qifHtNA2D26IZcxhv+yjwlDrDuzuxmjAsw+Mbtb9uUqf2oj85N/O3Au7pQahzsnDrYPr6AfRzW7Xzl6qjcqj7R9FDOjAJWGyRNZXYAxnbqQIAQIGv0hXcEc0/Kav1tkYBrl+pbg6S/0xlO3XdYC45SqapiQAAIABJREFUojT+EG76PhypifsqXNmDvOXTC0zffcrKHHm8I48vPCWbqCS0M9a8KbL5nzVAQWCoN2eo11PiTIux0tMKqQ0z2bu+UYe+ofzjVGZqr58KTTQS3hWvurjYxOKlgLfUEelYEMYVeHyTtmCjenYv+qvMmvJ2HPgR8BOegRsH0nP60c9u1s76QJlSTM8bRZd7jvd3ENwWNehwBkvx93WHtsVag/DfNWTuCvUXldSLp7FpcXCgElMbN4QavgvFaoXsgZ7C8dkDZ5dk6vwmg1AM5e/r9vd197mgMFYrNG0Mb3+1ThXVnKG+3BE+6x/j52LhiZOZ80vJnOXq2b2oR8YyLhYAILpX4LNsmXSAWrt8fd0nXAcla/DcVu2P67VTC6jl57FdnjaRjjL87e8C/WzYOf3oRzdqJ72t/E9f+vfDmaxj7AnbsiWaPSiTBwh1/ipXzYcHgUBUgd+sVj+pJa+ewZxywsl4ZtMU0rwxXP9tMLwrnjXQXTwpx9/P3XNG/ozjKra7iu0lU/PiB8XG70JbF+8FgNxR/rzRfotPOTmjiNowk735K3XkW8q/JjFj86jgtgwfINR5ShxCi3ycYcK3dmt3fqOVe+DfU5kRnZ8qke4wCNvntsH/jaDn9qf/sE7tv0T+/XDm+gFHjh8RjQR3RPvMPGrpVsaxBzibm139XeyKHfyp+dR3M1mPtSdRRmuFg6uaG74LuYsdeSf5+13Rm+Gx/Zd8zny+ZGpeydS8SE2iYW1ww4JdzgI+f2xWzjCfZS84vDZ4fhLz5m7t/KXKDQOZaRtClT2gClM05a9ytWyN5o0+cjR0XSO5bZUakuCv45mzii36rhkNg/B48h3w9ATmVwPpX69Sn9msPTbusNGF6N4E77NxnZyIlY5UApvyvZs/DP75ksILLbxGUJW0hnWhA181K3E1f4x/xG8qM28LFWvylDg8JY7y6QXNmyIHVrfsent/3kh/4YRsR57Z5+scw8wyelwe9Zv3EyOb1OwcpwlHM6Scv5/niCCsT8Dd36rv79XuG8XM7kv3gI0QjynzP8S7b1AW9dE57H/2kpu+Vgf5qUfHHRo4bNkSzeqf+QOE+2LkimVqjsPza2HPuNIis4vTvkSDtH9FU/26oK+Pq2xavr+vOy0mcWQYiqECQ72BoV6xRT6wqvmHp6odeXzRKYHswR4L7gle5KQe9Ea/qfCMfU95fBxzSYV1r/CSIqu/e8/7B4EAUKBo8LcftT99p/5PX3rLLFv6rv1NFgzCjprWm5pcxD66URvzzqHdm1q2RkvPyeSFEwDw1m7t+pXqLYOZ3w11fF9jC++Kt67Ctojg9ljdF42RmkT+2CxsAloEn2UrPSe/5Ky8ph/Ctcsaq989UHhqdsHYbMZurbBp+i409efFH7jZyz5XP6olfz2ZcWfunw/vt9lcTLQ28R1rv+ErtdABX56Hx6IdwgJAdXX1smXLQqHQiBEjJk2aBACbN2/euHFj652mTZvmdGbmMUOdwjNw1zD6FxXUrau0ca+Lf60TfX0y9teSUODXq9WlteSdKax+DnjOcF/DdyGLBCHRSON34drPGzSFFE3K6f8/JZYdlOqxKIbKGe7LGe6L7k3UftH4zSdb88dkFU/Ksch6lfgBURU1T6lzJAVrL2Bv/kod9bbyyunpvXn08dkr3c99GH7UZXs0I/ZBTCJ248aNkyZNOu+883Jych577LHJkyf/85//fOuttxYuXDhy5Ej9TmeccQYGYasSN/XmZObjzyLfOR3PLNf+Mra7J4Ja0MYWculn6tBsat0FP22ZljPcu2HBroqZheZ2c2kKqf+mZd9njZzPVnpuflZ/D/aCWpy7t6Pf5b3FFrn2i8Z1D23PGertdWau6atvG78L5Qzz6X88Lhb+OZF5dZd2zkfKncOYWwdbrye3ewjA89u0Vw86rgs1/3hVfgY3fLuGLS8vr6mpcblcAHDDDTdUVlb+8Y9/BIApU6Y8/fTTZhfPusoPxvqd5t1rh6Fvyg+MZmb3y5ya8/Rm7d616sNjmauqDrtmtGdz9mxbaEfM39eckVFNIQe+bq79vNFVZO97WS9vOV6cpRM+y9bn/MKSKXl1K5o2PL4rq7+79+RcE2fTNHwX6ndZr7a3XNyHHptL/WKZ+kmttmgim5cpF7hbQ2TuCjWhwDMXesTH99sVFWzpsP43hWiXy6WnIAC4XC6KoiiKAoC9e/e++OKLn332maIoppbQijRZa9oULhzhnT+G+eQc9rmt2mn/VrYEidnl6q4mEWZ+oj63VVvxM/aIFNTlj8nav6I59QUjKtm/ovnbP20LbY8NmF0y8JpSTME0xbqYkql5o+/u68jjv3+yetvL+4RGKfXFCG6LUjR19MlEZR7qi2nsiAA18m3l49q0r9GiCvPWaae+p1xUTn89nR2Rz/j7uxs3hM0ul+Uc1ll/5513XnTRRQUFBXa7XVXVzz//fNWqVSzLLlu2LDs7u93HC4Lw3HPPLV26tO2NJ5988uTJk9u9vyiKLMvqWZu+mr+PuHrxhFNFUe3rgs/Ogr9vo059T7uuH/x2MEmvXTclSQIAQsgXB6k5K6iLysgLEwhHq6LYzp19Qx27PzwY2hexp+pQAqKRpvWR2k+anAV85RUFrmI7AIjtFi5T6CfU22wZ3XtFQd6p3sAY98EVLd89vjNroKvozADnS93Y4d7P6gtO9YtS+39I9wyBiXnUnC/IrDIybwThaBBFkeMsvVfA0b6qp25YRfX1wuppWpFTkyUAgKxhrgNfNmeNSPvryI6/IxzHnTBxfrrqf+CBB7755hu9O/TXv/71Bx98sHDhwh9++CEQCMyfP/94T0HT1FE6Ur701fRdJDDip6VHNAXX9SNrziPft1Bj/0N/VZ9mP76swT3rqatXUM+MJw+OIsfZeoy20XnjfAeWt6SiWARaNkY3PV7TuC5ceWlh1VVFegqijMHwdNGZgSG/KWVd7MYn9uz9oFGJqyl43fh+MXFQyh52vB7+Sflkzc/Izgg16UN6aygFhUqmoAS/WkVdsZz6w3DttdO0ojap56tyJuolsUU2r3RWdOgS7LHHHnv++eeXLVsWCATafpthmLPPPnvFihXHerzdbp89e7Y+17QjZFnmeT6tk1KOqdHdwoCrSo84kbych3emwlu7tStXaOee6Ag069gWpv5nBdXbTX83k8ntQND0npS/9s/b+pzHdPZIl04JV8er3zugyVqf84t6wmLNtgghiqLwvEWXoicdz0PFdGfv0/NqPqrf+FhN8WmBolMDhh4DsmdlQ9HEHLvzBH/uhTy8MxX+sUWb/DHcNch56wg+LT62Xq/WblulzSilNl3E+Nr7CMod4QttTPSenN7VSpKkJNYRGgD+8Y9/PPHEEx9//HFR0aHl0oT81Dm+YsWKioqKZL1eBmj8LpQ90MMco910QRm96SLWRsOgJcrLO7UUl61TCMATm7QzP6JmV5G3p3QoBQHA5mJyR/rqljcZVKpEg7R5Uc3Wl/YWTsge8evKnpaCPRbnYSsvKhp6U3m0JrH2we31a4NgzAidGJSbt0QKTs468V0BAODa/vSKn7GLq5lzPlRqY5YeNayOkHM/Uu5bp71+JvO38e2nIADkjfI3rA2mtmhWx65fv37u3Lnjxo2744479JvmzZt3ww03VFRUZGVlrV69es+ePc8++6y5pbSU+rXBkrPyjnMHrw2eHM9cWUXmrlCf36Y9OZ6x4D7u1RFyzXJVVOGLc0hFJ7cNL56Us+HxXYWnBJK7IEyJqTUf///27js+jupeFPiZ2d67drVFvVq92LKsZlvG2MiOSwKEfHLDhYRAgA83vGseL5eEFLiBm4QEXgjlURNawLTYGAcsGSzbkmVbtoolW1bXStu1vc7szLw/1hGKg23JWmlX2vP98MfuzO7OQeOZ35xzfucci/WsS7NenvtvOjguMAFxFKy8f09xj/pH95kMR6bTd6hEmVEetzrxd4uqSjqvZVNyRMihRuzZEV75x+GnqmjfzYq7EXgYCX7fSz7VS/xnEe0/i65RnRakcCmS8k4G+NqVkhe7YHStVvvuu+/O3iSXy5955pkTJ0643e4HHnigqamJw4F/r0t8xmDIjolzrn1xrlEgp3bQ/9RP1u4P352P/qSExo2LYcSApMCf+snHzhL/u4T2YCFK4PPO2WPLmMnrpMMfGPLvSIlKkSiCMhybnmyxyUuF5Q9nM3jLKt0IijZhOrfkgQxrl2vwnSmehp22TcWJUnKW86LXNeQreyhrvl+ko+BnZei2FOSOI8S7I+QLNTQNL14e1FoM1P1tRI4IObWDnnbVZeMuQUBShdjU7si6Gd7YL0Fmt4Jeh4aGhl/96ldz7yP0er2RQRoLOWgMDbyh52s5mg3yuX/F6Ad7OohjZuo3a9BbM2L8LNlrp+45TtAQ8HLdpXpqJGt0vhlxFEGdfWo4ZbNCXipaUIEoYOtxjX1i5qpYadtUXGWidIxdBUVRPp+Pz4dtwoAMU4bW6akvbYoyUcrmJPrCnpAIjDz7m6HMb11Pr7PH4xEIBAAAnARPdJPP9hGPltPuzY/xRNWjHuqhk+RZG/VMNW1byjyKgvuIzicvlu/JYoqWa3LyzBmJiqWupATN2NDfzAABdA4NoSE0Jkrn0hg8GoNPZwrpTBGDJWGwRIz4nCskYAk5B31Zt8xvMetkLnhrA+2oifqPduJP/eTvqmhrFDH43/Pg4LGzxJ8HyccraD/IW+j1i9CQ7Fs1518bF2Xzr7sC5xr2jX1iogiQfasmTmZug+IKSke0G+XKKsnEZ5bOJy9q1svVdbLrXlF5/IBZlMlbYK8zAwWPlqE3pyM/Ok78ZZB8roZWGYsp2VwY+J9u4qUB8seFtDfX0+a7PjaDR1Oulkx+YcvYGadLUIWDRMiOh5w45sIxTzjsI3A/QWIkESKzv62J+pTCSx0ImVJG5s3JCIWEgwRFUARGRv4Pg3bMPebHXHjQjod9YZaUyVWxuEksnprN07A5clY8hEZ9i1VdJ7u+le3qVMjpnfTXB8lvNhM1SuSxCjR7qToOSQq8Pkj+7DR5oxbp3c2I1nwZglSOco2k/+WxgrvT5rtQvXcyMPF3i98USt2apCgXx8PJheIWg0fL3J2srpONf2rufHJQt0mhrJIgtPn9ozEem7b3e0r/V3Ty/vLFyBdN9L8Mkjs+JzZrkV9X0pKXamBeiAAvXiB/3UU06dCuXfTrbqHVbJCf+Z9BXaNiUdO/5yjsJ7xTAZ8h6DeF/OZQwBKiCIotZbIkDKaQzhAy2DImX0ejsVAaC6Vzot91stR/ApSB8LScqzeNkmEqaMP85lDAHLR2ucYOmHFfWJDCFaRyhOk8YTo3JousBu2Yo9+T8cj1P0ChCLgzB/12BvpMH1mzP9yUgv60FM0ULmIQoADYN04+2kkKGOBvm6P/6Jq6VUmEyL4X5xELvRMBfYvVMxHQNcrz70iZ7+0MSlgcBTPvdp13IjD+mWXysE3bKE9aLZljRpXh6LShdbro3vQo3kMRAG7PRnenob/uIoo+wO/KQx8qpkkXs2kfI8ErA+QTXWSpDGm+iV4oWdC1wxTQkyrEU0dsadtU0SrhPFDAbwq6hv2ecb973I97wzw1h69hC1K5yjUSjoLJ4C9pbFoefYRhP+EZ93vGA64Rn1cf4CazJTl8cS5fkMpZsgmgL741yZIxU7dcLV907lwYePoc+ad+YqMafagYrYh2iCIp8PE4+WQ3iRHgsUp0e8oVHx2ur4/wKxQY+ZvRNezL+pZakHrFp2KKoOz9HkPrdNCOaRrkqmrJog4UW9ZgH+E1ucf8k81Wjz6QvE6qrJJcpaEsHCT0n1mnz7mL7k1nSa6/Pe3qPVJ6H/X4WfLDMfKefPT+VdGfhd+FgRfOk3/sJ0ul4BcVUXuiDTnxs08Nlf2vrIX8ZeZ7RMcFr/Oi1zXopXPpokyuII0rSOVyk+bd5hfdPsLlEQhnI8OUZ8zvGPA6B7xBBybJFcgKBZJ8waJWE6d73WOfmMv+M/O6+ye+lhcHrwyQfzhHJnPB3XnoLRnowjNLbUHw5hD5p34yiQMeKkZ3pF7jSWGhgRAAQAHLaef4QbMwnauqlvJ1nJlzgXnCkZM13e3iJrNV1VJ5iTAO12iNKzAQzpHfHDIenbZ1u3hqtqxIJMzgclWsS/+6KOC3hJwDXn2LVbpKkHaTcoENgHO57Y54qN/3ku8Mk7vS0Lty0cjiZQt0ykq9eIH8cIzcnoLuKUKLpFG+diZbrM6LvsJ70haxe4IC3snAdK/b3u/BXLg4TyDJ5Ytz+AscfJXogXA2zIXb+z3T5zyeUb8wkysvFckKhVGPiLg3fPZ3Q/l3pFylxrMQBAU+1ZP/7wJ5zETdpENvzkA2qdH5rpPiwsBBPfn+GHXYQG5PQe/JR6vndh1GIRACAAAgMXKqddrR7/EZggw+HSAg7CcAAoRpXGEmT1EmgkvmzhEMhPNChinHeY+93+MZ8wcdOINLQ2hI2E/QuTRhOlezXs5TR2Fmvrnfdm1B8MoA+epFkomC2zLRHalIwTzbMCkAeuzUh6Pk3lEKI8EPctE7ctBFWuuNIqneZ0flZSJ1nezan57fTwPPRMDW5bT1uFEGKisSSgsEgpSoteHBQPg1wkHC0eexdrncw35xDk9RIZbkC6IzIpsC51+b4KpYqTct+mL01iD4cIx8f5Q8YaHKZEhDMlIiRcpkSAof+dd2RJwEox7qvJM6YaGOmameaaohGd2ZhnwzDb3SjBJfK1qBcAZFUCEnDgCgsWkMLg0mwswXDITXLZJ8R5EUjYVGt5NpvrddCoBjJur9UfJv4xQNAbUqpFaJlMiQbCEi+bp+RFsQ9Dmormmqw0odNpACBrIjFbklA12tWPShZkEb1v1/R4rvT4/Wklh+U8jS6bSddaEMRF4qkpeIuKro953CQHg14QAx3eO2dDp9xqC8WJRUKRamca//XkyBoQ8MfkOw6L70pUzr8IfBMTN13ER220GPnZryUVIWSOJcCof+MHBilCMENDwkVwTWKNAaJVKtRHjXdeFHPRBCCwQDYRxayG23z0EdN1NtZuqcg7roolAEyFiI+B+dYvYQMPopLh3ki5ESGVIpRzaqkVT+kj4/mk859J9ZC+9JY8uv/z6AecLWM07raSfuIxRlIkWFOCp18StZ3uMIFxudQ1NWSZRVkpATt55xDe01kBipqBAlVYjn+7xDkdTQu4agHSu4O22Jkxu5dLBZg2zWXMpwowAw+YE1SOEkAABw6EDMBHI2EtX+SgiCFkWBBCmQID/Mu/TWHgKOEOX8x4ROEhZI5iCcmN6JlaslFAF6nxstuDttvpNaECFyutdtPeP0jAdkRcL0HcmiTN6yawdaaYFwBkvM0G6UazfKfYagpdPZ+/wYg09TlIpkxUKO4tpn2m8KjXxkBAgouCs1ugky1wEBIJkLkrnL7R8XBEH/QsoC0vhbx0K1VoIykN7nRlO3KlVVkmtGMiJEOi54bV1O54BPmMlVrpHk35GyfFPBV2wgnMFTs9PVqvRtKteIz9btOvf8GMpAJPkCUQZPkMa9LHOJIijPuN9yxjXd49bdoEheJ4UD3SAISgRJFWKukjXyscl4fFq7Xi7O5V/Wz0qGKe9kwDPqdwx4PRN+YRpXXiLKukWzGCPcl9jKD4SXIECUyRNl8jJ3Aa8h6LzgMZ90DO2dAgjCEjMYAjpFUOEgEbRgbDlTskpQ8X+y6dxlf3YhCILmjq/lFN+fbutxWzudwx8amUI6g09HmSgRIjEnjnnDXCVLmMZNrpHm35ESk4lNFknCBMIZCOBr2Px/rHWOe8IhF457wggNobFpbBkTLn0AQVAikxcL5cVCiqT8plDYT5A4iTJRlpjBFDFW6uJoiRcI/xlDQI+HqfYgCILiCoIii5r2GVdWTt0WgiAIgq4DDIQQBEFQQoOBEIIgCEpoMBBCEARBCQ0GQgiCICihwYRJCFpqwXAQJ8MAAB/uJynyHxtDOIkDACiKCgQCXOzSUicC5leTjnLoHDpKAwDwmTxk2U1jBUHxCgZCCLpOOIG7MI875HGH3B7M68a8HszrxXw+3O/H/YFw0Iv5AuEARuB+PBAIB8IkEYl8bDqLgTIAAFwGl4ZcapWZ2QgAIEkSRS9t92DemSP68QBBETMbuQwODaHxGFwmjcGhc7gMDofO5tDZXAaXz+TxGTw+k8dn8oRMgYDFF7GEIpaQx1iUpcQgaFmDgRCCrsgd8tgCdrPP6gg6rX6bI+iyBezOoNMRdNmDDpzAhSyBkCUUMvlClkDA5AuYfD6TJ2EncxlcDp0tYPI5dDaTxoyEKDpKj4Suqx907qtPROKiF/PhJB4IB32YP0gEA3jQHw54MK8P81n8Ni/mdWNeN+ZxBd1uzIMTuIglFLPFco5UzBbKOFIZRyLjSOUcmYIrk3GkkRonBCUUGAghCDiCLqPXZPSajV6z0Wc2+6wWv83sszBRpowrVXIVUo5YwZWlirQVqmIxWyxiCaVsMZ/Ji22xuQwO+Oe202vCCdwZckUCuTPotgXsUx5Tt6VvOmA3+2zOoFPEFql4CiUvSclVJPOVkf9UPCUMkNAKBgMhlFiC4eCEe0rvnhpz6Sc9hkmPYcpjpKE0DT85ctPPlWY16NYl8RQqnoJNX2kzazBoDAVXruDKv3YvSZHTAbvJZzX5zGavdcA+9OXEcaPXbPVPK7gyrUCtFapThJoUoTZVqJVzo72mOQTFCAyE0EoWJolxt37EOTbiGB91jY+59I6gUyvQpAg1OqGmRrtGK1BrBMnzqlStYCiCRsJkkSJ/9vYwSZh85kmPYcI9NeqcODLRNu6eDIVDqSJdhjg1XZSaIU7NkqQLWVFbKBWClhIMhNCKEgyHBh0jF+3Dg46RIfvIhHsyma/MEKdlitO2Z21JF6eoeEoUgfmW80NHaVqBWitQr1VXzmz0YN4x18Soc2LYOXZE3zbsGOUyuNmS9GxpRrYkM0+WdaV6JwTFGxgIoeWNpKgx10S/baDfNnBhenDKa0oT6XKkmQXy3J3ZW9PFqSwaM9ZlXJkETH6RYlWRYtXMFqPXPOQYuWgfOTD8+e9PPQ8AyJNm58uzC+R5ebJsmLAKxS0YCKHlJxgO9tsu9lj7zlkv9NsGZBxJvjx3lSxnR87WTHE6TOuIlUgna52uOvLW4rddmB7stw38ufevg44RFS+pUJFfpMgvTipQ8ZJiW1QImg0GQmh5CIaDvdbzXebeLsu5IcdYtiS9SLFqd+62R2v2wK6p+JTElSdx5fW6agAAQRFD9tFe6/njkyefP/s6A2WUKgtLkwrLVcUwKEIxBwMhFL9Iijw/ffG0sfu0qWvQPpwjzSxTFv+g5N9WyXNhg+fyQkNoubKsXFnWt/K2AwAm3FM9lnOdpu6Xut9g0ZgVqpJKVWmFqgQ+00AxAQMhFHdsAfuJqdOnjGfPmHqSePJKVem/Fd5crChg01mxLhoUHSlCTYpQsy3rRgDAmEvfaer6fPSL33Y8myLUrlGXVakr8mU5KAJnQoaWCAyEUFwgKWrAPtg2eap96pTZb12dXFajXfPj1XdL2OJYFw1aXGkiXZpI983c7WGSOGftP2k8+/uTz1v901Xq8nWaNWvU5TDLBlpsSx0I/eHApN2I/CN/nYEy2HQWi8Zi0hhwLFcCwgn8rLn36OSJ45MdAiZ/nWbNf6z+YYE8D9YGEhAdpZUqi0qVRT8s/Z7Vb2ufOv33kcO/7Xg2X55Tp11bo10Dx2MkoGA4hJGYF/MRJOEPBwAA6aIUZrR7RpY6ENoC9qe7Xpx5i5N4MBwKEiGMwLyYj0Nn8xhcAZMvYAkis1hJ2GIpRyznyORcqYIjk3IkS1xgaDFgBNZh6Dwy0X7CcDpdnFqrrXp28/+o+apYlwuKFwqu/BvZW76RvSUYDp4ynj022fFKz1tagbohZd36lBqYX7MyBMMhi99q9U/b/NP2oHM6YHcG3a6Q2xVyezCvB/P6cB8DZTBpTAGTjyJopG3gl3UPJ/OV0S0JQlHUQr7f0NDwq1/9qqGhYY6f93q9PB4PucKI5kA46MW8XsznxrzOoMsZcjmCTnvAafVP2wLTVr/NhwdUvKRkvlIjUKn5yZH5QVS8JFiBuG4YhgEAmMylyD3BCfyEofOLiWMdhs5caVZDyro67Vr4cHOZuU+6nVAIijhr6j2ibzuqb1fykjak1m5IqVXyFEtzdI/HIxDARJ7r58Y8evfUhGtyymua8hiNXrPJZ/bjASVPIefKFFy5lC2WcaQilkDMEolYQgGLH1lB5Ur39uiekfjqI4wsInOVBpAQgUVmRp7yGKe8hg5Dp9495Qy5dAJNmiglU5KWKU7LlKTL4L01npAU2Wnqbh470jZ5KlOStjG17j8qfyhiCWNdLmg5oSG0yuTSyuTSB1ff02U+d3j86F0HH9QJ1Y2pDRtSayVsUawLCH0FI7BR58SQc3TYMTbmmhh1TYTCIZ1QkyLUaAXqGu0aNV+l4iXFz0NwfAXCa2LRmJGu9dkbg+HQuEs/6poYdo6dOv/RoGOEgdJzpVm5sqx8WU6eLBvec2NlwD50aPTLlvGjSq5iU1rDD0tvh88o0AKhCFquKi5XFf949T2nTWdbxlpf6XmzQJ63OX19rXYtTC2OiTBJDDtHL0wPXpgeHJgemvQYdEJtliQtQ5y2TrM6TaSL8ynal1kg/FpsOisyRGlmi9lnuWgfvmAfevf8xxemB6UcSYE8t0iRX6RYlSLSwqW9F5vNP/352JefjRzGCHxz+oY/3vCEVqCOdaGglYaO0taqK9eqK4Ph0PHJjs9Hv3z61Iu12qotGRtLlIXwMl9sHszbaz1/znr+nPX8RftwMl+ZL8vJl+XszL4pQ5zKoDFiXcB5WAmB8F8peUlKXlJkqieSosZdE+dsF7ot/W/2vR8MB0v1svPeAAAT/klEQVSSCkuSCstVxZfVLKEFwgjs2GTHwZGWC7bBhpR1e6ruL1TkwfsRtNjYdFZjWn1jWr0j6Dw0duSPnS97Md+WjMYtGRujnlWR4LyY76y5t8vS22XuM3pNebLs4qRVtxd9e7nPJRtfyTJLwOq3nTWf6zL3njH3YARWoSpdnVy2Ork0YcerRSVZZsgxcmC4uWWsNVuasTWjsU5XDWd+uW4wWWbhBh0jB4dbWsZaMyVpN2VuqtdVLzDhPpGTZcIkcc52/pTx7Gljl949VSDPK1cVlyQV5soyaUjM5vWN7hlJuEA4m9FrPm3qOm3sOmPqUfIUVeryKnVFoSI/oXJQFxIIfbi/eezIgaFDrpD7psxNWzIalyyLbwWDgTBacDJ8fLLjwPChC9ODm9LqmzI3Z0nSr++nEjAQmn3WDkPnSeOZM6aeFKF2dXJZZXLpKnkuA42LdkQYCKOPpMg+20CHobPD0GnyWVYnl9Voq6qSy/lMXqyLtuiuLxD2Wvs/Gfr82GTH6uSybZmby1UlcJG/aIGBMOosftunw4c+HW6RsEXbsjY3ptZzGZx5/UKCBEKSoi5MDx6f6mifPGUL2KvUFWvVFZXJpXGYbwgD4eKKTHTZNnWyy3wuT5Zdq11bo12zgis68wqErpD7s5HDnwwfAgBsy7zhxoyNcXiFLHcwEC4SkqJOG89+Mvz5GVNPva56W9bmVfLcOX53ZQdCjMA6TT3HJk8cnzwpYYvWadfUaNbkyXLi+el2JY8jjAdyjnRb1uZtWZtDBBaZ0uLPvX9V8hR1urX1uurUhMyvoQB1xtTzydDnJ41narVVD1XdN3s5VghaFlAEWaMuX6MudwSdfx85/Ov2pxkofVvW5s3pGxJzfkcf7j9h6Dyqbz9pOJMtzazVVn234ObETC+CNcJrIymyx9LXqj9xVN/OYXDqddX1uuocaWasyxUdV68R2vzTB0daPh1u5jK427M235C+flnnhi0LsEa4NChAdZvPfTJ0qN1waq26sinzhjJV0ZWSnFdSjdAd8hyf7GjVt3db+oqTVtXpqmu1VcuuaQc2jcYMBagL04OtE+1fThwHADSkrKvXVefLc5b1CIGvDYQ4gR+fOvnpcPN528UNqbVNWTfkSrOu8ANQlMFAuMQ8mPfz0S8/HT7kw/1bMxq/NudrBQRCR9B5VH/iiL7tvO3i6uSyel31Wk3l8n2uhYEwLgw5Rlr17Ucm2nx4oE63tiFlXbFi1XJMN70sEA5MD/199HDLWGuGOPWmzE0NKTVwIMQSg4EwVgbsQweHW1rGW7MlGVsyGut11TPz1CzfQGjx21on2lv1bSPO8Sp1RUPKujXJ5Stg/h0YCOPLhHuydaL9iL7N4rPV6qpqtWsrVMVRXyVk8UQC4TTmaBlr/XzsyzAZvjF9w40ZG+EE/7ECA2Fs4QR+bLLjs9HD56wX1mnX3JDWUKEq8Xl9yysQjrv0Ryc7jurbjV5zjXZNvW5dpapkec32cnUwEMYpk89yVH/imP7EoGOkMrm0WrO6WlMpZsX1XMAWv+3w6NEj+jajz9SQUrM5fX2BPC/WhUp0MBDGCUfQ1TLW2jJ+xOi11CSvuSGzoTipIJ5bfQiK6LX0t0+dPjbZgZN4jbaqTru2VFkYz2W+bjAQxjt3yNM2dap96tRpU1eKULtWXVmlLs+RZsVPLvK4S3986uRRffuUx1StrmzQrqvSVcRwkghoNhgI443Bazo40NxuPm0L2Gu1VbW6qnJlHLX6OILODsOZk4bOU8auZL6yWrO6VleVLcmIdbkWFwyEywZOhnssfScNZzoMndNBR7myuEJVUqosTBFql74wPtzfbTl30nC2w9BJUESNdk2tdm2pspDACbBU6xFCcwEDYRyK3HYNXtMx/YnjkycHHSMlSYVV6vIKVYlOqFn68ngxX4+1/6yp57Spy+qfrlCVVKkr1qjL5Rzp0hcmJmAgXJZsAXunseuMuafLfA4j8UJ5XqEiP1+WnS3N5NDZi3RQe8DRZ7twznqhx9o35tLny3JWJ5dVqSsyxKkzn1nKhXmhuYCBMA5ddtt1Y55OY3eH8UynsYsEVLmyqECeV6jITxenLFLLCklRkx7DhenBfttAr/W80WvKl+eUJRWVq0ryZFkrsvHz6uCA+mVJzpHemLHxxoyNAACzz3rOer7PNvDlxLER54SSp8gQp6aLUlNFWo0gWcNPnu/8TxFezGfwmibck2Mu/bBj7KJjGAtjBYrcVfLcu8v+vUCWu5K6yiEohoRMwYbU2g2ptQCAKY+x29LXa+3/8OIBs8+aIU7NlmSki1NShFqtQK3gyq4jSpEUafJZDB6T3jM16pwYdU0MOUZELGGeLDuy8mKONIuOwr6MqIGBMAaUPIWSp2hMqwcAEBQx5tKPOSdGXOMtY61THuOU18hAGXKuTM6RilhCIUvAprEErEv1AxaNFSJCAAA/HvDhPi/mcwRd0wG71T9NUmQyX5Ui1KSJdFsyNj4gvSsxJ4mAoKWkESRrBMk3ZW4CAPjxwJBjZMgxNuqcODLRNukxOoMuBVcm40jEbLGYJeQzeXwGD0VROkqnKIqgCABAMBzCCMwZcrtDblfIbfbZnEGnlCPVCFQ6gSZdnLI+pSZbmpGY098sDRgIY4yG0DLFaZnitMZZG90hjy1gt/mnXZjbFfSEiJAn5I3sshLTbBoLAMBhsCVsjYDJk7DFMo5UzpEKWcspvRuCVh4ug1OcVFCcVDCzBSdwi99mDzocQZcr5PZiPi/mJSkKJ3EUQSPtqCw6U8QS6oQaEUsoZomSeHIpWwIrfEsJBsJ4JGQJhCzB7J48CIKWIwaNEakyxrog0NUkXBcrBEEQBM0GAyEEQRCU0GAghCAIghIaDIQQBEFQQoOBEIIgCEpodADAyZMnX3/9dQDA7bffXlVVFdnx3nvvHTx4UKlU3n///VptDKYEgyAIgqAlgHZ3d2/atCk3Nzc/P3/z5s1nzpwBALzwwgsPP/xwY2NjIBCoqanx+/2xLicEQRAELQpaMBisra199NFHq6qqnE5nc3Pzzp07b7vttqeffvrmm2/eunXru+++y2KxysrKvvb7r7/++oYNG9LS0uZ4PAzDmExmAs41GrcIggAA0Ghw9G4cwXEczv4aVzAMY7GW/WK2K0l0zwja3t6+YcOGyJv169e3tbWZzebh4eGZjRs2bGhra4vW8SAIgiAortBNJpNcLo+8USgURqPRZDKxWKyZib0VCkVPT8+Vvu/z+R555JGZX4jYtm3bd77zna///OSIr/UDWCOMHyRJAgBQFKZNxQuKogiC8NPhrE9xJBwOB+AZiQ/c7T9AhdJAIDDHdiw2m33N+xudxWJFFuIBAIRCIQ6Hw2azcRwnSTLy5WAwyOFccTEEJpN5ww03FBQUzN6YnZ3NZn/90kK4RE5f3QgDYfwIh8MAADq8yOMGRVGhUOhKVxAUE4FA4Cq3QWgpsUQShMXGcXyO18hcnvLpWq1Wr9dH3uj1eo1Go1arKYoyGAyRZNHIxit9n8FgrF+/fu7rEdJ4Qm5pHQyE8QOuRxhvKIqifD4uXI8wnhAeDzd6q99BC4eiaBTbsdBdu3a99dZbFEVRFPX222/v3r1bKBQ2Nja++eabAAC3271///7du3dH63gQBEEQFFfo995773vvvVdbW4uiqNPpfO211wAATzzxRFNT0/HjxwcGBurr6+vr62NdTgiCIAhaFHSpVNrZ2dne3k5RVHV1daSJrLKycmBgoKOjIykp6UoDJyAIgiBoBUABAAwGo76+vqGhYXZHkVgsvvHGG6MeBfV6vdfrje5vQgthtVqtVmusSwF9xev1znTbQ/GAIIjBwcFYlwL6JyMjIzNpngu31EnzDzzwwLFjx5b4oNBVPP/88y+++GKsSwF9pbW19cc//nGsSwF9xel0btq0KdalgP7J7t27x8fHo/VrcPQYBEEQlNBgIIQgCIIS2kKHUYfD4e7u7rl/3uFw9Pb2crncBR4Xipbx8XEURY8cORLrgkCXnDt3zm63wzMSP9xudzgchmckrgSDwZMnTxoMhmt+srKyksfjXf0zCEVRCynNO++88+yzz859XhKv18tms+E8JvEjGAwCAOA8JvEjHA4Hg0E+HFAfNyiKcrvdIpEo1gWBvuJ2u/l8/lzG1L/yyitZWVlX/8xCAyEEQRAELWuwjxCCIAhKaDAQQhAEQQkNBkIIgiAoocFACEEQBCW0xc3exHH88OHDbrd7w4YNly3eCwA4ceLE2NhYeXl5Tk7OohYDmuF2u9vb2z0eT1FRUW5u7uxdnZ2dM68VCkVKSsqSly4R9fX1RRJ3AQACgeCya2F4ePjUqVM6na6mpiYWpUtEs88IAEAkEs3kHE5NTZlMppldpaWlc1wbFro+ZrN5cnIyMzNTLBZHtng8npaWFhRFGxsbLxsUQRDEF198Ybfb6+rqkpOT53WgRcwaDYVCGzdupCgqNTW1ubm5paWluLh4Zu+DDz64b9+++vr6Tz/99Le//e33vve9RSoGNOPChQtVVVWrV69WKBSff/75vffe+9hjj83spdFodXV1DAYDANDU1ARn+Voaq1at4vF4keu8pKTkd7/73cyu995777777mtqampra9u4ceMLL7wQu2ImkDvuuGNycjLy+syZM7t3737ppZcibx955JG33norOzs78vajjz6Co1wWT2Fh4ejoKI7je/fu3bFjBwDAaDRWV1cXFRWFw+GhoaG2tjaFQhH5MEmSW7dutdlsq1atOnjw4L59+9atWzePg1GL5o033igpKcFxnKKon/70p9/61rdmdo2MjHC53KmpKYqiPvvsM7VaHfkYtKgcDofBYIi8Pnv2LIIgVqt1Zi+KorPfQksjPz//2LFj/7qdIIj09PSPPvqIoiiLxSIQCPr7+5e8dAnN7/eLxeKjR4/ObPmv//qvPXv2xLBICeX8+fM4jufm5n788ceRLQ8//PCtt94aeb1jx46f//znMx/+5JNPMjIy/H4/RVFPPfVUpA42d4vYR7hv375du3ZFxs7ffPPN+/fvp/5R+zxw4EB1dbVarQYAbNq0KRgMnj59evFKAkWIxeKZFgOdTkdR1OwmIABAZ2dnW1uby+WKRekS1/nz51tbWy9bA6S3t9dsNjc1NQEAFArF+vXr9+/fH6MCJqi9e/cmJSVd1ihts9laWlouXrwYq1Iljry8vMumXtm3b98tt9wSeR2JKTO79u/fv337dg6HE9n1xRdfeDyeuR9rEQPh1NSUVquNvNbpdKFQyGaz/esuFEXVavXU1NTilQT6V0888cSmTZtmzgIAQCqV/uY3v3nwwQfT0tI++OCDGJYtoXC53DfeeONnP/tZenr6H/7wh5ntBoNBpVJFWqoBAFqtdi6zSUFR9Oqrr37/+99HEGRmC4Ig3d3dTz75ZE1NzbZt20KhUAyLl4Auiymzo8bsXRqNBkGQeV0vi5gsEw6HZ3qSIy9wHI+8xXF89tQ4dDp9Zhe0BF599dW9e/deth6WwWCI3HbffvvtO++8c8uWLdecoA9auPb29sif/cSJE+vXr29qaorky+A4PjsRg06nR3H1NeiaRkZG2tra3nnnndkbf/7znz/++OMAALfbvW7duj/+8Y979uyJUQET0WUxZfYVMXsXgiAois4rpixijVClUs209lgsFhqNlpSUFHmbnJw8uyHIYrHMN8kHum5vvfXWo48+eujQIZ1ON3v7TOXj29/+djAYhCuRLo2ZP/vatWszMzO7uroibyPXyExvgsViiXQlQEvj5Zdf3rp162X3pZmTJRQKd+7cOTvRGloCswPHZVfE7HBjt9vD4fC8YsoiBsL6+vrm5ubI6+bm5pqaGjqdHgqFgsFgfX39sWPHIh1UfX19Ho+nvLx88UoCzXj//fcfeuihzz77bCZNH8OwQCAw+zORPmqNRhOLAiYuu92u1+s1Gk04HPb5fIWFhSiKRvrOcRw/cuRIfX19rMuYKAiCePPNN++8886ZLR6PhyTJ2Z/p7u6e3bMALYH6+vpDhw5FXjc3Nzc0NAAAAoEAjuORcBN5cGxubi4sLJTJZHP/5UUcPmG324uKinbt2pWZmfn444+/+eabW7du/dGPfuTz+f7yl7/ceOONCILs2LHjueee+8Y3vvHf//3fi1QMaEZfX19paemWLVsKCgoiW+6+++433njj6NGj991334cfflheXu5yuV566aVvfvObzzzzTGxLmwi6u7t/+ctfVldXUxT15z//OTU19cCBA3/9619/8pOfjI2NRa6aBx544ODBgy6Xq7W1NdblTRQHDhz4/ve/r9frI1VADMNYLFZnZ+cvfvGLsrIyuVx++PDh48ePd3Z2XtasAkXRiy++ODo6+tJLL61fvz47O/vee+91uVy1tbV79uwJh8PPPPNMR0dHbm5uXV3d9u3b77///uLi4rq6urKysieeeOKpp576zne+M/djLe7qEwaD4bXXXnM6nbt27YqM6jhy5AiO45FM0ZdffnlkZKSqquqWW26Z3SMNLRKj0XhZ5uH27dsnJyeNRmNVVdXevXtHR0f5fH5dXd2mTZtiVciE4vF49u7de/78eQaDUV5evnv3bhRFBwcH29vbIyNrP/jgg+PHj6empt51111wFc8l097e7vf7GxsbI29JknzuueduvfXW06dPt7e3e73erKys2267TSKRxLacK9v+/fuNRuPM2507dyYlJfX19b399tsoin73u9+NTAny4Ycfpqenl5WVWa3WV1991WKxNDU1bdy4cV7HgsswQRAEQQkNzjUKQRAEJTQYCCEIgqCEBgMhBEEQlND+P9ZuXJ79A9s/AAAAAElFTkSuQmCC", "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": [ "ρ = real(scfres.ρ.real)[:, 1, 1] # converged density\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.5.3" }, "kernelspec": { "name": "julia-1.5", "display_name": "Julia 1.5.3", "language": "julia" } }, "nbformat": 4 }