{ "cells": [ { "cell_type": "markdown", "source": [ "# Monitoring self-consistent field calculations\n", "\n", "The `self_consistent_field` function takes as the `callback`\n", "keyword argument one function to be called after each iteration.\n", "This function gets passed the complete internal state of the SCF\n", "solver and can thus be used both to monitor and debug the iterations\n", "as well as to quickly patch it with additional functionality.\n", "\n", "This example discusses a few aspects of the `callback` function\n", "taking again our favourite silicon example." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "We setup silicon in an LDA model using the ASE interface\n", "to build a bulk silicon lattice,\n", "see Input and output formats for more details." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "using DFTK\n", "using ASEconvert\n", "\n", "system = pyconvert(AbstractSystem, ase.build.bulk(\"Si\"))\n", "model = model_LDA(attach_psp(system; Si=\"hgh/pbe/si-q4.hgh\"))\n", "basis = PlaneWaveBasis(model; Ecut=5, kgrid=[3, 3, 3]);" ], "metadata": {}, "execution_count": 1 }, { "cell_type": "markdown", "source": [ "DFTK already defines a few callback functions for standard\n", "tasks. One example is the usual convergence table,\n", "which is defined in the callback `ScfDefaultCallback`.\n", "Another example is `ScfPlotTrace`, which records the total\n", "energy at each iteration and uses it to plot the convergence\n", "of the SCF graphically once it is converged.\n", "For details and other callbacks\n", "see [`src/scf/scf_callbacks.jl`](https://dftk.org/blob/master/src/scf/scf_callbacks.jl).\n", "\n", "!!! note \"Callbacks are not exported\"\n", " Callbacks are not exported from the DFTK namespace as of now,\n", " so you will need to use them, e.g., as `DFTK.ScfDefaultCallback`\n", " and `DFTK.ScfPlotTrace`." ], "metadata": {} }, { "cell_type": "markdown", "source": [ "In this example we define a custom callback, which plots\n", "the change in density at each SCF iteration after the SCF\n", "has finished. For this we first define the empty plot canvas\n", "and an empty container for all the density differences:" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "using Plots\n", "p = plot(yaxis=:log)\n", "density_differences = Float64[];" ], "metadata": {}, "execution_count": 2 }, { "cell_type": "markdown", "source": [ "The callback function itself gets passed a named tuple\n", "similar to the one returned by `self_consistent_field`,\n", "which contains the input and output density of the SCF step\n", "as `ρin` and `ρout`. Since the callback gets called\n", "both during the SCF iterations as well as after convergence\n", "just before `self_consistent_field` finishes we can both\n", "collect the data and initiate the plotting in one function." ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "using LinearAlgebra\n", "\n", "function plot_callback(info)\n", " if info.stage == :finalize\n", " plot!(p, density_differences, label=\"|ρout - ρin|\", markershape=:x)\n", " else\n", " push!(density_differences, norm(info.ρout - info.ρin))\n", " end\n", " info\n", "end\n", "callback = DFTK.ScfDefaultCallback() ∘ plot_callback;" ], "metadata": {}, "execution_count": 3 }, { "cell_type": "markdown", "source": [ "Notice that for constructing the `callback` function we chained the `plot_callback`\n", "(which does the plotting) with the `ScfDefaultCallback`, such that when using\n", "the `plot_callback` function with `self_consistent_field` we still get the usual\n", "convergence table printed. We run the SCF with this callback …" ], "metadata": {} }, { "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "n Energy log10(ΔE) log10(Δρ) α Diag Δtime\n", "--- --------------- --------- --------- ---- ---- ------\n", " 1 -7.774389180708 -0.70 0.80 4.8 \n", " 2 -7.779019439385 -2.33 -1.51 0.80 1.0 25.0ms\n", " 3 -7.779313560012 -3.53 -2.63 0.80 1.5 26.3ms\n", " 4 -7.779350290351 -4.43 -2.89 0.80 2.8 31.2ms\n", " 5 -7.779350639634 -6.46 -3.09 0.80 1.0 24.9ms\n", " 6 -7.779350850604 -6.68 -4.90 0.80 1.0 25.4ms\n", " 7 -7.779350856126 -8.26 -5.14 0.80 3.2 34.1ms\n" ] } ], "cell_type": "code", "source": [ "scfres = self_consistent_field(basis; tol=1e-5, callback);" ], "metadata": {}, "execution_count": 4 }, { "cell_type": "markdown", "source": [ "… and show the plot" ], "metadata": {} }, { "outputs": [ { "output_type": "execute_result", "data": { "text/plain": "Plot{Plots.GRBackend() n=1}", "image/png": "iVBORw0KGgoAAAANSUhEUgAAAlgAAAGQCAIAAAD9V4nPAAAABmJLR0QA/wD/AP+gvaeTAAAgAElEQVR4nO3dd0AUZ94H8N8zM7vAUgTpRVQwSuxSbBAVxB4NtiSaqJeqyZ1nuuZyl3iXmKjJpRf1ziSa802CXRM1imIB7FIswYYFCwgiSIedmfePJSv2BZad3Z3v569lnJ35DeB+eZ55nmeYLMsEAACgVpzSBQAAACgJQQgAAKqGIAQAAFVDEAIAgKohCAEAQNUQhAAAoGoIQgAAUDUEIQAAqBqCEAAAVA1BCAAAqmZFQVhWVvb222+bvn9tbW3zFWP9cPlKl6AkvV6vdAlKwuWreWlMURTNfvlWFISFhYU//PCD6ftXVVU1XzHWD5evdAlKqqqqUvNHYXV1tSRJSlehmJqaGjX/9Gtqasz+07eiIAQAALA8BCEAAKgaghAAAFTNQkEoSVJWVlZ2drZlTgcAAGAiwTKnefLJJ3mev379eocOHebPn9/o41yvpRqRvBxv3X6mVG7ryppUIgAAqJIlWoSHDh26cOHCDz/8sHr16vXr11+6dKnRh/rxtDRwg77w5gGDM/eJU3aITa0SAABUyUJB2KdPHyLiOC4yMjIzM7PRh5oaxj3SmvX7RZ9XWbflrQPi5ovy6kEWatoCAICdsUR+lJSUODs7G167uLgUFxc35Wj/iuCJxLhf9ev6s7nHxQ25ctJwwdPBHIUCAID6NGMQlpaWurq6EpGPj8/+/fsNGwsKCnx9fZt45H9F8KIs9tqk9XaSU0ciBQEAoPGa1DW6bdu2UaNGBQcHDx8+vP72vXv3tmvXLiQkJDg4ODk5OS4ubsuWLZWVlUVFRenp6b169WpazUREeokkmU5flzflqneBCQAAaLomtQg1Gs24ceO6du26Y8cO40ZZlidNmvTaa69NmzYtMTHxiSeeOHfu3MsvvxwbGytJ0r///W9jN2mjGe4LZoyoefeY09M7xTXn5KX9eSfcJQQARW3btm3BggXNfRZRFDmOY0x14+SHDBnyzDPPNMeRWdPXrPvPf/6zdOnSXbt2Gb5MS0sbOXJkfn6+IAhE1LZt2y+//HLEiBH3Pc6xY8ciIyN79+5t3DJ58uRx48bdsts/s/jNl/n1sbXamlIXF5e3M4XFp3hfJ/pfdE3HFipaf6+srMzFxUXpKhSj8ssvLy/X6XQq/Cg0qKiocHR05DirWw/k448/Tk9PnzhxotKF2KG0tLQLFy589913lZWVWq2W53kT3+jo6GgIo3swfzPq9OnT7du3N544LCzs9OnTprzRycmpRYsWb775pnFL+/btdTpd/X0WZstJ+XLScM7TQSgt1et0uo/6kKNWWnlWHpGsfasb99dOavloEEXxlm+Oqqj88iVJUnMQyrLs5ORkhUGo0WjCwsLGjx+vdCF2SJKkVatWGX7tGxSEpvyemD8IS0pK6n9Cubq6Xrt2zZQ3MsYcHR0HDRp0j30mPUCPhVJLByIijuMMV/h+FDf1QblST49vE1Ov0H8e4t21TboEm2C8fHXC5auzc8yA+4PShdxKtT8Ry2CMcfWY8cjm/03y9vYuKSkxfllcXOzj42Oug7to6lLwFq1dWJg72/OIEKCjHqv1afkq6iMFAICmMH8QduzY8ffff6+srCQiURQzMjI6depk9rPckSNPn/XhP+nNjU7Szz4kSkhDAAC4nyYFYWFhYVJS0u+//15cXJyUlJSVlUVEXbp06d69+6xZsy5cuDB79mwfH5+HHnrITNWaJKE1dyBB2HZJHrxRf7nCkmcGAADb06QgPHXq1Lx58w4fPuzn5zdv3rw1a9YYtv/888+XLl0aMGDA4cOH165da/l+81bOLHmEEOPHwlfXbsxFwxAArFqlnu79OVWht1Al6tSkwTK9e/fesmXL7duDgoKWL1/elCM3Hc9odjg/wJ+btF0c04Z92IvXWt2ddQAAIqLnU0QngRbG8HdsNHyQKa0/J6WNwlzp5mLn4TDAn6WPFs6UUsx6/anraBoCgDWaE8ltvSg/u/MOIxs+OizNzxQ/7WPqbIFbhISEXLlypan1ERGRJElJSUmi2NRH/dTW1np5eVVU3OvGVY8ePU6ePNnEE5nOzoOQiLwcae1g/olQru86/f+dxnpsAGB1gl1Y8gh++2X5uV03ZeFHh6U56eJvw4Se3o28wVRdXd30VVMM9Hr9oEGDDAMhm4Ln+bffflurvdcst5qaGkmy3Me1KtrajGhGZy4+kD2+TdyUK38Twzur4roBwGYYsjD2V/H5FHFRDM8xM6Tg7RITEw8dOhQSEjJlyhQHh7q5aFu3bk1JSfHy8powYULLli2JaPv27Z6enl26dCGiixcv7t27d8yYMYZRIN99952Dg8PIkSP9/f1NOWNFRcWiRYsKCgp69+49cuRIw0ZHx7qnq69bt65Lly67du06efJkXFxcbGysua60Qey/RWjUyYPte0TwcKDINfqMq+gmBQDrYsjC5Evy8ynih1nmT8EZM2Z8/fXX7dq1W79+/cMPP2zY+P7777/44ove3t5Hjhzp0aPH1atXiWjJkiVbt2417JCdnf3+++837oxVVVX9+/c/duxYcHDwrFmz5syZQ0R6vX7q1KnV1dVENH/+/NGjRx8/ftzd3X38+PHbtm0zw3U2nLpaRk4CfdaHX3FGGrZJP6sb/9fOHNaBAIDmk3RRnpfVsJtq/jp56SlZkliEN3vrQMPeO/kBblK7OzdvcnNzly5deu7cOTc3tylTpoSGhu7cuTM8PPy9997bv3+/YbZ3QkLCV1999fbbb9/xCAkJCUT01FNPmb7M74oVK1xdXRctWkRE/fr1i4qKevXVV29ZFCYhIWH27NlEdO3atTVr1sTFxZl6teajriA0GNeWi/RiE5PF5Mvy4n48HmcIAM0kwovN7NqwcS4/5kgZVyUHLXk7sBkN/GP9Qfe7/tPhw4fDwsLc3NyISKPRREZGHj582N3d3cHBwbjmSUxMzL59+xpULREVFBS88sorRMRx3JIlS+r/U1ZWlvE5Cg8++KAgCGfOnAkNDa2/T/fu3Q0vgoKC6j/IyJLUGIRE1MaV7XxYeC9DDF+t/2EA388PLUMAMD8PB4oPbMDHy0eHpVVnpG0jBD8niv1V/DlHMtwvbLry8vKqqirjl5WVlTqdTqfTGUbTGGZ7GzYSEc/zxtGh9x7eSUTOzs5jx46lOy21Wl5ebhwUI0lSdXX17Wvlm758dvNR0T3CWwgczQ7nFz/EP5Esztov1mI8KQAoqv7omPr3C821WuThw4ePHDlCRBcuXEhNTY2JiWnTpo2Pj8+KFSuIqLKyMjExccCAAUTUunXrQ4cOEZEsyytXrjS8XavV6nS6oqKiWw6r0+kSEhISEhIeeeSR20+6bt268vJyIlq9enVAQECrVq3MczFmpd4gNIgPZAcThKwiud8v+jOlGEEDAMq4fYyoMQtvmVPRaF26dHn22WdHjRoVFRX15ptvPvDAA4IgfPvtty+99NKwYcO6dOnSpUuXSZMmEdHTTz+9a9eu/v37R0ZG1p/GMH369KioqMjIyPT0dBNP2qFDh+jo6IcffviFF15YuHChFT4zhFTbNVqfjxP9OkT4/IjUa63+y778oyHW+HMCADs2L1OalyUmDRcivG7qXQx2YVuH87EbxGmp4qKYpnYhenp6btiwISMjIygoKCAgwLAxLi7u1KlTv//+u7e3t7G5FhgYmJ2dffTo0aCgIG9vb2Pv6Ny5c999992ysjJXV1cTT9q3b98///nPJ06cePDBBw13KLVabX5+vrOzMxFt3LjROJXiqaeeMsSw5SEIif6YaNjXl01IFn85Ly+I4XX4xgCApYS60ZZht6agQRtXljyc33TBPP1VWq22Z8+et2x0cnIKDw+/ZaNOp4uKijK8NgSYgUaj8fDwaNBJPTw8evXqVX+L8dl89QPVwcHBOLXRwtD6uSHKmx0aLYgyRa3RHy5CNykAWMi4ttwdU9CgjSub9mAjP6uDgoIMo1HCwsIef/zxRtbXWIMGDerbt28j3hgQEKDRaMxez92g4XMTNw0ti+WXnpTiNuj/3p2f0Rl/KACADdu7d6/hheEWoIXPbph62Ah3fJxD88EH/R1MfoDb9bDw/Ulp3FbxWrXS1QAAQHNCEN5ZmDvbPUoI1FGP1frUfHSTAgDYLQThXTny9Fkf/rM+3Ngk/exDZpvKAwAAVgVBeB+PtOb2JwjJl+RBG/WXKhCGAAD2BoNl7q+VM9s2Qng3XYxYrV/cTxjeCuuxAcAdCILw/fffN/dAD+OKaKpy7do147KlZocgNAnPaHY4H+vPTdoujm7DPuzFa9GWBoCbTZ061QJP1KuqqtJqtda5REuzCgoKaqYjIwgboL8/yxorPL9LjF6n/zGOb+emuj/KAOAeXF1dIyIimvssFRUVjo6OKgzC5oNvZcO4aylxID+9ExezXr/sFBbqBgCweQjCxpj8ALd1uDA3U5q8XSyrVboaAABoAgRhI3XyYPseETwcKHKNPv0qRpMCANgqBGHjOQn0WR9+TiQ3fJP+syOYZwgAYJMQhE01ti235xEh8YyUsEW8ivXYAABsDYLQDFq7sB0jhB6e1GVl7eaLaBkCANgSBKF5CBzNDud/GCA8vVOcsVusxXhSAAAbgSA0p4EBLH20cOq6/NAv+jOlaBoCANgABKGZeTvSL0OECSFcr7X6n3PQMAQAsHYIQvNjRDM6cxuGCn8/IE3eLlbolS4IAADuDkHYXCK92KHRgkQUtUZ/uAjdpAAAVgpB2IxcNfS/AfzMblzcBv1nR9BNCgBgjRCEzW7yA1zKSOH7k9LYJPEaJhoCAFgZBKEldGjB9owSgpypx2p9aj66SQEArAiC0EIcePqsD/95H25skn72IRELsgEAWAkEoUWNas0dSBCSL8nxG/SXKhCGAADKQxBaWpAzSx4hPNKai1oj/pqLLAQAUBiCUAEcoxmduVXx/F/TxBm7xRqMJwUAUA6CUDG9fNih0UJ+JfVdpz9ZgqYhAIAyEIRKaqGln+L4v3biHvpF/79TaBgCACgAQai8yQ9w20YI87OkydvFslqlqwEAUBkEoVXo6M72jhI8HChijT79KrpJAQAsB0FoLZwE+qwP/0EUN3QTJhoCAFgOgtC6jGnD7X9ESLooD92kz6tUuhoAABVAEFqdYBe2fYTQ15eFr6797QIahgAAzQtBaI0EjmaH8/8bIDyzS5yxW6zFeFIAgGaDILRecQEsfbRw+rocs16fU4qmIQBAs0AQWjVvR1o/RJgYyvVeq//ptKSX6I7zK4prLF4ZAIC9QBBaO0Y0ozO3cajw9iGp3y/6QRv1pTdn4X+ypbhf9WgwAgA0DoLQNkR4sYMJQqgbO14i9/tFf/2PLFx8XHovQ0ocyDNFywMAsF0IQpvhqqEfBvCf9Oazi+Xuq/Slerb4uPSvdGnrcL6dG3IQAKCRBKULgIaZ8gDXy4c9tE4ftlbrqpV2PIwUBABoErQIbU9YC/bPCL5WZtdrZBcBKQgA0CQIQtuz+Lg0L0tKG1rd1o21X157qULpggAAbBmC0MbcuC/oImeMEVq7sjBkIQBAEyAIbcnyM9I/D90YHcMRZY4Rgl1YxxW11/H8JgCARkEQ2pK4AC55xE2jYziirLHCoEBu+KZb5xcCAIApEIS2xNOBQm8bI8oRJQ7ku3myYchCAICGQxDaA0b0ZV++a0s2/Dc9nnEPANAgCEI7wYi+iua7eLBhyEIAgIZAENoPQxZ2RhYCADSEhYJwwYIFY8aMef755y1zOtViRF8jCwEAGsJCQdi6devnnntuz549ljmdmhmzEPcLAQBMYaEgHDZsWKdOnSxzLjBkYUd3ZCEAwP2ZedHtEydO6PX6+lvat28vCFja29IY0Tcx/Asp4vDf9BuGCC4apQsCALBWZo6ozz//vKSkpP6WTz/91NPT07xnAVMY2oVP7xRH/Kb/FVkIAHAXpgZhZmbm3r17c3JynnjiiS5duhi3Hzp0aMGCBVVVVY899tiIESO+/PLL5qkTGoNj9G2/uizcMFRwRsscAOA2pt4jfOGFF5KSkhYvXpydnW3cmJOTExsb26FDh8GDB0+ZMmXDhg13e/uyZcvmzp2bl5c3a9asgwcPNrVqMJkhC9u6suGb9OX6++8PAKA2prYR0tLSiKh+W5CIvvnmm4SEhFdffZWIioqKPv744+HDh9/x7V27dvX19R0zZgwR+fv733EfSZIqKioSExONW3r06BEaGnq3kiRJkiTJxPrtT4Mu/78x7OldcsLm2jXxvJNdtAvx05ckiTGVPo0SP31cvum//Bx3//Zekz4U9+3bN3nyZMPr/v37v/XWW3fbs0uXLreE6O2qq6vLy8t//vln45aamprAwMB77K/RqPfGV0Mv/6somraHH7VZSuwn2kEWqvynX1VVxfO8aoOwqqqKMWbKB5xdqqqqItM+3+1SVVWVJEk8z5u4v1arve+AzSZ9Iubl5RkHwnh5eZWVlZWWlrq6ujbuaE5OTt7e3itXrjRxf1EUdTpd485lBxpx+Uvj6Kkd4oRUfu0gwdazUOU/fUmSdDqdaoNQlmUnJyfVJgEROTo6qvbyGWNardb0IDRFk76VTk5O1dXVhtdVVVUcxzk6OpqjKmgWPKPv+vN+TuyRLfpK3C8EACCiJgZhq1atzp8/b3h97tw5X19fNfdW2QRjFiYgCwEAiKiJQTh27Ngff/yxpqaGiJYuXTpu3DgzVQXNyJCFPk4sYYu+SlS6GgAApZkahJMnTw4NDT1x4sT06dNDQ0P37dtHRBMmTPD09OzevXt0dPTu3btnzZrVnKWC2fCMvu/P+zixRzYjCwFA7UwdMjFv3rzKykrjlwEBAUTk4OCwefPmrKysysrKiIgI9IvaEEMW/mmHmLBFv2aQ4GjOG88AALbE1CC82+Q/xli3bt3MVw9YjiELpyALAUDdVDoAFwx4Rkv6816OuF8IAOqFIFQ7YxaORhYCgCohCKEuC1s6IAsBQI0QhEBExDNaOgBZCABqhCCEOoYs9EAWAoDKIAjhBp7RDwN4Dwc2JklfjSwEAHVAEMJNDFnormWjkYUAoA4IQriVIQtbIAsBQB0QhHAHPKP/DeBbaNFHCgD2D0EId8YzWtqfd+CRhQBg5xCEcFcajn6OQxYCgJ1DEMK9GLJQy7GxyEIAsFMIQrgPDUeJA3kNxyYmi7WS0tUAAJgbghDuz5CFkkyPb0MWAoC9QRCCSTQc/TyQF5GFAGB3EIRgKi1HichCALA7CEJoAGMWTsD9QgCwFwhCaBhDFtZKyEIAsBMIQmgwLUfLkYUAYC8QhNAYxiycmCzqkYUAYMsQhNBIhiyskWgCshAAbBmCEBrPkIXVIrIQAGwYghCaRMvRini+WkQfKQDYKgQhNJUhCytFGVkIALYIQQhmoOVoZbyALAQAW4QgBPMwZuET25GFAGBLEIRgNoYsrNAjCwHAliAIwZy0HK0YKJTXIgsBwGYgCMHMHHhaGY8sBACbgSAE8zNm4ZPIQgCweghCaBaGLCxDFgKA1UMQQnNx4Gn5QOFqtTxphyjKSlcDAHAXCEJoRk4CrRskFFbJT25HFgKAlUIQQvMyZGFBJbIQAKwUghCanZNA6wYLBZXys7tECVkIAFYGQQiWoBNo3WAht0x+BlkIAFYGQQgWYsjC82VoFwKAdUEQguXoBFo/WDiHLAQAa4IgBItCFgKAtUEQgqUZs/A5ZCEAWAEEISjAkIVnSpGFAKA8BCEoQyfQL0OQhQCgPAQhKMaYhc+nIAsBQDEIQlCSIQtPX0cWAoBiEISgMGQhACgLQQjKc/4jC6ciCwHA4hCEYBUMWXgKWQgAFocgBGthyMKT1+VpqchCALAcBCFYEWeBfh0inCiRp6XikU0AYCEIQrAuxiycmoIsBABLQBCC1TFk4fFiZCEAWAKCEKyRs0AbhgrHi+VpyEIAaGYIQrBSzgKtHyJkFSELAaB5IQjBerlp6LdhyEIAaF4IQrBqbhraNFTIKpJfQBYCQPNAEIK1a6GlTUOFjCL55T3IQgAwPwQh2IAWWvptqLAxV4pao7/ln7KL5TFJIhISABoNQQi2oYWWtgzTHC+Ww1ffyMLjJfKgjeIjrRnPFCwNAGwbghBsRrALHR2nOVFSl4XHS+T4DeJ7kdyUB/BrDACNJyhdAEADGLKw04raXpu0pXqkIACYgSWCUJbl/fv3p6amCoIwduzYgIAAC5wU7FVrF1oTLwzfrCeSV52VC6qkGF8W4cU0CEQAaBRLfHhUV1fPmzevZcuWjLE+ffpcuXLFAicFe3W8RH5qlzgvXN/ZgxHJhVXyX9LElj/UxqzXz9ovrj8vldQoXSIA2BRLtAgdHR1XrlxpeL19+/aDBw8OGzbMAucF+2O8LzjGT/9UR6chG/UhrnQgQSirpT1X5JR86fOj0oRtYmsXFuPH4gNZrD/n5ah00QBg3Sx6j7CgoCAzMzMiIsKSJwW7UX90TGkpuWvpt2HCkI36l/eIn/Tm4wNZfCBPRLUSZRXJSRflpSelqSmih5ZF+7IYPzYokLV1xehSALiVmYNw4sSJZ86cqb9lyZIl7du3J6KKiopHH330/fff9/HxMe9JQSVKamhuFPdEuxv9+e5a2jRU+OKYJMnE/ZFxGo4ivFiEF5tJnChTdrGcmi+n5MnvpUt6WY7x5Qy5GO7FkIoAQERMlk2dilxbW5ubm+vl5eXm5mbcKMvywYMHS0tLe/bs6ezsXFtbe8sBNRoNY6yysjIhIWHChAl/+tOf7nb8s2fPxsbG3pKj91BaWurq6mrizvYHl9+Iy88plVPy6nLxQrnc04fFB3DRvqynD9Pa1FibsrIyZ2dnptYoLy8vd3Jy4jib+pmZT0VFhaOjo2ovv7KyUqvV8jxvxmOa2iIcP378unXrJEn64osvpk2bZtgoimJCQsKpU6eCg4OPHDmydevWsLCw29+r1+vHjBnTo0ePfv365eTkeHt7q/kTHBQU4spCXNnkB4iILldQSr6Ukie/tEc6USJ3bcli/Opy0QmzigDUxNT/8a+//vrChQvHjh1bf+P69euPHz+ekZGh0+neeOON2bNn//TTT7e/t6KiwtXVNScnZ9asWUT04osvDhgwoMmVAzSJv47Gt+XGtyUiKq2lvVfklHxpXpa4O18Oc6+7rRgXwHk6KF0oADQzU4OwZ8+et29ctWrV2LFjdTodET355JM9e/aUJOn2Brubm1tiYuJ9T1FRUZGbm+vh4WHc8tprr82YMeNu+5eXl6u2a4hw+Wa9fEbUuwX1bkGvtacKPWVe4/YUct/+Ts/u5FpqKdZP6u0lxfjIwc7WsqRpeXm5LMuq/QWoqKgQRVG1fYMVFRV6vV61l9/QrlFHR0dBuE/SNakP6MKFC3369DG8bt26dXV19ZUrV/z8/Bp3NJ1OFxgYmJGRYdzi6up6jwuQZdnFxaVx57IDuPxmunwXokHuNKgtEZFxrE3SRfmtTEnDkfWMtVHzPULGmJrvEXIcp+Z7hDzPK3aP8I6qq6s1Go3htVarJaLKysqmHJDjuPotQgBl8Yw6ebBOHuz5MCLijWNtPjsqVejlKG9myMVePljXBsCGNSkI/fz8rl69anhdWFjIGPP39zdHVQDWqP5Ym0sVdaNPDWNtevmwaN+6XMRYGwDb0qT/sr17905OTp45cyYRbd++vXv37o6OWMYDVCFAx8a3ZYaxNtdrad8VOemSNPuQmH5VfvCPsTYDA7iWGGsDYPVMDcK1a9dmZ2efP39+8+bNJSUlo0ePbt++/VNPPfXhhx++9dZbYWFhM2fO/Pe//92stQJYJzcNGde1KddTeqGcmi8vypae3ikGO7MYPxbtywb4s2AXld7SA7ByDWsRPv/88/W/9PLySktL++abb5KTkxcuXDhy5Eiz1gZge5wFivFjMX5sZjdOL1FmkZySJ/9yXn5lj9jij8Xeon1ZRw+1DnQBsD4NWFmmuWFlmQbB5dvc5eeUykkX5ZQ8eVe+XKWXI71ZjC8XH8h6eDKugamIlWXUPGoUK8tY16hRADBdiCt7PswwAPXGWJupKdLJErnnH2NtYvyY413+g+eUyr8X04hWN4VftUhLTkrPhTU0SQHgBgQhgALqj7W5Ukl7C6TUfNk41iY+kEX7sof8OHftjbcUV9MzO/ULYviE1nVNgWqRxibpHQX2dHsSVNo8ADADBCGAwnycaGQwNzKYiKj+gxUnJteNtYkPZAP8uXAvtmmoMGyTXpZpkBfVSPToNtFRYD/G8khBgKZAEAJYEZd6A1CNY22W58jTUmrdtSzalz0bxj2/S/w4kltxQdJw9GMsj7n8AE2EIASwUsIfD1ac0ZlEmc8qknflybvyZInoT2kad638545s7TmpW0sW4sZ43CQEaCwEIYAN4Bn18GQ9PNnUMBqbJBdUisevc9kldLiIZhZJ+ZXyg+6sa0vWpSXr0pJ1bcm8sbIFgMkQhAA2w3hfcFNc7dka3fDfxK+judFt+BqJTpbIBwvlY8Xy1iPS/gK5RqKO7izCi3XyYB3dWbgX0+H/OsBd4D8HgG2okWj8VtFwX7C6grp7so1DhWGb9EQ0ug1nWBzcuPO1ajp6TT5YKB8slJeelDKLZB9H1tGDDLkY4cUedG/w5EUAe4UgBLANv12QtBz9Xyyv4aiaiIi6e7Jfhwh/SRMfDuZuGTLj4VC3wI3hS71E58tlQzT+cl6elynllMohrn80GT2opzfn62TpKwKwEghCANswMph7OJhuacWFe7HUUcJ9m3YCV/foDMMkDSK6XksnS+qiMemSnHlV1Ms39aZGeOExGqAW+E0HsBl3DLzGdXC6aeqGpBqeKkVElyrkY9fo6DU5JU9elC0dL6nrTTVGI9ZHBXuFIAQAIqIAHQvQUXxgXdjVSnSiRD5WLB+9Ji89KR27Rler5XZuzNhq7OaJsalgJxCEAHAHGo4MA3AM68ARUUkNHS6qi7evi3cAABrUSURBVMblZ6TMItlVc2P0jWHnuy2UCmDNEIQAYJIW2psG4BDRpQr5YKF87BolXbwxAMcw+sYQjW1d0ZsKNgBBCACNFKBjAcE3BuDU701dlC0du0ZF1bKhpWhoNXb3ZC4aRSsGuBMEIQCYx+29qYbpjMbe1IyrspuW6g9MDXM3dXG4gio6ek0e4H/T3pJMq89KY9tiuVVoEgQhADSXW6Yzkgm9qSGudw7Gs6Xyo1v1SwcIQ4PqdpBkej5FPFEiPxzMOeDeJDQBghAALOeW3tT6i8MtypYOFMjVEnV0v9Fk7OHFnAUioihvtnawMGqz/tt+fJwnyUTT08TsYnnjUAEpCE2EIAQAxWj/6E01bjEuDnesWF5+Rkq/Kreo15v6WR/+6Z3iVz1Z8hXp8DXaOFRwxU1HaDIEIQBYkdsXhztxXT5cJGcVySvOyEeuyeV6emKXprWLnDFGg6E3YBYIQgCwXgJHHd1ZR3f2WAgRkSTTUzvFrRfFc2Vse578cCvMzgAzwGgrALANMtFf0sTT1+WDI/Q9fejxrfr15yWliwJ7gCAEABsgE/05VcwqkjcOFVwE+X/9OQ1PU3aIyEJoOgQhANiATbnykWvypj9Gx7R2Ya914Tt7cLP2SVWi0sWBjUMQAoANGNaKbRsu1B8d83pXrrBKfi+CwwKn0EQIQgCwDcLNH1dajhbE8C/vlcpqFSoI7AWCEABsVT8/9pAfey8DfaPQJAhCALBhH/Xivz0uZRbJShcCNgxBCAA2zNeJ/hXB/yVNRBJCoyEIAcC2PR/GiRJ9fwLzKKCREIQAYNs4Rgti+Df3i4VVSpcCtglBCAA2r2tL9ngo9+Z+jJqBxkAQAoA9eDeC33RB3nEZ9wqhwRCEAGAPXDX0cW/uL2liLe4VQgMhCAHAToxvy7V2oc+OIgmhYRCEAGA/PuvDz80Qz5aigxQaAEEIAPYj1I3N6My/sheNQmgABCEA2JWZ3bjsYhmPZwLTIQgBwK5oOfommp+eJpXrlS4FbASCEADsTX9/Fu3L5qRjWiGYBEEIAHbo4978f49LWViMG0yAIAQAO+TrRLPDsRg3mARBCAD2adqDXK1ES09i1AzcB4IQAOwTx2hhDD9zHxbjhvtAEAKA3erakj0awr11AKNm4F4QhABgz96L5Dfkymn5uFcId4UgBAB75qahj3px01KxGDfcFYIQAOzcYyFcK2f6Aotxw10gCAHA/n3Wh38/QzxXhg5SuAMEIQDYv3ZubHon/lUsxg13giAEAFV4szt39Jr8y3k0CuFWCEIAUAUtRwui+b+kiViMG26BIAQAtejvz/r6svczMK0QboIgBAAV+aQ3/9/j0rFidJDCDQhCAFARXyf6Rw9+WgoW44YbEIQAoC4vPsjVSPS/UxhBCnUQhACgLhyjr/ryr+3FYtxQB0EIAKoT4cUebcv9HYtxAxEhCAFAneZE8b/myruv4F4hIAgBQJXcNDS/Jzc1BYtxg6WCUJbl/Pz8yspKy5wOAOC+JoRyQc701TEkodoJFjhHdXV1r169AgMDL1y40Lt37wULFjDGLHBeAIB7+7wP33edflxbFuSMDyX1skQQOjg4ZGRkEJEsy126dDlx4kSHDh0scF4AgHtr58Ze7Mi9vEdaPpBXuhZQjOW6RnNycnbu3CkIQmBgoGVOCgBwX2924w8XYTFuVTNni1AUxe++++6Wjc8++ywR1dbWzpo16+TJk/Hx8c7OzmY8KQBAUzjw9EVf/rldYmyA4GyJPjKwOmZuEdbexrBdq9UmJiamp6cfOXIkKSnJvCcFAGiKQYGstw+bm4lphSpl0t8/kiR9+umnBw8evHLlyhdffBEWFmb8p3nz5i1evJjjuGnTpr300ksvvPDC7W8XRZHneeOhzFI3AIAZfdqH77aqdmIo96A7Rs2ojkktQkmSjhw5Eh0dnZqaWlJSYtyemJj4zTffrFmzZvny5R999NEvv/xyx7cfOnQoMjJy/PjxkZGRAQEBAwcONE/tAABm4udEb3XHYtwqxWS5AT/3Fi1abN68uVevXoYvBw0aNGLEiJdeeomIPvjggz179qxdu/aOb6yqqrp06ZKXl5ebm9vdDn7ixImePXsa7ikaDBw4MDY29m77l5aWurq6ml68ncHlq/nyy8rKnJ2dVTsNqby83MnJiePMfGdHkqn/Ju7FMHlCW6tOw4qKCkdHR7Nfvq2orKzUarXGXsb70mg09925SbeGjx49+o9//MPwOjw8/Ntvv73bno6OjiEhIfc+GsdxPM+3bNnSuMXT0/MeF8DzvOnfC/uDy8flqzYIDZdv9iTgib7oQ2O2sRGtyMPBvMc2p2a6fFvB/8HE/U35b9KkILx69aqxhdeiRYuCgoKmHE0QBDc3t7/97W8m7q/RaDQaTVPOaNNw+bh81Qah4fKbIwl6+9HoNuLsTPo62nr/zGq+y7cJer3elEZegzTpW+nh4VFaWmp4ff369fqNOQAAW/RBFL/uPBbjVpcmBWFoaOjvv/9ueJ2dnR0aGmqOkgAAFOOmoXlR3LQUUY8R7qphahDm5ubm5ORIknTp0qWcnBy9Xk9EU6ZM+eqrr0pKSq5du7ZgwYIpU6Y0Z6kAAJbwRDvO2xGLcauIqfcI//SnP509e9bHx+e1114joh07dgQFBT399NOHDh0KDg427DBx4sRmrBQAwFK+ieaj1+vHtWWBWIxbBRo2feKODHPkm37n9uzZs7GxsWfOnDFxf5UPoMflq/nyMX2iOaZP3OLtg+LxEvo5zupGzWD6RIOmT5jCDN9KjuNU+yMBAHv1Vnc+46r8ay5Gzdg/BBgAwB0YFuOesVuswhKk9g5BCABwZ4MDWbgn+yADSWjnEIQAAHf1eV9+QbaUXYwOUnuGIAQAuCs/J/pbN35aKhbjtmcIQgCAe/lLJ+56Df14GtMK7RaCEADgXnhGC2P41/dKxTVKlwLNA0EIAHAfUd5sVGv29wMYNWOfEIQAAPc3N4pfc07eg8W47RGCEADg/lpoaW4U95c0EcNm7A+CEADAJE+241po6Wssxm13EIQAAKb6Jpp/N128WI5WoV1BEAIAmKp9C/Z8GPfaPjQK7QqCEACgAd7qzu8vkDdgMW47giAEAGgAJ4G+wmLc9gVBCADQMEOCWHdPNjcTSWgnEIQAAA32eR/+62PS8RJ0kNoDBCEAQIP56+jN7vz0NDQK7QGCEACgMf7aiSusop+wGLftQxACADSGYTHuV7EYt+1DEAIANFKUNxsZzP6BxbhtHIIQAKDxPojiV5yR9mIxbluGIAQAaDwPB5rXk/8zFuO2ZQhCAIAmmfQA10JLC37HqBlbhSAEAGgSRvRNNP/PQ+KlCrQKbRKCEACgqdq3YM+Fca/tRaPQJiEIAQDM4O/d+b1X5KSLaBTaHgQhAIAZOAn0VTT/YhoW47Y9CEIAAPMYGsS6eLD5WeggtTEIQgAAs/msD/f5ERGLcdsWBCEAgNkEOTMsxm1zEIQAAOY0oxNXUEWJOeggtRkIQgAAcxI4WhjDv7wHi3HbDAQhAICZ9fRmI4LZOwfRQWobEIQAAOY3N4pPzJH2FWDUjA1AEAIAmF9LB/ogiv9zKhbjtgEIQgCAZjGlPeeqoUXZGDVj7RCEAADNghF9Gc2/fRCLcVs7BCEAQHPp6M6e7cC9sQ+NQquGIAQAaEZv9+B358tbL6FRaL0QhAAAzchJoC/78i+kYjFu64UgBABoXsNasU7u7CMsxm2tEIQAAM3u877cp0fEE1iM2yohCAEAml0rZzazGxbjtlIIQgAAS3i5M5dfScvPoIPU6iAIAQAswbgYdwkW47YyCEIAAAvp5cOGBrHZh9BBal0QhAAAljO/J//TaSn9KkbNWBEEIQCA5bR0oDlR/NQUUUIUWg0EIQCART3VnnMWsBi3FUEQAgBYFCP6si//j4Pi5QqlSwEiQhACAFheJw/2dHtu5j6MmrEKCEIAAAW8E86n5MvbsBi3FUAQAgAoQCfQV335F1LFajQLlYYgBABQxrBWLMydfXQYo2YUhiAEAFDMF324z46IOaXoIFUSghAAQDHBLuy1rvyLqegeVRKCEABASa905vIqaCUW41YOghAAQEkCR1/25f+6G4txKwZBCACgsBg/NiSI/SsdHaTKQBACACjvw178/53CYtzKsFwQVlZWTp8+fffu3RY7IwCArfB0oPcisRi3MiwXhO+8805KSsqJEycsdkYAABvydAfOWaD/HseoGUuzUBDu3bv3+vXrDz30kGVOBwBgc4yLcV+pVLoUlTFzEMo3M2ysrq5+6623PvjgA/OeCwDAznTyYJPbcW9gMW7LEsx4rP3797/wwgv1t3Tv3v2///3vnDlzpk2b5uHhYcZzAQDYpX9G8J1X6rddkuMCmNK1qIWpQXjhwoWDBw8WFBQ8++yz9bfv2bNn06ZN/v7+EydOjIqKOnDgwO3vraysXLhw4cKFC48fP75t2zZ/f//BgweboXYAALujE+jj3twLqWLWGMGBV7oadTCpa3T79u0dO3Z85513bmnwJSYmjho1ShCEjRs39uvXr7a29o5v//DDD7ds2bJly5aEhITXX38dKQgAcA8Jrbkwd/bxEYyasRCTWoTR0dHFxcVHjx4NDw+vv/3dd9/99NNPJ06cKIpi165d165dO27cuHsc57HHHvP29r7HDnq9Picnx/ilj4+Pi4uLKRUCANiTL/pwkWv0j4WwEFd0kDY7k4JQo9HcvjE/P//IkSPDhw8nIp7nhwwZsnXr1nsHYXR09D3+taKiIj8/f+DAgcYt06dPf+655+62f1lZ2f1Lt1+4fKVLUFJ5ebkkSYyp9COyoqJCr9dznD2vB+JB9Of2wtQd4qr+t/a0VVZW1tbW2vfl30NlZaVWq+V5U3uNHR0d7xhh9TV+sExeXp6Dg4O7u7vhSz8/v9TU1EYfjYh0Ol1gYOCZM2dMf4urq2tTzmjrcPlKl6AYxpizs7Nqg5DjOCcnJ7tPgrciKWKNPumqdnSbm66U53lHR0e7v/y7EQShQUFoisZ/KzmOkyTJOEdCFEVBMOcYVAAANRM4+qov/9IeqezOoy/AbBofhP7+/rW1tVevXjV8efnyZX9/fzNVBQAAFOPH4gLY7EOYVti8Gh+EXl5ePXv2XL16NRHV1NRs2LDBcL8QAADM5aNe/LJTUgYW425OJnVmFhcXjx8/vqysTBTFQYMGeXl5/fjjj0T0z3/+84knnsjMzMzIyAgODh46dGgzVwsAoC6eDvSvCH5qirh7lMCp9KZwszMpCHU63cyZM41fOjg4GF4MHTp07969O3bsGDRo0PDhw1V78xYAoPk804H7/qT07Qnp2Q74jG0WJgWhVquNj4+/4z+1a9euXbt2Zi0JAABu4BgtjOHjN+hHBXM+TkpXY4/w9wUAgLXr5MHGtuFm7r911ExhlSLl2BsEIQCAtdt+Wd6YKyddlJIv3xg1s+acFL5aX65XsC47gSAEALB2sf5sZneuRmJTU8RaiYjo11x5Woq4Mp53xvztJsO3EADABkwN44jolT3Sv7L42CB6dpd+/WAhyhsDSc0AQQgAYBumhnEXy+U56dx/T4or45GCZoMgBACwGeFezEmgGpHGbdV7OrBePqyXN+vpw8I9GR5e2GgIQgAA22C4L7hxoP5omXZelrwslsstp5Q8efkZKeOq3NqFxfixaF8W4cU6eqh1RfZGQRACANiAX3PlZ3bq1w8WOjnXRgcyxtiEZGnbcH58W46IaiXKKpJT8uSki/K8TOlihdzZoy4X+/hwXo5KV2/dEIQAANYuJU9+Zqd+wxAh3ItVVBARTQ3jqvQ0aKOYOUZwFkjDUYQXi/CqawfmVdL+AulgobwoW5qyQ3TkKcaXMzQWo7zRiXorBCEAgLXr6cO2jxDC3G/q75zRmYsNYHecPuHnRCODuZHBRESiTNnF8sFC+WChvPyMlH5VftC9rgc1wot18kAfKoIQAMDqaTm6JQUNura8f4zxjDp5sE4ebPIDRERltZRxVT5YKCddlN85JFXp5UjvulCM8eU8HMxeuw1AEAIAqIiLhmL8WIxfXYJeqqhrLH5+VHoiWWzlXBeKEV6slw/TqGPNFQQhAIB6BehYQDAzdKLqJTpeIh8slFPz5UXZUm653LVlXSj282NtXO22ExVBCAAAREQCd1MnamktZV6VU/Pl5Wek1/bKwh/jcWJ8ub6+TGdH6WFHlwIAAObj+kcn6kzi6I9O1NR8efYhMbNIDnZmEV51MzQedGc2/dBgBCEAANzf7Z2oqflySp782RHpQrncpWVdKPb24bxtbdoighAAABrG2In6fBgRUUkN7S+QU/KlRdnSUztELU+GHtRoXxbpzRytftoighAAAJqkhZbiA1l8YF3i5ZTKKXnywUJ51n7bWPsNQQgAAOYU4spCXOtG3JTrKb2wbtqiYe23KO+6UOzry3lax7RFBCEAADQXZ+HO0xYXZUuTt4stHW6scdPTh2nvMm0x/ar89TFpQQzP12tOFlbRtFTx2368m6apRSIIAQDAQuqPuKm/9tsPp6QTJTemLcb4sZB60xY7urO8SnnCNvH/Yut6XwuqKH6DfkgQa3oKEoIQAAAUccvab4ZpiwcL5V/Oy2/sEzlGxjVuHvLjVgwUxm3VT0wWF/emkioa8pt+SBCb39M843AQhAAAoDzjtMUZnYmIP1Mq77ki7y2QP8iQJhSJD7RgkV5sf6H08Db+Wo00tBVnrhQkBCEAAFihtq6srSubEEpEVCNRxlV57xX5ajW3/rzUxoV9EGXOORnqWFEVAABslpajnt7s8VDu9HX5xQ5ShxbyhG2iXjLb8RGEAABg7YyjY+aGS4mxXKUoT0w2WxYiCAEAwKoZU9BwX9CBpxUDBTNmIYIQAACsWnG1PDH0ptExDjwtHyiEuFG1yoNw0aJFxcXFSlehmMWLFxcWFipdhWKWLFmSl5endBWKWbZsWW5urtJVKCYxMTEnJ0fpKhSzatWq48ePK12FRT3Qgs3sVpdW69evP3r0KBE58jQ3inc2x4hPGw7C77///vz580pXoZhly5adPn1a6SoU8/PPP6vts6C+lStXGj4L1GnNmjVZWVlKV6GY9evXp6enK12FYjZs2LBv3z7zHtOGgxAAAKDpEIQAAKBqCEIAAFA1Jsuy0jXUOXnyZOfOnYOCgkzc/+LFiz4+PhqNOZZctUGXLl3y8vLSarVKF6KMy5cvt2zZ0sHBOh7iYnF5eXktWrRwcnJSuhBl5Ofnu7q66nQ6pQtRxpUrV5ydnZ2dnZUuRBmFhYWOjo4uLi4m7j9x4sR333333vtYURASUYNGglVXV6v2c5Bw+bh8FV9+TU2NRqNh1viEV0tQ+eXX1tbyPM9xpnZn+vv73/dPRusKQgAAAAvDPUIAAFA1BCEAAKgaghAAAFQNQQgAAKpmkw/mzc/PP3DgwMWLFwcOHBgaGqp0OZZ27NixjRs3Xrp0KSQkZNKkSW5ubkpXZFHbtm1LTU29du1aq1atJk2a5OXlpXRFCpBl+YcffvD19R0yZIjStVjUjh07jEvr8Tz/zDPPKFuP5RUVFX3//fe5ubnBwcGTJ0/29PRUuiLLWbZsWXl5ufHL0NDQgQMHmuXINtki7Nev3/vvvz9z5sz9+/crXYsCBg8efPbs2eDg4I0bN4aHh5eUlChdkUUlJiaKohgSErJnz55u3boVFBQoXZECvvvuu+nTp3/99ddKF2JpS5cu/emnn3JycnJycs6cOaN0OZaWm5vbrVu3AwcOtGnT5uzZs2pbefzcuXM5f/j73/+emppqriPb5PQJSZI4juvevfusWbMef/xxpcuxtKqqKkdHRyISRbFdu3bz588fP3680kUpo0OHDu++++6jjz6qdCEWdfny5fj4+FGjRh07dmzt2rVKl2NRzzzzTFhY2Ouvv650Icp47LHHfHx8vvjiC6ULUdiFCxdCQkJOnjzZunVrsxzQJluEpk+ltEuGFCQixlh1dbWrq6uy9Sjl1KlThYWFHTt2VLoQS3vxxRfnzJnj7u6udCHK2Ldv3/z585cvX15TU6N0LZa2cePGsWPHLlmy5Ouvv1Zhg9ho8eLFcXFx5kpBstEgBIM5c+YEBATEx8crXYil/e1vfwsKCurUqdMHH3zQuXNnpcuxqGXLlgmCkJCQoHQhymjVqpWPj09xcfHcuXOjoqLq3zGye4WFhaWlpS+//PLx48dPnDgRERFx4MABpYtSgCzLS5cuffrpp818UBvVrVu3H3/8UekqFLNkyZKgoKDTp08rXYgCysvLL1++vGLFCk9Pz7S0NKXLsZyCgoJ27dpduHBBluW5c+eOGjVK6YoUU1tb27Vr108//VTpQizn6tWrRPTJJ58Yvnz11VfHjRunbEmK2LJli6enZ1VVlRmPaZOjRuGnn3568803t27dGhISonQtCtDpdDqdbuzYsevXr1+1alWfPn2UrshCNmzYUFRU9MgjjxBRXl5eeXl5XFzctm3blK5LAYIg9OrVS1Xdgx4eHk5OTsZ7AZ06ddq5c6eyJSni22+/ffLJJ8271i6C0PasWrXqlVde2bx5c1hYmNK1WJper5ckyfDMjdra2szMTDP3kFi3ESNGdOrUyfB6yZIlmZmZn3zyibIlWVhlZaVhAeWysrLk5OQ33nhD6YoshzGWkJCwZ8+ewYMHE9GePXtUeIO8uLh4zZo1aWlp5j2sTY4anT59+u7du48dO+bn59eyZcsFCxZERkYqXZSF1NTUuLi4eHl5BQQEGLbMmDFj0qRJylZlMXl5ed26devTp4+rq2tKSkqbNm02bNigzqcRzZs3Ly0tTW2jRn18fHr37u3q6rp9+/Zu3bqtXbtWVQ9iy87OjouLi42NraioSE9PT05Obtu2rdJFWdSXX365ZMkSs0+cs8kgPHny5PXr141ftm/fXj0jJ2VZPnToUP0tQUFBvr6+StVjebm5uenp6VVVVe3atQsPD1e6HMUYukbVtqDEuXPn0tPTq6ur27dv36NHD6XLUUBxcfHWrVudnZ1jYmJMfyaf3Thz5owgCK1atTLvYW0yCAEAAMwF0ycAAEDVEIQAAKBqCEIAAFA1BCEAAKgaghAAAFQNQQgAAKqGIAQAAFVDEAIAgKohCAEAQNUQhAAAoGoIQgAAULX/B6d3KhSZlSoyAAAAAElFTkSuQmCC", "text/html": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ], "image/svg+xml": [ "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", " \n", " \n", " \n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n", "\n" ] }, "metadata": {}, "execution_count": 5 } ], "cell_type": "code", "source": [ "p" ], "metadata": {}, "execution_count": 5 }, { "cell_type": "markdown", "source": [ "The `info` object passed to the callback contains not just the densities\n", "but also the complete Bloch wave (in `ψ`), the `occupation`, band `eigenvalues`\n", "and so on.\n", "See [`src/scf/self_consistent_field.jl`](https://dftk.org/blob/master/src/scf/self_consistent_field.jl#L101)\n", "for all currently available keys.\n", "\n", "!!! tip \"Debugging with callbacks\"\n", " Very handy for debugging SCF algorithms is to employ callbacks\n", " with an `@infiltrate` from [Infiltrator.jl](https://github.com/JuliaDebug/Infiltrator.jl)\n", " to interactively monitor what is happening each SCF step." ], "metadata": {} } ], "nbformat_minor": 3, "metadata": { "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.9.0" }, "kernelspec": { "name": "julia-1.9", "display_name": "Julia 1.9.0", "language": "julia" } }, "nbformat": 4 }