{ "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.774274701735 -0.70 0.80 4.8 \n", " 2 -7.779026246327 -2.32 -1.52 0.80 1.0 30.5ms\n", " 3 -7.779318850934 -3.53 -2.59 0.80 1.5 32.5ms\n", " 4 -7.779350256427 -4.50 -2.93 0.80 2.2 36.0ms\n", " 5 -7.779350722965 -6.33 -3.25 0.80 1.2 31.2ms\n", " 6 -7.779350849834 -6.90 -4.58 0.80 1.0 30.0ms\n", " 7 -7.779350856116 -8.20 -5.31 0.80 2.5 37.5ms\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+gvaeTAAAgAElEQVR4nO3dd0AUZ94H8OeZmaXsgkiVrmIPxUizYaUpRVGj5pJoYsxpyuXVxLtokkte85rkjMmd8VIuIcVozsRYURCVoti7KDYUBRUQBUSlw055/1hD7C6w7Owy389fu+PuzG+WZL/7zFOGSpJEAAAAlIqRuwAAAAA5IQgBAEDREIQAAKBoCEIAAFA0BCEAACgaghAAABQNQQgAAIqGIAQAAEVDEAIAgKIhCAEAQNFMKAirq6vff/99/V+v1WrbrhjTh9OXuwQ58TwvdwlywukreWlMQRAMfvomFITl5eU///yz/q+vr69vu2JMH05f7hLkVF9fr+SvwoaGBlEU5a5CNo2NjUr+6zc2Nhr8r29CQQgAAGB8CEIAAFA0BCEAACiakYJQEIQDBw4cO3ZMyZe2AQDABHHGOcykSZMcHBwqKytdXV2XLFnS4v1UaYmKIVbsvdvL6omzVasqBAAAZTJGEB46dOjGjRtr166VJMnPz6+oqMjT07Nlu9p2RUzMFddGcHdm4Xe54rYS6dcR98UjAADA4xjj0uixY8dCQ0MJIZTSoKCgnJycFu9qbGfmqa7MuHS+7vd5RD+cFZMuiUuHIgUBAKAljBGElZWVarVa91itVldWVrZmb9N6MpN8mPEZfB1Pfjgrrrt4bwMRAABAf211aVQQhBs3bnTs2JHjuE6dOh04cEC3vbS0tFOnTq3c+bSejERI/zRLH1sxNZqzRAoCAEBLtapFmJ6eHhkZ6erqGhMTc+f23bt3d+nSJTAw0NPTc+vWreHh4WlpaTU1NWVlZcePHx8wYEDraiaEEEEkokj2XJPWX1TuAhMAANB6rWoRajSaGTNmnD59OiMjo2mjKIrPP//8ggULXnjhhY0bN06dOrWwsPDdd9+NioqilH755ZfW1tatLFrXL3h4dMOX+eqXdgnJl6Rvh7A2qlbuFQCgVbZt2/bNN9+09VEEQWAYhlLa1gcyNdHR0dOnT2+LPdPWT+z77rvvli9fvmvXLt3T3bt3jxs37urVqyzLEkK6dev2+eefx8fHP3Y/p0+fDg4OvrO9OHny5Gefffaely3PZzcWMf8drOXrqm1sbH48z/7zNEsI+WEQP8BJQa3D6upqGxsbuauQjcJPv6amRq1WK/CrUKe2ttbKyophTG49kH/961/Z2dnPPPOM3IW0Q3v37i0qKlq6dGldXZ2FhYUuX/RhZWXFcY9p8hm+jzA/P7979+5NVfbq1Ss/P1+fN1pbW9vZ2b3zzju6p5RSf3//e77ski6Jm6+KG6I5S9ayipFsbGz+50lioxZ/PCs+t0f1ch/mvX4sq4wvB0mSlJwECj99QohGo1FsEFJKra2tTTAILSwsevfuPXHiRLkLaYdEUVy/fr2NjQ3Lss0KQn0YPgjvHCNKCLGxsbl586Y+b6SUWllZRUREPOI1Ee7MaE/mntExL/ZkxnZmeJG8uJMPS+b/O5zt1kGhXxAAANBchv9J5ezsfGfy3bhxo/XDRJvYqMgDx4g6WpJO1iQlmnvahxm4kf/5vIKukQIAQGsYPgj9/PzOnDlTU1NDCOF5Pjs729/f3+BHeSBKyCw/Znss91mOOClTuNFgnMMCAIAZa1UQlpeXZ2RknDlz5ubNmxkZGbolY3x9fUNDQ+fMmXP+/Pl3333X09Nz8ODBBqpWL7729MBYzk1N+q3nd17FGt8AAPAorQrCS5cuJSYmFhUV9enTJzExMT09Xbd95cqV1dXV48aNKygo2LBhgyHqbB4rliwZyCYOYZ/ZLszaJzTiQikAmLAanjzia0oi5GajEatRnlYNlgkKClq1atX9211dXf/73/+2Zs8GEeVBs8dx03cKYcn8iuFsDzuMoAEAUzQ1SxAksiqctbivbSIRMmufsPuqdDiBY/Ad1jZMbvyxYTlbkY1R7F+eYAYl80tOomEIAKbofwOZPdfEcel8g3DXdomQ2fuExFxxQTDbshT08fEpLS01SJGiKGZkZAiC8PiXPpJWq3VycqqtrX3Ea/r165eXl9fKA+mvnQehztQezM44blmeOCFDuI4RNABgYgIcaGYMd7BMGp/xRxbqUvDbXHFtBBfr1cLGYENDg6Fuh87zfGRkZF1dXSv3w7Ls+++/b2Fh8YjXNDY2iqLxmi6KCEJCSJ+O9MBYzt+BBKzltxRhBA0AmJZ7stAgKXi/VatWzZs3LzExsaHhjzZBZmbmBx988NVXX1VUVOi2ZGVlnThxQve4uLh43bp1hJCkpCRCyNKlSxMTE0tKSvQ8Ym1t7eeff/7uu+8mJyc3bbSyun0j9Y0bNxYUFCxfvvy9997bvn17q8+vhZQShIQQFUPmB7LLh7Mzdgmz9gkNrW3fAwAY0p1Z+Ppew6fgrFmzvv766+7duycnJ8fFxek2fvzxx6+++qqzs/PJkyf79et3/fp1QsiyZcsyMzN1L8jNzf34449bdsT6+vphw4adPn3a29t73rx5H330ESGE5/mZM2fqknjRokXjxo07e/Zsx44dJ06cuG3bNgOcZ/MZ4w71JiXcneZM4F7ZIwQn8b+MYP0d0PsMAG0lo1j6JKd5P7q7daBbiyVJkgIc6ecnhc9PNuO9U3swU7o/uHlTWFi4fPnyS5cudejQ4fnnn+/WrdvOnTsDAwM//PDDQ4cO+fr6EkISEhK++uqr999//4F7SEhIIIRMmzZN/9UN16xZY2trm5iYSAgZOnRoSEjInDlz7lkbLyEhYf78+YSQGzduJCUljRw5Ut+zNRzFBSEhpKMF+XUEuzxPDE/l332S/R8/DMUCgDYR5ETnBjRjVUyJkK9PiwyRVCzRsHSOP8M15+upT8eH/tOJEyd69+7doUMHQohKpQoODj5x4kTHjh0tLS11KUgICQsLO3jwYDOORwghpKys7M033ySEMAyzbNmyO/8pJyen6T4Kffr04TiuoKCgW7dud77mySef1D3w9PTcsWNHc49uEEoMQp2pPZghrvS5LGFzkbh0KOemfvxbAACaxd6SRHjoG2W6fsHNReL6SM5LQ8JT+S9OSesiDHPv8Zqamvr6+qandXV1arVarVbrRtPoFnDXbSSEsCzbNDr00cM7CSEajWbChAmEkPtXga+pqWkaFCOKYkNDw50rUesYdvnsllFQH+H9utrSHbHcABcalKRNuYwRNAAgm3tGxzxwHGkrnThx4uTJk4SQoqKiPXv2hIWFdenSxcXFZc2aNYSQurq6VatWDR8+nBDSuXPno0ePEkIkSVq7dq3u7RYWFmq1umlATRO1Wp2QkJCQkDB27Nj7D7px40bdipvr1693d3f38vIyzMkYlKKDkBDCMWR+ILsugpu9X5i5W6jl5S4IAJTngWNEDZ6F/v7+L7300pgxY0JCQt5+++0ePXpwHPfjjz/Onj179OjR/v7+/v7+U6ZMIYS8+OKLu3btGjZsWHBw8J3TGF5//fWQkJDg4ODs7Gw9D9qrV6/BgwfHxcW98sor3377rQnePIso+dLonQa40KPjuNf2CCFJ/IoR7JOO6DQEAON5c7/w3VkxOYqLvPs6qi4Lw1P5yduEdREtnFPfxNHRMTU19dixY56enu7u7rqNI0eOPH/+/JkzZ5ydnZuaax4eHrm5uadOnfL09HR2dm66Orpw4cIFCxZUV1fb2trqedBBgwa99tpr586d69Onj66H0sLC4tq1axqNhhCyefPmpqkU06ZN08Ww8SEIb+ugIj8PZ1cXiKO38LP92L8FMBhCAwDGMdCFxnjdm4I6uizcUSIZ5BvJwsIiNDT0no3W1taBgYH3bFSr1SEhIbrHugDTUalU9vb2zTqovb19//7979zi4uKie3BnoFpaWlpaWjZrz4aCILzLxK5MqDOdkiVkXBF/Gsp6aBCGANDmJvk86oJhgAMNaOlEL09PT91olN69ez/99NMt20mLRUZG3hmi+nN3d1epVAav52EQhPfqbEN1dzQM3SB8PZgZ29kUr2gDAOjjwIEDuge6LkAjH1039bAFmu5lZBz4ln8AlpK5fZn1kezfDopTs4RqrdwFAQBAm0EQPlSoMz2SwFlzJDiJP1KOyRUAAO0TgvBRbFXk2zD2o2Amdis//6ggIA0BANodBOHjTejKHErgdpRIQ5L5giqEIQBAu4LBMnrx0tBtsdy/T4r9N/CLB7DPPmRZWwBQMo7jfvrpp7Ye6NG0Ipqi3Lhxo2nZUoNDEOqLEjLLjxnpTp/NEjZckhKHsB0fdV9JAFCcmTNnjhgxoq2PUl9fb2FhYZpLtLQpT0/PNtozgrB5/B3ogTHcvENCv/X8z8PYMFfF/S4DgIextbUNCgpq66PU1tZaWVkpMAjbDj7KZrPmyJKB7DeD2T9tF+YdErTi498CAAAmC0HYQtGe9EgCd7JCGpLCn6/ECBoAAHOFIGw5F2uSHM292JMJS+YTc9EwBAAwSwjCVqGEzOjNbIvlvj4tTswUKhrkLggAAJoJQWgAT3Sk+8dy3TqQwPX8jhJcJgUAMCcIQsOwYsnCEPb7IexzWcKsfUIjLpQCAJgJBKEhRXjQo+O4i9Vk8Eb+3C00DQEAzACC0MCcrciGSPZ1X2ZwMr/kJBqGAACmDkHYJqb2YHbFccvPi+MzhOsYQQMAYMIQhG2ld0e6fwwX4EAC1vKphbhMCgBgohCEbUjFkPmB7H9HsC/vFmbtExoEuQsCAID7IAjb3Ag3emICV15PgpL4nAo0DQEATAuC0BjsLMiKEexbAUxEKr/kpIgwBAAwHQhC45nagzk4lltdII7azJfUyl0NAAAQQhCERtbFlmbFcgM70aAkbfJlTK4AAJAfgtDYOIbMD2TXRXBv7BenZgk1vNwFAQAoG4JQHgNc6NFxnDVHQpL47OvoNAQAkA2CUDYdVOTbMPaDIGbUFn7+UQFDaAAAZIEglNnErsyhsdz2K1LUZr64BmEIAGBsCEL5edvQ7bFcvDfTbz2/8gJG0AAAGBWC0CQwlMzyYzZFc+8fFadmCdVauQsCAFAMBKEJCXGmx8dx9pYkYB2/9xoukwIAGAOC0LRYc2TJQPaz/sz4DH7+UUFAGgIAtDEEoSka34U5Nl51sEwaksznVyEMAQDaEILQRLlak03R3GQfZsAG/r/nMYIGAKCtIAhNFyVklh+zLZb7NEeclCncwA1+AQDaAILQ1PnZ0wNjOTc1CUziF+WI35+9t3V4pVb6y17c6hAAoIU4uQuAx7NiyZKBbIyXNH2XoOZIHU9e9739C+ZqHZmQISwewMpbIQCA+UKL0GxEe9KjCVwPO/K/R8X5RwRCyNU6Mi6dXzyAHeBC5a4OAMBcoUVoTlysSUoU958z4pv7haNlFmWNSEEAgNZCi9DMUEJe7cNsiGS3XGF4ifRHCgIAtA6C0PxcrSPzj4q/DdEWVJL+G3gt5lYAALQCgtDMNPULRrgKBU9zRTUkcD1fh7v7AgC0FILQnNwzOsZWRS5M5iq1JGAdf6tR7uIAAMwTgtCcXKqS/nX36BhrlpydyHW2JWHJ/JVaLMYGANBsCEJz0t+FDrxvdIwVSzJGc891Z4YkCxcqkYUAAM2DIGwn5vZl5vgzwzYJORXIQgCAZkAQth+vPsEsCmUiN/P7SpGFAAD6MtKE+s2bN2dlZWk0mvfff984R1SmZ7oxdhZ0bBr/83Au2hNTDAEAHs9ILcJLly516dJlzZo1xjmcksV60eQobtpOflU+JhgCADyekYLw5Zdfjo2NNc6xoL8LTR/N/fWAmJiLLAQAeAwDXxoVBEGS7uqg4jgsZyoDX3u6I46N3iJcqZXmB+LeFAAAD2XglIqLi6uoqLhzy8aNGzt16mTYo4A+utrSnXHcqM38jQbh84EsOgwBAB5I3yCsrKw8evRoXl5eeHi4j49P0/YrV66sWLGirq5uwoQJvr6+mzdvbps6oSVcrcn2WC4+jX9hh/DDEJbDGGEAgPvo+9Xo5+c3e/bsOXPmHDlypGljaWlpUFDQhQsXRFEcNGjQwYMHH/b2jIyMX375paKiIjEx8dy5c62tGvRmb0nSRnOlddKETKEe97EHALiPvi3CCxcuqFQqf3//OzcmJiaGhoZ+8803hBCWZT/55JO1a9c+bA8ODg6PnTuh1WrvDNqePXva2trqWSE8jJojG6K4qVnC6C38hiiug0ruggAATIm+QahSPeDrc/v27RMnTtQ9HjVq1D//+c+HvT0iIuKxh6itrS0vL3/ppZeatjz//PPTp09/2Ourq6sfu892rLmn/20wefOIathGYd0wrZOl2c+4V/hfv6amRhRFShXa81tbW8vzPMMo9Fp/XV2dVqtV8ulbWFiwrL5jAK2srB6YX3dq1WCZkpISFxcX3WNXV9dbt27V1taq1eqW7U2tVru5uWVnZ+v/FoW3F5t7+j+MIJ8cF2Oy2LTRrJfG7L9DlfzXp5RqNBrFBiHDMNbW1opNApZlraysFHv6HMc1Kwj10aqPkmEYUbw9U00QBEqpYv825mJuX+bVPszQFOHsLbNvFAIAGESrWoTu7u4lJSW6x1euXHFwcLCysjJEVdCGXvdl7C3JsBR+QyTX/757WQAAKE2rGnAxMTHr16/XzaBfv3491o4xF891Z74bwsan8enFaBcCgNLpG4Rz586NjIy8ePHiggULIiMjT5w4QQiZNm3a1atX4+LiXnjhhZ9++mnevHltWSoYUrw3szqce3Y7v7YAy7ABgKLpe2l08uTJkZGRTU+9vLwIIXZ2docOHUpNTa2vr1+0aFHTwBkwC8Pc6OZRXHwaf7ORTO+Fzl0AUCh9gzAwMPCB2zUaTdMMCjA7QU50dzwXtVm40Uj+6o8sBAAlwnef0vnY0l1x3PI8cd4hLDwDAEqEIATipiZZsdzOEumVPYKI0TMAoDAIQiCEEAdLkh7D5VdKEzKwJCkAKAuCEG7TcCQ5muMYEruVr9LKXQ0AgLEgCOEPFgxZOZLt1oGGp/Ll9XJXAwBgFAhCuAtLybdh7Ag3OiyFL6pBhyEAtH8IQrgXJeSTUHZGb2ZoipCHJUkBoL1r1Vqj0I7N8mPsLcmIVCElin3SEUuSAkC7hRYhPNTUHsxXg5jRW/jdV9EuBIB2C0EIjzK2M/PLSO6pTD61EFkIAO0TghAeY4QbTYnmpu/kfzqH5bkBoB1CHyE8XrAT3RbLjdoi3Ggkb/jhxxMAtCv4UgO99OlId8ax35zBkqQA0N4gCEFfnW3ozjgurUj6y14sSQoA7QeCEJqhkzXZGcedvSVNyRK06DEEgHYBQQjNY6MiKVFcvUDGpfN1vNzVAAC0GoIQms2SJavC2U7WdEQqX9EgdzUAAK2DIISWYCn5figb1okOS+Gv1KLDEADMGIIQWogS8ll/9rnuzJBk4XwlshAAzBWCEFplbl/mrwHMsBQhpwJZCABmCUEIrfVKH+az/kzkZn7vNWQhAJgfBCEYwJ+6Mb+M4MZn8FuKkIUAYGYQhGAY4e50YxT34k5+VT4mGAKAOcFao2Awoc40fTQXs1W40Uhm9sZvLAAwD/i2AkPytadZsexnOViSFADMBoIQDKyrLd0Vz20plGbtw5KkAGAGEIRgeK7WJCuOO1IuvbBD4NFjCACmDUEIbaKjBUkbzZU3SOMzBCxJCgCmDEEIbUXNkaRITqMio7fylVq5qwEAeAgEIbQhC4b8MoLt50hHbuLL6uWuBgDgQRCE0LYoIYsHsBO7MkNT+MvVGDwDACYHQQjGMLcv89oTzKBk4eQNZCEAmBYEIRjJX55gPglhIlL5/aXIQgAwIQhCMJ5nuzPfD+HGpPHpxchCADAVCEIwqjhvuiaCe3Y7v6YAEwwBwCQgCMHYhrrSLaO4/9knfH8WWQgA8kMQggwCnejueO6T4+KnOchCAJAZghDk4WNLd8ZxP58X5x0S0GEIADJCEIJs3NRkRyy3s0SahiVJAUA+CEKQk70lSY/hSmqliZlCPW7cBAByQBCCzDQcSY7mVAyJ3cpXYUlSADA6BCHIz4Ihv45ku3Wg4al8OZYkBQDjQhCCSWAp+TaMjfGiQ1P4ohqMngEA40EQgqmghMwPZGf2ZoamCOduIQsBwEg4uQsAuMssP8bekoxMFVKi2CcdqdzlAED7hxYhmJypPZivBjGjt/C7rqJdCABtDkEIpmhsZ+bXkdyEDH79RUwwBIC2hSAEEzXcjaaO4l7dIyw9hywEgDaEPkIwXcFOdHssF71FuNFA3vTHjzYAaBP4cgGT1rsj3RnHJuaK8w5h4RkAaBMIQjB1nW3o3jFcVon0yh6h5kFLz9xsNHpNANCOIAjBDDhYkozR3MEyqfsq7a27Y2/jJfHpbbxMdQFAe4AgBPNgoyJ747nOtrTbb9rSutsbUwulT3LEX0egqxsAWg5BCGbDkiV74rkAB9pjlbawlkktlD46JqREcfaWclcGAOYMP6XBnLCUZMZyMVv44FSVTwdhWwxSEABaCy1CMDOUkFf6MPaWpKhG6rFaOz5D+C5XLMY63QDQUmgRgpnR9Qvuj27cd0v9xSlxXBe6+6r0/hHBwZLGd6YR7swwN6rCDzwA0JuRglCSpMLCQpVK5ebmZpwjQrvU1C/INUqTfBhCyDdnxI1R3NeD2WPXpeTL4rxDwqVqaYQbE+FBx3RmXK3lrhgATJ4xgrChoSEoKKhr165VVVV2dnZr1qxRqVRGOC60MztKpIXHhdRozs6CVDUSQsgkH6ZRJH/azidHcUFONMiJnR9ISuvIliIx5bI096C2Wwca503jvZlAJ4o7WQDAAxkjCFUq1Z49e+zs7AghkZGRO3fuDA8PN8JxoZ0JcaYp0VyHu39EPdedCXe/60qoizWZ2oOZ2oPwIru/VEopFKdkCRUNUpQHE9+ZRnsyHfAzDADuYIwgZBhGl4KiKF6/ft3FxcUIB4X2R/2Q/1rd1A/ezjEkzJWGubILQ0h+lZRRLC3PE/+8SwhxphHuzJjOtE9HtBIBwNBBuGzZsuvXr9+5ZcqUKc7OzrrHH3zwQWhoqL+/v2EPCvBYPrZ0Rm86ozdTy5PMK2LKZSlqs2jBkAgPGudNozwYS1buEgFAJgYOQicnp3v6/zju9iH+9a9/nTp1auXKlYY9IkCzqDkS783EexNCyKkbUspl6d+nxOe2C0PdaLw3E+tFPTRoJgIoi75B+N1336Wnp+fn57/99tsTJkxo2v7jjz8uWLCgqqpqzJgxX331VWxs7MPenpmZuW7duqZcBJCdrz31tadz+zLl9WR7iZh8SZp3SHCzvj0NY7gb5TANA0AB9P0fvbS0NCIiguf5srKypo05OTlvvvnmqlWrLly4cPHixQ8//PCB762oqJg5c+bly5cHDx4cHBycnJxsgMIBDMfJikzsyiwfzpY9p1o+nLViybxDgtsv2kmZwvI88UaD3PUBQFuiktSMJTlGjBgxefLkl19+Wff0jTfeqKys/OGHHwghmZmZzz33XElJSYtLyc3NDQ0NnTRpUtOWmJiYmJiYh72+qqrK1ta2xYczdzj9tj79yzUk/QrdVkIzS4iPLYnxlGI9yZMOkilcOa2urtZoNFSpU0Jqamqsra0ZRqEN9traWisrK8Wefl1dnYWFBcvq26uvUqke++JWXag8e/ZsdHS07nFAQMDVq1dv3rzZsWPHlu1NpVJZWloGBwc3benWrdsjZhyqVColz0fE6bf16XfrSLp1JC8/QeoFsucaySyRpu0h9QKJ8iDh7iTag9rK9/HrTl+xQag7fcUmgcJPn+d5fbKtiT4fVKuCsKKioulXeYcOHQgh169fb3EQsixrY2PT1NzU5/X6fxbtD07faKevYUmUF4nyIp+EkvwqKfmS9P05ccZuMdSFRrgzCV1oLztjB5Lu9BUbhLrTV2wS4PQN/r9/q4LQwcGhsrJS9/jWrVuEECcnJwMUBWCqfGzpLD86y4+p4cm2K2LKZWnkJtGKJbr1a4a6UQuFfjsBmLFWBWHPnj1PnTqle3zy5EkXFxfdxHmAdk9z3zSM+UeFnAppmBuN92bivKm7WqHNNQCzo28QlpaWVldX19XVlZeX5+fnd+rUSaPRTJs2bfjw4a+99lr37t0XLlz44osvtmmtAKapaRpGWT3JKhGTL0lvHxK62tIIDxrnxQzqRBlkIoAJ0/c6zuLFiydNmsTzfFJS0qRJkw4ePEgI6du376JFi8aMGePl5eXq6vree++1ZakAps7592kYpc+pPh/AEkJm7xdcV9yehnGzUe76AOBBmjd9ok1dvHhxxIgRBQUFer4e8wdw+nJXoZeCKim9WMooltKKRT97Gu/NRHjQIKdWNRIxfQLTJxR7+s2dPqEPrPMC0La63l7mlNTx7J5rUvJlcXyGyFIS6UEjPOhoT8ZGubNgAEwCghDASKw5EuFBIzzYJQNvT8NIzBVf2imEutA4L2ZcF+pto9AWHoC8EIQAMmiahlHRQDKviBnF0icbBWv29m2Eh7lR1X3XvW41kmqtdP+a4KdvSk/gflIAraDQq8wAJsLBkkzsynwbxhb9SbUqnHVX0/lHBbcV2kmZQmKuWFL7xyvzq6Sx6UJRzV2d+v85I757SDSVfn4A84QWIYBJYCgJcqJBTnRuX6a0jmwpElMuS/MOaX1+n4Yx2JUmhrFj0oR1EaxuqM13ueKmy+LaCA7tQYDWQBACmBwXazK1BzO1B+FFdn+plFIozt4vFNZI0R7M0z50bLrw62C6u1DceFlaG8HhlsIArYQgBDBdHEPCXGmYK0tCSN4taVOhlFoonr8lhaSqvGyk1eEsUhCg9RCEAOahhx2dbUdn+zGfnxQTz/DlDVLMVp5j6ChPOsqThrtjGgZAC2GwDIA5+S5XzCgWd0U1bo5mna3osmFMP0f6/VnRdYU2OImff1Q4Um4ya2QAmAm0CAHMxne54oZL4toITltHAp3o90PYl3YJ6yLYGb25Op7suSZlXBGnZgnXG6ShroxuJoa9pdxFA83SruwAABz3SURBVJg8tAgBzMOOEmlTobQu8o/RMYFO9D+D2em7BOn32foLQ9hTT3F7x3ARHjTlsuTzmzY4iZ93SNh9VcIcC4CHQYsQwDwMdaODOrH3TLTv70I3j7p3+oTP74u61Qvs7qtSxhVx9n7hcrU03I2J8KBjOjOu1kasG8DkIQgBzAMl5P7lZgghj7gVsBV7e1E3EkIuVklpxVJG8R9zEyPcmeFulMNVIVA8BCGAInT5vZnYNDdx3iHhUrU0wo2J8KCxXvT+xdsAFAJBCKAsd85NvFZHthaJKZeluQcFdzWN70wj3B+80ilAO4YgBFCuTr8vYSNI7LHrUkaxNP+okFMhDXOj8d7MaC/qhWYiKACCEAAIe8dKp2X1JKtEzCiW/veIoOZu3xBjqBt9RGckgFlDEALAXZytyMSuzMSu5D+D2ezrUkax9EmO8FSmFOJM47yYhC60M+6bCO0LghAAHuzOG2JcbyDbrogZxdKijaJuMGqcN430YKyw2CmYPwQhADyeo+XtZiIh5NQNKeWy9O9T4pQsIcSZRrgzYzrTPrg5MJgtBCEANI+vPfW1p3P7MtVasr1ETLksRW8WVYxuziId5cnYYvlvMCsIQgBoIRsVifdm4r0JISS/Skq+JCXmin/edbuZGOFBg5zQTAQzgCAEAAPwsaWz/OgsP6aGJ/uuScmXxQkZIkNJpAeN8KBRHoydhdwlAjwEghAADEnD3V7XbclAkl8lZRRLy/PEaTuEJx1pvDcT4UEDnSjaiWBSEIQA0FZ+X/6bqeXJ3mtSxhVxSpZY0SBFeTDxnWmkB9MRzUQwAQhCAGhz6t+biQtDbjcTV+dLM3Zpu3W4PWEfzUSQEYIQAIyq6S5RdTyru5nwzN237xIV503jvBmHB91MuEEglVribHXv9qIayRPrwEHrYNEkAJBH082EDydwB8bevplw91UPvpnwmZtS7Fb+Wt1de1h/UZyaJQi45zC0DlqEACC/rvfdJWr2/j/uEhXnTZ90pIsHsHFb+Q1RrB0hhJBNhdK/T4kbojgWDUJoHQQhAJiQR98l6rnuzJg0YWUYyb1O/nlC2BDFdcDkfWg1BCEAmKimu0RpRXbvNWlLkfhTnphfKQVtUvl0EHbEqZCCYBDoIwQAU6diyDA3+o8QNnsct3gAq+HImZvkcrXcZUF7gSAEALOxqVD6+bx4NI5/tjszYIP2fCXGyYABIAgBwDxsKpQ+yxE2RHG2nPTdYBrpyQQn8YU1yEJoLQQhAJiBQ2XSpzl3jY5ZNZLt2YGGbuBFRCG0DoIQAMxAkBPdFH3XGFEVQzJjOWdLuuSUKF9d0B4gCAHADDCUaO4b5G6rIltGs0tOimsLkIXQcghCADBj7mqaEs3+Za+w9xqukEILIQgBwLz52dOlw7gJGfy5W8hCaAkEIQCYvVGe9MNgNmarUFr3+BcD3ANBCADtwfRezGQfGp/G1/JylwLmBkEIAO3Eh8FsLzv6wg4BEyqgWRCEANBOUEJ+GMpWNEjzDgly1wLmBEEIAO2HiiFrI7gtRdIXmFwIekMQAkC7YmdBNkayC4+LSZeQhaAXBCEAtDddbOnGKHbGLmF/KXoL4fEQhADQDgU50aXDuHHpPO5QAY+FIASA9inWi34QxI5JE240yF0KmDYEIQC0WzN6MzFeNCGdb8AwUng4BCEAtGef9mc9NfSFnQKukMLDIAgBoD3TTS4srJbeP4JWITwYghAA2jkrlmyI4lblS/85gwkV8AAIQgBo/xwtyeZR7EfHxK1FuEQK90IQAoAi+NjS9RHs8zv4Y9eRhXAXBCEAKEWIM/1uCDsmTbhcjSyEPyAIAUBB4r2ZWX5MzFbhZqPcpYDJQBACgLLM8WfC3en4dL4RQ2eAEIIgBAAFWjyA7WhJX9mNCRVACCGcEY4hCMLXX3+dl5dna2s7derUXr16GeGgAAAPw1CyYjg7MpX/MFv8ez+0B5TOGP8FiKJoZWX1pz/9KSAgIDIy8tatW0Y4KADAI1hzJDmKW35eXJ6HK6RKZ4wgVKlUf/7znwcOHDh58mQnJ6eysjIjHBQA4NGcrMjGSPatg0JGMQaRKpqRrglotdqZM2dGRETExMR0797dOAcFAHi03h3p6nDume18TgWyULkM2UfY2Ng4cODAezYeOXKEEMKy7IwZMy5evPj++++/+uqr7u7uBjwuAECLDXGlXwxi47YKe8ewnhoqdzkgg2YE4YULFwoKCkJCQuzs7Jo23rp1a8uWLQzDjBo1ytbWVhd792MYJigoKCgoaNOmTXv37n3qqadaWzgAgIFM9mEuVJKx6cKOWM5GJXc1YHR6XRrVarX29vahoaGjRo3Kzc1t2l5cXOzr67t69eoVK1YEBASUlpY+8O35+fkbN248e/ZsSkpKVlZWSEiIYWoHADCQd55kBrjQSdt4HkNnlEevIOQ47ujRo9evX9doNHdu//zzz4cPH75mzZqkpKSgoKCvvvrqgW/XaDQHDhz44IMPtm3blpSU1Llz54cdSJKkG3eQJFy1BwAj+fdAVsXQV/dicqHi0GaFjZ2dXVpaWv/+/XVP+/Tp849//CMhIYEQ8uuvv3722WcPuzSqj9OnTwcEBNja2jZtmTNnzqxZsx72+urqahsbmxYfztzh9JV8+jU1NWq1mlKFdmjV1tZaWVkxjOHH+lXzdFSmamJncVZv3uA7N5S6ujpLS8u2OH2zUFdXZ2FhwbKsnq+3srJSqR5zvbtVg2WKi4s9PDx0jz08PIqLi1uzN7Va7eXlVVBQoP9b7kxNBcLpy12CbCilGo1GsUHIMIy1tXVbJIEtIVtiyMCNfFd7y2e7m2jSsCzbRr8DzALHcc0KQn206qMUBKHpf0WWZXnedH9DAQDow01NUkexbx4Qtpega0YpWhWEbm5uTbPjr127hkkRANAOPNGR/jaSe3Y7f/YWslARWhWEw4cP37p1q+5xWlra8OHDDVARAIDchrvRz/qzMVuEa3VylwJtT98+wvfee6+0tLS+vn7hwoUuLi4ff/yxo6PjG2+8MWjQII1GIwjCqlWrDh482Ka1AgAYzTPdmLM3pbitfFYcpzHG7QlANvr+ef38/G7duhUUFKR7amFhQQjx9fU9dOjQb7/9Rik9fPiwj49PW5UJAGB084PYyzXC09v4pEiOVejIJEVo3vSJNnXx4sURI0boP2q0qqpKyeMGcfpKPv3q6moljxqtqalpo1Gj99OKJHYr72tPFw8w5DDF1mi72SNmobnTJ/Sh0I8SAEAfKoasieC2XZGWnMSSM+0WghAA4FE6qEhqNLv4pLjuIrKwfUIQAgA8hoeGJkexr+0R9pWaSl8SGBCCEADg8fwd6I9DufHpfB4mF7Y7CEIAAL2M9qILgtnRW4WyerlLAYNCEAIA6OulXsxTXen4DL4e96hoRxCEAADN8I8QtqsNnZoliLhE2l4gCAEAmoES8v1Q9nqD9M5htArbCQQhAEDzWDBkdTi34ZL01WlMqGgPEIQAAM3mYElSo9mPj4kbLiELzR6CEACgJbra0g1R7J93CQfL0Fto3hCEAAAtFOxElw7jJmQIl6qRhWYMQQgA0HKxXvTv/ZiYLcKNBrlLgZZCEAIAtMrM3kyUJx2XwTdgGKl5QhACALTWP/uzjpb0hZ0CrpCaIwQhAEBrMZT8dzh7uVqafwStQvODIAQAMABrjmyM4lbmS9+cwYQKM4MgBAAwDEdLkhrNLsgW04pxidScIAgBAAymWwe6PpKdmsUfr0AWmg0EIQCAIYU60y8GsXFbhcIaZKF5QBACABjYxK7M//gyMVuEW41ylwJ6QBACABje3wKY4W508jaex9AZk4cgBABoE58PZK1Y+vIeTKgwdQhCAIA2wVLyywj21A3p42NoFZo0BCEAQFtRc2RjFLf0nPjzeWSh6UIQAgC0IWcrsjGK/dsBIfMKBpGaKAQhAEDb6tORrgrnntnOn8DkQpOEIAQAaHNDXemSAezYdOFqndylwH0QhAAAxvB0N+bFnkzcVr5aK3cpcDcEIQCAkfy9HxPqjMmFJgdBCABgPEsGsoJEXtuLyYUmBEEIAGA8KoasDucOlkn/PIFWoalAEAIAGJWtimyKZr84Jf5yAVloEhCEAADG5q6mm6LZOfuFPdcwoUJ+CEIAABn42tNfRnITM/lzt5CFMkMQAgDIY4Qb/TCYHb1FKMXkQlkhCAEAZPNiT+aZ7jQuja/l5S5FwRCEAABy+r8gtk9H+vQ2QcAlUpkgCAEA5EQJ+X4IWydIbx3E5EJ5IAgBAGSmYsjaCC69WPr3KUyokAGCEABAfh1UJDWa/ecJcf1FZKGxIQgBAEyCp4auj2Bn7hb2l6K30KgQhAAApiLQiS4bxo1L589XIguNB0EIAGBCRnvR/wtix6QJFQ1yl6IYCEIAANPy595MnDdNSOcbMIzUKBCEAAAm55NQ1tuGPr9DEHGJtO0hCAEATA4l5MehbGmd9PfDaBW2OQQhAIApsmDI6ghu3UXp69OYUNG2EIQAACbK0ZJsHsV+dEzceAlZ2IYQhAAApqurLd0Qxc7YLWRfR29hW0EQAgCYtGAn+r+BbESqcLn6rizUiuT9IwKPtmKrIQgBAEzd9F6Mh5r03yDcbLy9RSuSp7cJzlaUw7d4q+EjBAAwdRYMOTyO03AkNIlvFG+n4HA3+rovvsMNAB8iAIAZsGDIyae4Gp6MTOP+tF1EChoQPkcAAPNgxZJj47ncW+RwOXmqK769DQYfJQCAedCK5OXdwt8DRA1Heq3WLjkpYt0Zg0AQAgCYgaZ+wdl9xOxxTIgzs+SUODSFP30TYdhaCEIAAFN3z+gYC4ZsHsX2daDeNnRoMj/vkIDluVvDeEGYn5/v7u6+cuVKox0RAKB9OHdLivS4a3SMBUN+C2e9bUj2ONWFSuK/jt92BU3DFuKMcxhJkt56663+/fs3NOAWWwAAzeNrT33t6T0bLRiyMIQlhKwOZ5Mvi9N2CsNc6eKBrKOlHCWaMyO1CL/99tvw8HAvLy/jHA4AQFHivZmc8Zy9JQlYyy/Pw2IzzWPIFmF5efmePXvu3OLo6BgWFlZcXJyUlJSamjp79mwDHg4AAJrYWZAlA9nJPtKM3cLqAvGrQay3zb2NSHggQwZhdXV1Tk7OnVu8vLzCwsLeeuutQYMGrV27Ni8vTxTFYcOGdenSxYDHBQAAnUGdaPY47l8nxH7r+bcC2L8GMCzS8HH0CsKKioqVK1cePXq0pqbm119/bdre2Ng4f/78tLQ0Jyen+fPnDxgw4L333rv/7fHx8ZcuXcrPz6+srCwrK6uqqjJY+QAAcDcVQ+b2ZSZ0pS/vFlYViIlhbJATwvBR9ArCgoKC3bt3d+rUadmyZXcG4QcffLBjx46lS5ceOHBg9OjR58+fd3R0vP/tTz/9tO5BcXFxUFCQv7+/QUoHAICH6d6BpsdwP+eJcVv5ST7MxyGsxkiDI80PlSR9R9yeOHEiMDBQq9Xqnmq1Wnd39/Xr14eFhRFCoqKiYmJiHt0LWF5ebmlpaWtr+8B/PX/+fGho6F//+temLcOGDQsNDX3Y3qqqqh62KyXA6Sv59KurqzUaDaUK/ZlfU1NjbW3NMAqdBl1bW2tlZaX/6V+rI/MOk31l5MsBJMK9TUszhrq6OgsLC5Zl9Xw9y7KP/axa/guhpKSkvLy8KahCQ0OPHz/+6Lc4OTk94l9FURRF8caNG01bampqRPGhw590r29Oye0KTh+nr9ggxF+/WafvbEl+GEy2FNOX99IQJ/Lv/pKjpRnPOGzuf/z6/GJoeRCWlpaq1WoLCwvdU3t7++zs7BbvjRBiYWFhb2//6aef6vn6xsZGS0vlzpfB6Sv59LVaraWlpWKDkOd5S0tLxbYIBUFowemP9SGR3uT/soUnN4ofBnMzepvrpyeKYrNahPpo+WdhZ2dXX18vCLcX9qmuru7YsaOBqgIAAANTc2RhCLt1NPddrjhiE593y4zbhYbV8iB0d3dnWfbChQu6p3l5eZgUAQBg4vo50n1juITOzOBk/pPjIq/ca8x/aHkQajSahISEL774ghBSUFCQkpLyzDPPGK4wAABoExxDZvkx+8dy266IwUn8oTKlNw316iMsLy/v2bMnIcTW1tbBwcHNze3UqVOEkE8//TQhIcHb27uqquqdd97x9fVt22IBAMBAfGzp1tHc6gJxbDo/sSvzUTBro5K7JpnoFYROTk4VFRX3b+/cuXN2dnZJSYm9vb2VlZWhawMAgLY1sSsT4c7MOyQErOO/HsyO8lTiCCwDTLB0c3Nr/U4AAEAW9pbk2zB2R4k0Y7fQ14F+OYh1sZa7JuMy1xG0AABgQMPc6LFx3BP2xH+dNjFXVFS3IYIQAAAIIcSaI/MD2fTR3A9nxeEpfO5NpaQhghAAAP4Q4ED3jeGe7c4MTeHnHxUaFTC/AkEIAAB3YSiZ0Zs5Np7LqSDBSfz+0nbeNEQQAgDAA7ir6boI9qNgZlKmMHO3UKmVu6A2gyAEAICHivdmciZwVix5Yg2/tqB9XidFEAIAwKN0tCBLBrK/jmDfOyLGp/FFNe3tSimCEAAAHm+IKz0+ngvrxAQl8UtOiu1pggWCEAAA9KJiyNy+zK44LumSODSFP91e5lcgCAEAoBl62tFtsdzUHszQZH7eIaFBkLugVkMQAgBA81BCZvRmTkxQ5VcSv7X8tivm3TREEAIAQEu4qcmqcPZfA5hpO4WpWcL1BrkLaikEIQAAtFy8N5MznrO3JAFr+eV5Zjm/AkEIAACtYmdBlgxkV4ezi3LEuK385Wozu1KKIAQAAAMY1Ilmj+OGuDL91vOfHBcF80lDBCEAABiGbn7FwbFcerEYksQfKTePMEQQAgCAIXXrQNNjuNl+TNxWftY+oYaXu6DHQRACAICBUUKm9mCOjVfdaCB91/HpxSbdNEQQAgBAm+hkTZYPZ/89kP3zLmFSplBeL3dBD4EgBACANhTjRU8/xfl0IE+s0SbmmuL8CgQhAAC0LTVHFoawaaO573LFEZv4c7dM60opghAAAIzhSUe6bwyX0JkJS+bnHxW0JtM4RBACAICRcAyZ5cccGMvtL5VCkviDZSbRNEQQAgCAUXW1pVtGce/2Y8ak8TN3C9VametBEAIAgAwmdmXOPKUihASs47cUydk05GQ8NgAAKJm9Jfk2jN1RIs3cLQQ40C8HsS7WMpSBFiEAAMhpmBvNHsc9YU/812kTc0Xjtw0RhAAAIDNrjswPZNNHcz+eE4el8Lk3jZqGCEIAADAJAQ50bzz3XHdmaAo//6jQaKz5FQhCAAAwFQwlM3ozx8ZzJypI0Hp+f6kxmoYIQgAAMC3uaro2gv04hJmUKczcLRy5Li3Ivrd5WC+Q1/YKDYIBDocgBAAAUxTvzeRM4KxYMjZN2HNNfOvgH6HXKJJJmYK/PbVkDXAgBCEAAJiojhZkyUB25Qj2cjVZXSC+skcghDSK5KkMIcaLvtzHMBGGeYQAAGDSwlzp0XHcx8eEz06IZ2+yGk6K9WYMlYIELUIAADB9Viz5vyB2Vzx3oIzeaJAMmIIEQQgAAGahUSQfHBH/ESgGONI7+wtbD0EIAACmrqlf8M89xH8PoDU8MWAWIggBAMCk3TM6hhLy5SDWgFmIIAQAAJNWXCNN9rlrjKguCztZU6XPI0xMTLx586bcVcjmhx9+KC8vl7sK2Sxbtuzq1atyVyGbFStWFBYWyl2FbFatWpWfny93FbJZt27d2bNn5a7CqLra0me7306r5OTkU6dOEUIoIXP8GaXPI/zpp58uX74sdxWyWbFixYULF+SuQja//fab0r4L7rR27Vrdd4EyJSUl5eTkyF2FbJKTk7Ozs+WuQjapqakHDx407D7NOAgBAABaD0EIAACKhiAEAABFo5Jk/LsBP1heXp6fn5+np6eery8uLnZxcVGpVG1alcm6cuWKk5OThYWF3IXIo6SkxMHBwdLSUu5C5HH16lU7Oztra2u5C5HHtWvXbG1t1Wq13IXIo7S0VKPRaDQauQuRR3l5uZWVlY2NjZ6vf+aZZxYsWPDo15hQEBJCmjUSrKGhQbHfgwSnj9NX8Ok3NjaqVCpKqdyFyEPhp6/ValmWZRh9L2e6ubk99iejaQUhAACAkaGPEAAAFA1BCAAAioYgBAAARUMQAgCAopnlHeqvXbt2+PDh4uLi8PDwbt26yV2OsZ0+fXrz5s1Xrlzx8fGZMmVKhw4d5K7IqLZt27Znz54bN254eXlNmTLFyclJ7opkIEnSzz//3KlTp+joaLlrMaodO3Y0La3Hsuz06dPlrcf4Kioqfvrpp8LCQm9v76lTpzo6OspdkfGsWLGipqam6Wm3bt3Cw8MNsmezbBEOHTr0448/njt37qFDh+SuRQZRUVEXL1709vbevHlzYGDgrVu35K7IqFatWiUIgo+Pz/79+/v27VtWViZ3RTJYunTp66+//vXXX8tdiLEtX7585cqV+fn5+fn5BQUFcpdjbIWFhX379j18+HCXLl0uXryotJXHL126lP+7v//973v27DHUns1y+oQoigzDPPnkk/PmzXv66aflLsfY6uvrraysCCGCIHTv3n3RokUTJ06Uuyh59OrVa8GCBZMmTZK7EKMqKSmJiIgYM2bM6dOnN2zYIHc5RjV9+vTevXv/7W9/k7sQeUyePNnFxeWLL76QuxCZFRUV+fj45OXlde7c2SA7NMsWof5TKdslXQoSQiilDQ0Ntra28tYjl/Pnz5eXlz/xxBNyF2Jsr7766kcffdSxY0e5C5HHwYMHFy1atHr16sbGRrlrMbbNmzdPmDBh2bJlX3/9tQIbxE1++OGHkSNHGioFiZkGIeh89NFH7u7uERERchdibO+8846np6evr+8//vEPPz8/ucsxqhUrVnAcl5CQIHch8vDy8nJxcbl58+bChQtDQkLu7DFq98rLy6uqqt54442zZ8+eO3cuKCjo8OHDchclA0mSli9f/uKLLxp4p2aqb9++v/76q9xVyGbZsmWenp4XLlyQuxAZ1NTUlJSUrFmzxtHRce/evXKXYzxlZWXdu3cvKiqSJGnhwoVjxoyRuyLZaLXagICAzz//XO5CjOf69euEkMWLF+uezpkz56mnnpK3JFmkp6c7OjrW19cbcJ9mOWoUVq5c+fbbb2dmZvr4+MhdiwzUarVarZ4wYUJycvK6desGDhwod0VGkpqaWlFRMXbsWELI1atXa2pqRo4cuW3bNrnrkgHHcf3791fU5UF7e3tra+umvgBfX9+dO3fKW5Isfvzxx+eee86wa+0iCM3PunXr3nzzzbS0tN69e8tdi7HxPC+Kou6eG1qt9vjx4wa+QmLaYmNjfX19dY+XLVt2/PjxxYsXy1uSkdXV1ekWUK6urt6+fftbb70ld0XGQylNSEjYv39/VFQUIWT//v0K7CC/efNmUlLS3r17Dbtbsxw1+vrrr+/bt+/06dOurq4ODg7ffPNNcHCw3EUZSWNjo42NjZOTk7u7u27LrFmzpkyZIm9VRnP16tW+ffsOHDjQ1tZ29+7dXbp0SU1NVebdiD755JO9e/cqbdSoi4vLgAEDbG1ts7Ky+vbtu2HDBkXdiC03N3fkyJEjRoyora3Nzs7evn17165d5S7KqL788stly5YZfOKcWQZhXl5eZWVl09OePXsqZ+SkJElHjx69c4unp2enTp3kqsf4CgsLs7Oz6+vru3fvHhgYKHc5stFdGlXaghKXLl3Kzs5uaGjo2bNnv3795C5HBjdv3szMzNRoNGFhYfrfk6/dKCgo4DjOy8vLsLs1yyAEAAAwFEyfAAAARUMQAgCAoiEIAQBA0RCEAACgaAhCAABQNAQhAAAoGoIQAAAUDUEIAACKhiAEAABFQxACAICiIQgBAEDR/h+xpLhg76bIgAAAAABJRU5ErkJggg==", "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.8.5" }, "kernelspec": { "name": "julia-1.8", "display_name": "Julia 1.8.5", "language": "julia" } }, "nbformat": 4 }