{
"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 the silicon lattice,\n",
"see Creating slabs with ASE for more details."
],
"metadata": {}
},
{
"outputs": [],
"cell_type": "code",
"source": [
"using DFTK\n",
"using PyCall\n",
"\n",
"silicon = pyimport(\"ase.build\").bulk(\"Si\")\n",
"atoms = load_atoms(silicon)\n",
"atoms = [ElementPsp(el.symbol, psp=load_psp(el.symbol, functional=\"lda\")) => position\n",
" for (el, position) in atoms]\n",
"lattice = load_lattice(silicon);\n",
"\n",
"model = model_LDA(lattice, atoms)\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 Eₙ-Eₙ₋₁ ρout-ρin α Diag\n",
"--- --------------- --------- -------- ---- ----\n",
" 1 -7.844843852024 NaN 1.97e-01 0.80 4.0\n",
" 2 -7.850295129198 -5.45e-03 2.92e-02 0.80 1.0\n",
" 3 -7.850646529252 -3.51e-04 3.02e-03 0.80 3.0\n",
" 4 -7.850647500614 -9.71e-07 4.76e-04 0.80 3.0\n",
" 5 -7.850647510869 -1.03e-08 2.10e-05 0.80 1.0\n",
" 6 -7.850647511582 -7.13e-10 4.98e-06 0.80 3.8\n"
]
}
],
"cell_type": "code",
"source": [
"scfres = self_consistent_field(basis, tol=1e-8, callback=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+gvaeTAAAgAElEQVR4nO3deUAUdf8H8O8cy8ICIoqCgCgoKHiDeCQeCAKCKKh4ZNJTlplpPP26fJ6sx8rM7FdP9utSy4qy8kIUFBUQUhRPNMUjD0ARFUVRjuXYOX5/bJEp6i677Ozsvl9/zQ67Mx+njTff73xmhhJFkQAAAFgrWuoCAAAApIQgBAAAq4YgBAAAq4YgBAAAq4YgBAAAq4YgBAAAq4YgBAAAq4YgBAAAq4YgBAAAq4YgBAAAq2aiIMzOzo6IiIiJiSkoKHjQexobG0+cOKH7NnmeN0ZpVgRHTF+CIEhdgszgiOkLR0xfgiAY/c6glAnuNVpZWTlw4MC9e/dWV1fHxMScPHlSoVDc/7aSkpLQ0NDi4mIdN1tdXe3o6GjUSi0cjpi+ampqHBwcpK5CTmpra1UqFUVRUhciG2q12tbWlqYxOaeruro6GxsbhmGMuE1THP2dO3eGh4e7ubn5+vr26NHjwIEDJtgpAACALkwRhFeuXHF3d9cuu7u7l5WVmWCnAAAAujBFENrZ2TU0NGiX6+vr7e3tTbBTAAAAXRgUhNXV1VlZWf/973/XrFlz9/rGxsbly5c/++yzH374YV1dnb+/f1MXTGFhYc+ePQ3ZKQAAgBGxhnx41apVP/30k0KhYFl2xowZTetnz55dVFQ0Z86c5OTkvXv3btq0aeHChe+9997t27d9fHy6d+/e4j1Wa4iCJrb3nSW9UU862LZ4qwAAYL2M0DW6atWq5OTkPXv2aF+WlZV169atpKTEzc2turrazc3tyJEjnTt33rp1q1KpjI6ObrZllBBSUlISEhLy3XffNa0ZPny4Uqm8+z2bLworzwgbw1lb5q8eyFVnhF1XxZ9DjdlEZJHQNaovdI3qC12j+kLXqL5ao2vUoBFhsw4ePOjr6+vm5kYIcXR01F44MWvWrClTpjzys2q1esmSJU0vv//+e2dn57vfENaeXHVnxm/n1wzT8PW1FEUlFzFbLtM/DtPU1Bj9n2Jpamtr8RtKL7W1tVKXIDNqtVoQBHzNdKdWqzmOQxDqTt8gVKlUjzy8xg/Ca9euubi4NL3s0KHD1atXdfysk5PTrl27Hv6eOX2IUikk7me+HyyuLVNtvSqkRrK2jPLhnwJCiCiKGN/oC0dMLxRFYUSoF5qmMSLUC8MwMhgRKpVKjUbT9LKhocHW1sin757yo0VCwncp29sK26NYJeZEAQCgpYz/Z4iHh0dpaWnTy9LSUk9PT6PvhRdIHUf2lYvvHuPrOKNvHgAArIXxR4QjRoyoqanZs2fP8OHDCwsLz58/HxkZadxdfPO7kHpR2BfZ8ONl1f8eF77/XXwvmE70xdwCAJhOVlbWypUrDdwIz/M0TWMyWReRkZGzZs1qjS0bFIR79ux56aWXKioqbt68OXDgwNGjRy9btszOzm7ZsmWTJk0aMWJEXl7e4sWL72l4MdA3vwspJcLGcFajJi/1ZtraUF+dFpb+JvxwXvi/oUzPtvg+AYApHDx4kOO46dOnS12IVdi3b19mZqY5BmFQUNC6deuaXjbdMmbWrFnh4eEnT55cunSpIVcN3i/1opB6UUgJZ5UM0Z6H1J4v3F4qDnOlQtK4Gd3pJcGMvfEHugAA9/L3909ISJC6CqsgCMKmTZtaaeMGJYZKpfLx8Wn2R126dOnSpYshG29WuDs91pO+pzvmaT96QhfSXkkSfKgFBwX/9dzigZgpBQAAncgsLRwUpNke0fZKQghxV1HJo5jvRzHLjgtjMrgzt1v9CVMAACB3MgtCXYR2oo7Gs+M60yFpXFI+X4ueUgAAeDALDEJCiIImSb3p45PYygbiv55LPodnQAMAQPMsMwi1MFMKAObm5wvCI/80zywT/1uIP99Nx5KDUAszpQBgPiZ2pTcUi98/OAszy8S3jvAzu1v+L2fzYRXHGjOlAGAmlAxZH8ZsfEAWalMwLYJ1adGNKbt06XLz5k1DSySEECIIQlZWFs/zBm6noaHBxcWl6dnszerXr9+FCxcM3JEhrCIItTBTCgDm4EFZaGAKEkIaGhoMf7KeVmNj45gxYx4eYLpgWfatt95i2YddqtfY2CgIUo5PrO7Kc+1M6RenBFx9DwBS0WZhQjZPiPCkL02MkYL3EEVx7dq1x44d69atW2JiYtOzXbOysvLy8jp27Dh9+nTtbb9ycnI6dOjQu3dvQsjly5cPHToUHx+fmppKCFm9erWNjc348eO1T9Z7JLVavWLFioqKiqFDh44bN067sum5C6mpqYGBgbm5uefOnQsLCxs1apRx/qkGs6IRYRPMlAKA5O4eFxo9BQkh8+fPX7lyZffu3Tdv3jx+/HjtysWLF8+bN69jx46//fZb//79b926RQj59ttvm55/d+rUqffff79le6yrqxs+fPiZM2e8vLxee+21pUuXEkI0Gs1zzz2nfSTRBx98MG3atLNnzzo5OU2aNCk3N9fwf6ZRWO9oSDtTmnNVnL+Px31KAcAoYnZwjfr8aS2IJClf4AQysAM1PUe/Xr71YWxbm+Z/dPHixZ9++unixYuOjo5PPvmkt7d3Xl5ev3793nvvvYKCAn9/f0JIbGzsF198sXDhwma3EBcXRwh5+umnVSqVjvWsW7fO2dl5xYoVhJCQkJAhQ4a89NJL97xn7Nixb775JiHk1q1bqampZjIotN4g1MJMKQAY0ct9GEGfk3RHKsQrtYItS4Z0pMLd9ZuiUz34l9WJEyf8/f0dHR0JIQqFIigoqLCw0MHBwd7eXpuChJDhw4cXFBTotUdCyPXr119++WVCCMuy33777d0/On78+JAhQ7TLvXr1oiiqpKTknntt9u/fX7vg4eGxb98+fffeSvBb/4+ZUtynFAAMN9pdj4mlzDIx9aKwJ5a1Z0l8FuffVnzSSL98amtr6+vrm17W1dWpVCqVSlVfXy+KovapT2q1WjvaYximqTu0rq7u4Vt2cHCYNGkSIYSm7y21tra2aVkQhMbGxvtHk8Z9sryx4Df+H9BTCgCmdPd5QTuWbB7DPuiaipb57bffTp06RQgpLS3dv3//sGHDvL2927Vrl5KSQgipq6tbv369dmaya9eu2qGhKIobN27UflypVNrZ2WlPIt5NpVLFxcXFxcU1nXe825YtW9RqNSFk48aNnTt39vDwMNY/p1VhRPg3mCkFABO4vzvm/j5SA/Xt2/epp55ydXU9ePDgG2+80a1bN0LIt99+m5iY+PXXX589e3bo0KEzZswghMyaNWvo0KEjR46sqakJCAjQfpyiqHnz5gUFBXXu3Pnrr79umtJ8OF9f38cee8zT0/PgwYO//PLL/aNG84Rf8/fCTCkAtKoH9YgaNwtdXFzS09OPHTvm6enp7u6uXRkWFnb+/PnTp0936NChc+fO2pWenp5nzpw5deqUp6eni4tL05zqsmXLFi9eXFtbqz3XqIuQkJDnn3/+7Nmz/v7+bdq0IYTY2tqWl5drr6DYuXNn06UUzzzzzD/+8Q8D/43GgiBsHnpKAaA17LoivnWET49k2yub+amSIWtHM/FZHEuRGQbfZc3GxmbQoEH3rLSzswsMDLxnpb29fXBwsHZZoVDcvQUbmwd0pj6As7Pz4MGD717TsWNH7cLdgapUKpsubZQcxjoPg/uUAoBxdWtDHpSCWnYs2RTOBrq05C/vzp07a7tR/P39p06d2uIiWyYyMnLo0KEt+KCHh8fd6Wt6GBE+AmZKAcCIujg8OuHsWOLfoimoQ4cOaRf69u3bt2/fFmzBEPHx8S37YFZWlnEr0Rd+p+sEPaUAAJYKQagHzJQCAFgeBKF+cJ9SAAALgyBsCcyUAgBYDDTLtByuvgewZt7e3q+88sqOHTsM2UjTDc/g4SorK5tuZGp0+M1tEPSUAlit6dOn9+jRw8AH4dbX19vY2MjlDizS8vT0bKUtIwiNAFffA1in+69M15darba1tUUQSgtH32jQUwoAIEcIQmNCTykAgOwgCI0PPaUAADKCIGwtmCkFAJAFBGErwkwpAID5QxC2OsyUAgCYMwShiWCmFADAPCEITQczpQAAZghBaGqYKQUAMCsIQmlgphQAwEwgCCWDmVIAAHOAIJQYZkoBAKSFIDQLmCkFAJAKgtBcYKYUAEASCELzgplSAAATQxCaI8yUAgCYDILQTGGmFADANBCEZg0zpQAArQ1BKAOYKQUAaD0IQnnATCkAQCtBEMoJZkoBAIwOQSg/mCkFADAiBKEsYaYUAMBYEIQyhplSAADDIQhlDzOlAACGQBBaAsyUAgC0GILQcjQ7UyoS0uyEKWZRAQC0EISW5p6Z0rVFwkv7+Xtir6BCnJTFS1MfAICZQRBaoLtnSl87IFyoEp/P+ysLj90Un9nDfzQY/+kBAAhBEFqwppnS4mqyo0yYvosXCTl2U3x6N78xnPF2pKQuEADALLBSFwCtK7QTVRDP/veE8HYBX3jLhqX5zRFMFwekIADAHxCEls+GJq/3o/u2Iwm7iJIRa3B9BQDAXTA1ahWO3RTfOCzsj2x4rCM1cBP3SaGArlEAAC0TBeHChQv79+8fHR1tmt3B3ZrOC3axF9Mi2Kf86A9+46MyuKtqqSsDADADJgrCKVOmrFmz5vLly6bZHTQ5elN8eje/acxf3TFfDGOm+tBVGhKUqtlaipEhAFg7E50j7Nu376VLl0yzL7hbOyXZNObe7phPhjK7r4kKmszM5bdcpP47hFHhZDEAWCsj//7LzMxsaGi4e014eLitra1x9wK6e1CD6Ag3ihByJI59YR8/aDO3JpTp1w6tpABgjYwchMeOHaupqbl7zfDhwxGEZsvJhvw4illfLERlcAv6MS/2phGGAGBtdA3CtLS0ffv2Xbhw4X/+53+GDBnStD49Pf2TTz6pr6+fNm3avHnzXn311dapE1pRgjcd7EI9kctvKxW+G8l2UkldEACACenaLLNmzRpBEPbu3VtaWtq08vjx4zNmzHjhhReWLVv20Ucf/fDDDw/6+Mcffzx37txLly5NmTIlNzfXwKLB6Lo6Urkx7FBXKihVk34JHTQAYEUoUdTjt16fPn3eeuuthIQE7cu5c+dSFPX5558TQlavXr1y5cr9+/c3+8Hy8vLa2lrtcocOHRwdHe9/T0lJyZAhQ1588cWmNc8//7yDg8ODiqmurm52O/AguhyxAzfIk3vI6E7kf4MJOmhqamoe8g2E+9XW1qpUKorCFLuu1Gq1ra0tTeOSbl3V1dXZ2NgwDKPj+1mWfeQX0qBfdcePH3/22We1y4MGDZo3b96D3unq6qrLBgVBqKysbHrJ87wgPPDReoIgPOSncD9djlhwe7I/hiQdoIamU9+FCP3amaY0M4XvmL60RwxBqDt8x/TVGt8xg4Lw+vXrbdu21S63a9eurq7uzp07Tk5OLd6gvb39hx9+qOObGxsblUpli/dlhXQ8Yh2U5Kcwsr5YmLCL/LM382pf2mpbaDQaDb5jeuE4TqlUIgh1x/O8UqnEiFB3giDoNSLUhUFH39HRUa3+4/YkNTU1DMNgHsliJHjT+yew6aVC1HbuihpnDQHAYhkUhF27dj1//rx2+fz58507dzZuSoO0ujhQOdFsmDsdnMqnXcLsDQBYJoOCcMaMGcnJyVVVVYIgfPHFFzNmzDBWWWAmWJq83o9OCWde2i8k5vJqPLkCACyOrkE4ZswYiqIKCwunTJlCUdTu3bsJIXFxcaGhoT4+Pl5eXnV1da+99lprlgqSGdyRKohnCSEDU7ljNzFNCgAWRb/LJ5p1+/bthoYGHftCH6KkpCQ0NLS4uFjH9+PyCX0ZfsTWFwsv7uOtp4MGl0/oC5dP6AuXT+hL38sndGGEo9+2bVvDUxBkQdtBs7VUiEQHDQBYCvwZAvrp4kDlxLDh7nRwKr/lIjpoAED2EISgN4Yir/ejN41hXj4gJObyteigAQA5QxBCCw3qQB2JZ+1YEpzKHUUHDQDIFoIQWq6NgqwIYd4OoqO2c4sKeAFpCAAyhCAEQyV404cmsDlXxIgMrqwWYQgAMoMgBCPwcqB2xbBjPOhBm/nN6KABAFlBEIJxNHXQvIIOGgCQFQQhGNOgDlRBPGvHkoGpXEEFpkkBQAYQhGBkjgqyIoR5N4geuwMdNAAgAwhCaBWT/+ygGYMOGgAwbwhCaC1eDlRODDveix6wifvlAjpoAMBMIQihFdEUSepNb41k3ypABw0AmCkEIbS64A7UkTh00ACAmUIQgiloO2gWB9HR6KABADODIATTmeRNH4pjc6+KYzK4y+igAQDzgCAEk+ps/0cHTSA6aADAPCAIwdQoQpJ609ui/uigqdFIXRAAWDcEIUhjoAv1WzzrrCT9Urj865gmBQDJIAhBMnYsWT6U+XAwHZ/JLSrgeaQhAEgBQQgSm9iVPhTH/npVHLMNHTQAIAEEIUivsz21K4ad0IUO3MT9jA4aADAtBCGYhaYOmkXooAEA00IQghkZ6EIdi2edlaRvCrevHNOkAGAKCEIwL9oOmo8G0xOz0EEDAKaAIARzFN+VPjZRceC6ODyNK65GGAJAK0IQgplysyPbotipPvTgzdxP6KABgFaDIATzpe2gyYhi30EHDQC0GgQhmLsgF+ronx00e9FBAwDGhiAEGdB20Hw8hJ6EDhoAMDYEIchGXBf6t4mKgzfEkDSuCB00AGAkCEKQE1c7sjWSneZDD9nMrTmPDhoAMAIEIciMtoNmVwz7wXFhSjZ/u1HqggBA5hCEIEu9nakD49lOKjJgE5d3DdOkANByCEKQK20HzZfDmGk5PDpoAKDFEIQgb1Ge1JE49tANMSSNu1CFMAQAvSEIQfZc7Uh6JDvNhx66hfsRHTQAoCcEIVgCbQdNTgz74XFhSjZf2SB1QQAgHwhCsBy9nKkDE9hOKhKYig4aANAVghAsii1Dlg9lvhrGTM/hFxziNZgoBYBHQRCCBYr0pI7EsYW3xOHp6KABgEdAEIJl6mhH0iLZp/3oYWncD+igAYAHQxCCxaIImd2Tzo5m/xcdNADwYAhCsHBNHTQDNnF70EEDAPdBEILl03bQrBzOPI4OGgC4D4IQrEWExx8dNCFp3Hl00ADAnxCEYEW0HTSzetAhadzKMxgYAgAhCEKwNtoOml0x7BenhAR00AAAghCsU0Bbav8EtlsbEpjK7UYHDYB1QxCClbJlyNJgZmUI83gOn5SPDhoA64UgBKs2xoM6Gs8WVYvD0rhzdzA0BLBGCEKwdh1syZYI9pke9PB0buUZoYEn1Zpm3lZRb/LKAMAkEIQAf3XQfHlaiNzOj93OVf09C1NKhBk5nETVAUDrQhAC/CGgLZU/nh3SkZyrIiPT/8rClBLh05PChnBW0uoAoLUgCAH+ou2gWTOKKa0lvTZobjZS2hRMi2AdFVIXBwCtA3/kAtwr3IM6NYmN3qEJ2Kzo6SzkxiAFASwZRoQAzehoR/7Vn3G2IScqxVcP8rdw3T2A5UIQAjQjpUT4v5PC4RjNDyPZrZeEHus1ywsFDtcaAlgiEwWhIAjnzp0rLi42ze4ADNF0XtCBFaf6UF+FsD6OVHqp0CeF23EZ1xoCWBpTnCNsaGgYMGBAQEBAdXU1RVFbtmyxsbExwX4BWmD3NfH//uyOqWkghJCYzlQjT391Wlg2iJ67lw9wJp8OZbwdKakrBQDjMEUQKhSKgwcPOjg4EEIiIyN3794dHh5ugv0CtMDgDtSW+3pE47vSIzrR7ZUk0pP+8pQQnMrN6E4vHsigiQbAAphiapSmaW0KCoJQXl7eqVMnE+wUoGWUDGk23torCSHEhiZJvenjk9jKBuK/gUs+J2CqFEDujDwiXLVqVUVFxd1rnn76aVdXV+3ywoULhw8f3qtXL+PuFMDE3FVU8ijm0A3xxXz+s1PCp0OZIR0xUwogV0YOwi5dujg7O9+9RqlUahfef//94uLiH3/80bh7BJBKcAdq33j2h3PCxCwu3J3+cDDjaid1TQCgP12D8LPPPsvMzDx9+vSiRYsef/zxpvUrVqx45513ampqxo0bt2rVqoiIiGY//sknnxw+fHjt2rUMwxihagDzQBGS6EtP8qY/PM732qCZ14v+Vz9Gie84gKzoeo6wrq5u4sSJ9vb2VVVVTSuPHTu2YMGCbdu2lZWV3bhx49133232s5WVle+888758+eHDBkycODAtLQ0IxQOYDbsWbIokMkfzx6pEPukcOmXcN4QQE4oUdTjf9rQ0NCpU6fOmTNH+zIpKUmtVq9atYoQkpOTM23atPLy8haXUlJSMnDgwAkTJjStee+999q0afOg99fU1Gh7cEBHOGL6asERy7lGv3qEdrMTPwwS/J2sLhFra2tVKhVF4YyprtRqta2tLU3j3ia6qqurs7Gx0X1yUZfDa9A5wrNnz0ZFRWmXe/fuff369du3b7dt27bFG1QoFMHBwU0v7e3tm04x3q+xsfEhP4X74YjpS6PR6HvEorqQsM7ky9NiVDY93YdaNIBysqaLZjmOUyqVCELd8TyvVCoRhLoTBEGvINTl22hQEFZWVjo6OmqXtUO3mzdvGhKEtra2TcPNR2IYBmcc9YIjpq+WHTGGIS/1JYk9yDsFvH+K8O9+zLxeNGMd0aA9YghC3WmPGIJQd8yfjLhNg45++/btm04Z3rlzhxDi4uJihKIA5K+9kiwfymREshtKhEGbubxrVjdNCiAXBgWhn5/fiRMntMvHjx93c3NzcnIyRlUAFiLQhdozjl0USM/8lZ+SzV+sQRwCmB1dg7CoqOjIkSPV1dWXLl06cuTI7du3CSFPP/30xo0bDxw4UFlZuWTJklmzZrVmqQByFetFn5zEBjiToE3cogK+npe6IAC4i65BuGrVqueee44QsnPnzueee+7o0aOEkD59+nz66adPPPFEz549u3fvvnDhwlasFEDOVCxZFMgcm8gWVRG/dVzyOTzSCcBc6Hf5RKsqKSkJDQ3V/VFN1dXVTa06oAscMX210gUnuVfFpHy+vZJ8MpTp286i+kpw+YS+cPmEvvS9fEIXOPoApjaqE3U0nv2HHx2RwSXm8jfqpS4IwLohCAEkQFMk0Zc+PVnhrCS9NmiWFwoc5koBJIIgBJCMs5IsH8rsHsduvyz0SeG2XzaX8xQAVgVBCCCxnm2pjCh22SD6hb187E6uqBpxCGBSCEIAsxDrRZ9OYMPd6UGpXFI+X62RuiAAq4EgBDAXNjRJ6k0fn8TW8yRgA5d8TsDYEMAEEIQA5sVdRa0IYVLCma9OC4M3c/nXkYYArQtBCGCOgjtQe8ez8wLoSVlcYi5/rU7qggAsF4IQwExRhCT60uemKHzakN4bNIsK+Abcmw2gFSAIAcyaPUsWBTL549lTlaRPCre+GNcbAhgZghBABnydqHVhzBePMW8XCGMyuJOVOHEIYDQIQgDZCPegjsaz4zrTo7ZySfn8nUapCwKwCAhCADlR0CSpN30mQUEI8VuvWV4o8BgcAhgGQQggP+2VZPlQZnsUu7FEGLSZ23MNYQjQcghCALka0J7aPY5dFEjPzOVjd3IXaxCHAC2BIASQt1gv+tRkNsiFCtrELSrg6zipCwKQGwQhgOypWLIokDk2kS2qIj3Wc8nncIkFgB4QhAAWwtOeSh7FJI9iPj4hhG7ljt/CTCmAThCEABZlVCeqIJ59yo+OyOASc/kb9VIXBGD2EIQAloamSKIvfXqywt2e9NmoWV4ocJgrBXgwBCGAZXJWkqXBzK/j2B1lQu+NXEYpZkoBmocgBLBkPZyobZHsh4Ppefv42J1cUTXiEOBeCEIAyxfrRZ9OYMPd6UGpXFI+X62RuiAAc4IgBLAKNjRJ6k0fn8TW8yRgA7fyjCBgcAhACEEQAlgVdxW1IoTZNIb57qwwZAuXfx1hCIAgBLA+A12ovePZeQH0pCwuMZe/Vid1QQCSQhACWCOKkERf+twUhU8b0nuDZlEB38BLXROARBCEANbLniWLApn9E9hTlaRPCre+GNcbgjVCEAJYu+5tqHVhzJfDmHcKhPBt3MlKnDgE64IgBABCCAlzpwri2Sk+dPg2Limfv9ModUEApoIgBIA/KGgyuyd9crKCEOK3XrO8UOAxOAQrgCAEgL9ppyTLhzLbo9iNJUJwKrfnGsIQLByCEACaMaA9tXsc+3YQnfgrH7uTK8G92cByIQgB4IFiveiTk9gQV3rQZm5RAV/HSV0QQCtAEALAw6hY8no/uiCeLaoifuu55HO4NRtYGgQhADyapz2VPIr5cRTz8QkhdCv32y2kIVgOBCEA6GpkJ6ognn3aj47K4BJz+et/3put2VREVIJcIAgBQA80RRJ96dMJCnd70jdF88FvQu5V8end/D0TpuerxOjtHGZRQRYQhACgt7Y2ZGkw8+s4dvc1YU4eTxPyRC7fdNFhSbX4eA6/bBBDU5JWCaAbBCEAtFAPJ2prJPvhYPrXa+Kxm2J8Js+LpKRanLKL/2Y406cdYhDkgZW6AACQt1gvOsKD/rhQWHyUD96msGWFH0chBUFOMCIEAEMpGfKvfnTmWKasjrqiFrs6IgVBThCEAGAEJdXiP/cL2eEaXyeq21quEvfsBvlAEAKAoZrOC/ZuK+aNY/ycSPe1XEW91GUB6AZBCAAGuVAlTsz6qzuGImRPLOvXlvTcoLmBLAQ5QBACgEEcFNQPf++OoQjJj2UjPejIDIwLQQYQhABgEFc70su5me6YNaFMjBcVvo3DuBDMHIIQAFrLu0FMXFdqZDp3VS11KQAPhiAEgFa0KJB50pcO3cpdUeN+a2CmEIQA0Lpe70fP6kGHbuUv1yILwRwhCAGg1b3al57dkx6ezhfjSfdgfnCLNQAwhZf70PYsGb2Nz45mfHDrGTAnCEIAMJE5/jRNkdFb+axopnsbZCGYCwQhAJjO7J40TZGR6XxmNBPQFlkIZgFBCAAm9UwPWsWSiAx+RxTT7AWIACaGIHM1w4wAABfcSURBVAQAU3u8G00TEpHBb4/CA5tAeghCAJDAtG40TZGo7XxGFNMXWQiSQhACgDSm+NA0RcZkcOkRbHAHZCFIxhRBKIri0aNHf//99zZt2oSFhdna2ppgpwBg/iZ703YsFbuT2zyGHdwRWQjSMMUF9Y2NjZ9++umVK1d27tw5aNCghoYGE+wUAGQhpjP17Qh2QiaXfx3X2oM0TDEiVCqV3333nXZ58ODBFy5cCAgIMMF+AUAWxnamvh/Jjt/JrQ9jR3XCuBBMzUTnCDmOW716dVlZWdeuXXv06GGanQKAXER6UhvC2Km7uJ9C2dHuyEIwKWMGoUajSUxMvGflzz//3LRcU1Nz69at+vp6e3t7I+4XACzAyE7U+jA2IZv7cRQb7oEsBNMxQhCWlpbSNO3h4cGy7Ntvv938blh29uzZhJCZM2fu2LFj4sSJhu8XACzMcDdqYzgbl8mtHsHEeuGRAGAiOn3VOI4bN26cm5sbRVEHDhxoWq9WqyMiIoYNGzZ48ODx48c3Njb63YcQcvPmzRs3bhBCKioqjh075uPj00r/GACQu2GuVEYUO3sPv+WiIHUtYC10CkKKoqZPn56bm+vo6Hj3+i+++KK+vr6oqKioqKi8vPybb75p9uO3bt2aMmVKcHDwhAkTXnnllf79+xuhcACwUANdqPRI9rk8PhVZCCah09QowzAzZswghFDU3ybuf/7553/+858syxJCZs2a9cMPP8ydO/f+j/v6+ubk5Oiyo8rKysDAwKaXmzZtat++/YPeXFtbe0898HA4Yvqqra2VugSZUavVgiAY/jXrYUs2jKAm71bUqrkJnS05DtVqNcdxNI15YF3V1dXZ2NgwDKPj+1Uq1SMPr0HnCC9evNitWzftcvfu3S9evGjI1gghjo6Oq1atanrp7u6uUCge9GZRFB0cHAzco1XBEWsBHDG9UBSlUqmM8vfWMAeyfawYlUGJCuaJ7habEzRN29raIgh1xzCMXkGoC4OCsKamxs7OTrusUqmqqqoMrYZlg4KCDNwIAFiGfu2oXTFsRAYviCTRF1EBrcWgIHR1da2srNQu37p1y83NzRglAQD8wb8tlRXNjNnGCyL5hx+yEFqFQV+s/v3779+/X7u8f//+AQMGGKMkAIC/9HCisqKZt44IX5yy5JOFICFdR4QZGRk1NTUajSY7O/vSpUvR0dH29vbz5s2bPn36gAEDOI777LPPtmzZ0qq1AoB18nOido9jwrbxvEjm98K4EIxM1yDct2/f9evXZ86cefHixYsXL44ePdre3n7MmDGff/75J598QtP06tWrQ0JCWrVWALBaXR2pXTFM2DZeEElSb2QhGBMliuZyx/eSkpLQ0NDi4mId319dXX3PdY3wcDhi+qqpqUHXqF5qa2uN1TXarEs1Ytg2/ik/+t/9LSQL1Wo1ukb1ou/lE7rA0QcA2fByoPbEsmvOC4uP4nwhGA2CEADkxM2O7IphfykS3i5AFoJxIAgBQGZc7ciuaHZjibDgEC91LWAJEIQAID8d7Uh2NJtRKr52EFkIhkIQAoAsdbAluTFs7lXxlQPIQjAIghAA5MpZSXaOZfPKxef38ubS/g4yhCAEABlra0N2RLG/3RTn5PECwhBaBEEIAPLmZEMyo9mzd8TnkIXQIghCAJA9e5akR7JF1eLMXB6TpKAvBCEAWAJ7lqRFsNfrxSdyeQ5XGII+EIQAYCFULEmPYGs04oxcXoMsBJ0hCAHAcigZsjGcrefJ4znIQtAVghAALIoNTdaHMbxI4jO5BlxhCDpAEAKApbGhydrRjC1LxWdx9chCeBQEIQBYIAVN1o5m2iupuEyujpO6GjBvCEIAsEwMRb4bybjaUWN3cDUaqasBM4YgBACLxVBk9QjG25GK3sFVIwvhARCEAGDJtFnYpx01djtXhSyE5iAIAcDCUYR89hgT6EKFbeVuNUhdDZgfBCEAWD6KkOVDmWFu1JgM7iayEP4OQQgAVoEi5JMhTJQnNWYbV1EvdTVgThCEAGBF3hvIjPOiwrdxN5CF8CcEIQBYl3eCmLiu1Mh07qpa6lLAPCAIAcDqLApknvSlQ7dyV9R4aBMgCAHAKr3ej57Vgw7dyl+uRRZaOwQhAFipV/vSs3vSw9P54mpkoVVjpS4AAEAyL/eh7VkyehufHc34OFJSlwPSQBACgFWb40/TFBm9lc+KZrq3QRZaIwQhAFi72T1pmiIj0/nMaCagLbLQ6iAIAQDIMz1oFUsiMvgdUUwvZ2ShdUEQAgAQQsjj3WiakIgMfnsU06cdstCKIAgBAP4wrRtNUyRqO78tiumHLLQaCEIAgL9M8aFpikRkcOkRbHAHZKFVwHWEAAB/M9mbXj2Cjd3JHbiO6wutAoIQAOBeMZ2pb0ewEzK5fGShFUAQAgA0Y2xn6vuR7PidXM5VZKGFQxACADQv0pPaEMZO28XtuoIstGQIQgCABxrZidoQxj6ew2WVIQstFoIQAOBhhrtRG8PZ6Tlc2iVB6lqgVSAIAQAeYZgrlRHFzt7Db7mILLRAuI4QAODRBrpQWyPZmB0cL5L4rhhCWBQEIQCATgJdqG1RbMwOThDJJG9koeVAEAIA6GpAeyojio3K4Op48kR3ZKGFQBACAOihXztqVwwbkcELIkn0RRZaAgQhAIB+/NtSWdHMmG08L5Kn/JCFsocgBADQWw8nKiuaCd/G13FkbgCyUN4QhAAALeHnRO0ex4Rt43mRzO+FLJQxBCEAQAt1daR2xTBh23hBJEm9kYVyhf9yAAAt18WByolhPjslLDmGa+3lCkEIAGCQzvbUnlh2zXnh3aPIQllCEAIAGMrNjuyKYdcWCW8XIAvlB0EIAGAErnZkVzS7sURYcIiXuhbQD4IQAMA4OtqR7Gh2e6n42kFkoZwgCAEAjKaDLcmJYXOviq8cQBbKBoIQAMCYnJVk51g2r1x8fi+Ph/nKAoIQAMDI2tqQHVHsbzfFOXm8gDA0ewhCAADjc7IhmdHsuTvic8hCs4cgBABoFfYsSY9ki6rFmbk8h6sqzJjpglAQhG+//fbs2bMm2yMAgLRULEmLYK/XizN/RRaaL9MF4Wefffaf//wnPz/fZHsEAJCciiXpEWyNRpyRy2uQhWbJREFYUlKyZ8+eCRMmmGZ3AADmQ8mQjeFsA08ez0EWmiNTBKEoiv/85z8/+OADiqJMsDsAAHNjQ5N1YQwvkvhMrh5XGJoZYz6G6dChQ3PmzLl7Tf/+/b/55psVK1aMHDnSx8fHiPsCAJAXG5qsHc1Mz+EnZnEbw1k75t43iIRgrCAJPYKwqKjo9u3bgYGBd68sLy/ft2+fi4tLSEhIcHDwkSNH7v9gRkbGpUuX1qxZU1pampaW1q5du9jYWEMLBwCQGwVN1o5mxu/kfX7hziawjjZ//ehSjfjMHn5bJMuil9/kdArCw4cPR0RENDY2NjQ0aDSapvV5eXlxcXHh4eEnT5709fXduHFjs5Ofmzdv1i68+OKLQUFBSEEAsFoMRbZEMEO2cN3WcWcTWG0UltaK8Vn8l8MYpKAkdDrq3bp1O3z48P0Nn2+88cbChQt/+eWX/Pz8w4cP5+bmPnw7jz32mJ+fX8sKBQCwDAxF9o9nvR2J73rudiNVWivGZfJfDmMGdcDMqDQoUdT1ngcnTpwIDAxsGhFWVla2a9fu6tWrbm5uhJC5c+cqFIrly5e3uJSSkpI+ffrU1NQ0rTl37px2482qqalxcHBo8e6sEI6Yvmpra+3t7aWuQk7UarWdnR3a4nQhiCQsy+ZCFXG3F78cxA1oh9vP6KSurs7GxoZh7jvF+gAqlYqmHzHka3mzTFlZmUKhcHV11b709PQsKCho8da0XFxcqqurdXyzKIr4ta4XHLEWwBHTC0VRKpUKQaijDWPEwE2as1X0usu2ns60tyOO26MxDKNXEOqi5RPSGo2GYZimb7xCoWhsbDRSVQAAFq60VpyYxW8axb/Rn951RRyUysXu5PaVY1wogZaPCN3c3Orr66urqx0dHQkh169f79Spk/EKAwCwWE3nBXvbN4Z4UG1sqD1XxQhP6qndvKOCvNiLntGdZjA+NJWWjwjd3Nx8fX2zs7O1L3ft2jV8+HAjVQUAYLEu1YhxmfzKkL+6Y17qTQ9zo3KvioWT2P8E0l+dFnqs55YXCmpO2kqthU4jwpqamsWLF9+4cUMQhAULFrRp0+bf//43RVGvv/76/Pnzy8vLDx48WFVVNXny5NYuFwBA7uwV1OoRTL92fxvxvdyH3nVFZGkS60XHetF518RPTwpLf9M86Ucn9WI6qaQq1iroFIQ0TTs7Ozs7Oy9ZsoQQop0LJYTMmjWrU6dOGRkZPj4+S5cutbW1bcVKAQAsQnslaa9sZt5ztPtfK0PcqBA35nwV/X8nhYANmlgv+l/9af+2mC1tFXpcPtHaSkpKQkNDi4uLdXx/0+lJ0BGOmL5wwYm+amtr0TWqF7VabWtr+/D+/hv15PNT/BenhMEdqaReTLiHVR9efS+f0AVuYwAAYNY62JJFgUzpdEWCNz0/nw9K5ZLPCXi6oREhCAEAZEDJkERf+uQkdlEgvfKM4LeeW14o1KKbxhgQhAAAskFTJNaLzotl14cxRypEr581Sfl8Wa25nOGSKQQhAID8BLlQyaOYQ3EsIaRPCpeYy5+sRBy2EIIQAECufByp5UOZ4qmKIBcqajs/JoNLu4STh3pDEAIAyJuTDUnqTV+Yys7sTi84KAzYxCWfEzQIRJ0hCAEALIENTRJ96cLJ7IeDmPXFQpdfNIsK+Nu4A7QOEIQAAJaDIiTcg0qLYLdGskVVpNtaTVI+X4pumodCEAIAWKAB7ankUcyRONaOJUGbuCnZ/KEbiMPmIQgBACxWV0dqaTBzfqpimCs1MYsPSePSLgnIw3sgCAEALFwbxR/dNLN70m8cFvqncCvPCPW81GWZDQQhAIBV0HbTHJ/Ifv4Yk3ZJ8P5Fs6iAv9UgdVlmAEEIAGBdQtyotAh2+1i2qIp0X6dJyucv1lj1dCmCEADAGvVrRyWPYk5NVjgrSdAmLnYnt/+6lcYhghAAwHq52ZFFgUzxNEW4Oz11l5V20yAIAQCsnaOCJPWmi6eyr/ej3zsm9FjPLS8U6qzm0RYIQgAAIOTPR1vsH8+uHs5kXRG812oWFfA3raCbBkEIAAB/o+2m2TOOrWwg3ddqEnP53+9Y8nQpghAAAJrh60QtH8qcSVD4tCEhaVzsTm5fuWXGIYIQAAAeyFXbTTNVEetFP7WbH5jKJZ8TeMsKRAQhAAA8goOCzO5Jn57M/ieQ/ur0H900akvppkEQAgCATrTdNPvGs9+NYPaWiz5rNQsO8VfUsh8eIggBAEA/IW7UujBmbyxbx5FeG7jEXP70bRnHIYIQAABaolsbavlQ5uwURS9nKnwbH7uTyyqTZRwiCAEAoOU62JLX+9FFU9kEb/rFfD4olUs+J3CC1GXpA0EIAACGUjIk0ZcunMQuCqR/OC/4reeWFwq1MummQRACAIBxaLtpMsey68OYIxWi18+apHy+rNbc50sRhAAAYGRBLlTyKOZwHEsI6ZPCJebyhZXmG4cIQgAAaBXejtTyoUzxVEWQCxW13XwfbYEgBACAVuRkQ5J600VT2dk96X8dEgI3ccnnBI05ddMgCAEAoNXZ0CTRlz4xif1wELO+WOjyi2ZRAX+7UeqyCCEIQgAAMBmKkHAPKi2C3RrJFlWRbms1Sfl8qdTdNAhCAAAwtQHtqeRRzJE41o4lQZu4Kdn8oRuSxSGCEAAApNHVkVoazJyfqhjmSk3KkqybBkEIAABSaqMgSb3p81PZpN70u0eF/incyjNCPW+6AhCEAAAgPRuaJHjTByewnz/GpF0SvH/RLCrgbzUQQsjvd8R3j97bZlrPkxf28Q3GyEsEIQAAmJEQNyotgt0xlr2qJj3Wa5LyeSVDrtWJrx38K/QaBTIlm+/jTCkZI+wRQQgAAGanbztqRQhzYpLCWUmCU7lLNWJRNdFmYaNAJmfx0Z2pOf7GiTDWKFsBAAAwOjc7siiQebkPs/p34b+Fwr5ysfAWTVPiOC/aWClIMCIEAAAz56ggSb3pc1PYJcHM3nKaJqIRU5AgCAEAQBZEQlKKxXcHCJ4O1N3nCw2HIAQAAHPXdF7wWV/h0yFULUeMmIUIQgAAMGv3dMdQhHz2GGPELEQQAgCAWSurFaf6/K1HVJuFrnaUtV9HWFZW1thoHrcul4krV640NDRIXYWcXLt2ra6uTuoq5KS8vFytVktdhZzcuHGjpqZG6irMnbcjNaP7H2lVUVFRXV1NCKEIebkPbe3XESYkJJw5c0bqKuRk5syZR48elboKOXn22Wf37dsndRVyMn/+/OzsbKmrkJNXX31169atUlchJ2+++eaGDRuMu00ZByEAAIDhEIQAAGDVEIQAAGDVKFGU+NHATW7evDlu3Ljr16/r+P7KykpHR0eWxV3idIUjpq/bt2/b29srFAqpC5GNO3fu2NnZ2djYSF2IbFRVVSmVSqVSKXUhslFdXa1QKGxtbXV8f3p6ur+//8PfY0ZBSAi5ceOGth0IAADAcJ6eno/8y8y8ghAAAMDEcI4QAACsGoIQAACsGoIQAACsGoIQAACsmiw76cvLyw8fPlxWVhYWFtatWzepy5GBU6dObd++vayszNvbOzExsU2bNlJXZO5ycnLy8vIqKys9PT1nzpzZoUMHqSuSjR9//LF9+/Zjx46VuhBzl5eXd+rUqaaXs2fPlrAYubh9+/b3339fUlKi/R+zY8eORtmsLEeEI0aMWLJkyeuvv37o0CGpa5GHyMjIoqIiLy+vHTt2DBgw4M6dO1JXZO7WrVvHcZyPj8+hQ4f69u1bXl4udUXy8MMPP7zwwguffvqp1IXIwE8//bRmzZqiP0ldjgxcvXp1wIABeXl5Xbt2LSsrO3v2rLG2LMvLJwRBoGm6f//+CxYsmDZtmtTlyEB9fb32+lNBEHx9fd9///0pU6ZIXZRsBAQEvPnmm9OnT5e6EHN37dq1sLCwiRMnHj58OCMjQ+pyzN3cuXM9PDzeeOMNqQuRjSeffFKhUHz99ddG37IsR4Q0LcuyJXT3XRjq6+sdHBwkLEZeLly4cP369V69ekldiAzMmzfv7bffdnZ2lroQ2Thy5MiyZcvWrl2L56PpYtu2bQkJCcnJyZ9//vn58+eNuGUkinVZsmSJq6trRESE1IXIwFtvveXp6RkQEPDOO+/07dtX6nLM3bp16+rr6ydPnix1IbLh6enp5uZ2+/btjz76KDAwsKqqSuqKzFp1dXVFRcWrr7568uTJoqKi4ODgvXv3Gm3romz169fv559/lroKOUlOTu7UqdOZM2ekLkQeamtrr127lpKS0r59+z179khdjlmrqKjw9fUtLS0VRfGjjz6KioqSuiI54Thu4MCBy5Ytk7oQs6Z9fPGSJUu0LxcuXBgTE2OsjcuyaxRaYO3atQsWLMjOzu7Ro4fUtciDSqVSqVTx8fHp6ekpKSkhISFSV2S+du7cWVFRERcXRwgpLy+vqqoaMWLE7t27pa5LHhiGGTJkCPplHs7e3t7Z2TkgIED7slevXqmpqcbaOKZGrUJKSspLL720Y8eOR96FHQghPM83NjZqlzUazbFjx7y8vKQtycxFRkZmZmauWLFixYoVCQkJvXr1QuPoI9XV1WkX1Gp1VlYWzkM/Ulxc3P79+7XL+fn5RjxisuwanT9/fn5+/qlTp9zc3Nq1a/fVV18NHDhQ6qLMl0ajsbe3d3FxcXd3166ZP3/+k08+KW1V5qyioiIgIOCxxx5zdHTcu3evp6fn9u3bVSqV1HXJw8cff5yZmYmu0Ufy8PAIDAx0cnL69ddfe/bsmZ6ejocxPdyFCxdGjhwZEhLCcdyBAweys7P9/PyMsmVZBuG5c+fuPrHs5+fn6OgoYT1mThTFgoKCu9d4eHi4ublJVY8sXL58uaCgoK6urnv37kFBQVKXIyfaqVFfX1+pCzF3paWlBQUF9fX1+I7prqqqKjs729bWdtiwYUa8MYgsgxAAAMBYcI4QAACsGoIQAACsGoIQAACsGoIQAACsGoIQAACsGoIQAACsGoIQAACsGoIQAACsGoIQAACsGoIQAACsGoIQAACs2v8D6DJeyFIk+mUAAAAASUVORK5CYII=",
"text/html": [
"\n",
"\n"
],
"image/svg+xml": [
"\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.7.0"
},
"kernelspec": {
"name": "julia-1.7",
"display_name": "Julia 1.7.0",
"language": "julia"
}
},
"nbformat": 4
}