{
"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",
"kgrid = [3, 3, 3] # k-point grid\n",
"Ecut = 5 # kinetic energy cutoff in Hartree\n",
"basis = PlaneWaveBasis(model, Ecut; kgrid=kgrid);"
],
"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.real - info.ρin.real))\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.844829700757 NaN 1.97e-01 4.3 \n",
" 2 -7.850320818523 -5.49e-03 2.93e-02 1.0 \n",
" 3 -7.850616896981 -2.96e-04 2.80e-03 1.5 \n",
" 4 -7.850646836043 -2.99e-05 1.27e-03 2.8 \n",
" 5 -7.850647347985 -5.12e-07 6.02e-04 1.0 \n",
" 6 -7.850647506210 -1.58e-07 5.80e-05 1.0 \n",
" 7 -7.850647511527 -5.32e-09 4.72e-06 1.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+gvaeTAAAgAElEQVR4nO3de0BUZf4/8M9zzhnmAqgoolxEMTFveENLFOWOonjJFKvNLDWtdLPaSrPLWr9tW9vdynIrdetbmbqleQPviAgqmeI9NUVSkUSRi4AMMOfy+2OMUFEBZ+bM5f36Cw4zcz5ONG+e5/k85zBFUQgAAMBVcWoXAAAAoCYEIQAAuDQEIQAAuDQEIQAAuDQEIQAAuDQEIQAAuDQEIQAAuDQEIQAAuDQEIQAAuDQEIQAAuLTbBqEoiv/+97+nT5+empp6h+fn5ORcu3at4eeTZbkR1Tk1vBW1ZFnGpf7MFEXBW2GGt6KWoij4uKhljbfitkE4Z86cy5cvz5w587XXXsvOzr7dw6ZMmbJ///6Gn6+yshK/3GZVVVX45Tarrq6WJEntKuxCTU2NKIpqV2EXRFGsqalRuwq7IElSdXW12lXYBVmWq6qqLP6ywu1+sGzZspycHHd395kzZ37zzTehoaEWPzcAAIDq6h8RlpaWurm5ubu7E1HHjh3Pnj1r06IAAABspf4g5DiudgJTURTGmA1LAgAAsJ0/pkYvXLiQnJxsMBjGjBnTvHlzURQrKio8PDxycnKCgoJULBEAAMB6ro8Ijx071rNnz6NHjyYnJ4eGhpaUlEyaNOn111/ft2/ff/7zn0mTJqlbJQAAgJVcD8L58+dPnTr1008/XbVqVadOnb744ou//e1vXbp0Wb58+UcffdS7d+/bPV9RlOLi4oIblZeX131M+kVl+i7pplbR01eVxC0i+kcBAEBdzLwW6O3tvW7dukGDBhHRZ5999sMPP9x5+2Ct3r17nzt3zs3Nre7BqVOnzpkzp+6Rd48JF43s434mY+U1g8GQW8FNynJb/GBNt+auG4VGo1Gr1XIcrmlARqNRo9EIwm17mF1HdXU1x3EajUbtQtRnMplkWdZqtWoXoj5RFE0mk16vV7sQ9cmyXF1d3ai3wmAw8Dx/58cIRFRdXV1cXOzn52c+5Ovr+9tvvzXwHM2bN1+7dm1ERMSdH/aPMHorW/rLYc0HvVmBZHhqr7wsmu/h5XbnZzk3nud1Oh2CkIgEQUAQmrm5uSEIzRCEtRCEtWRZ1mg0BoPBsi/LEZF5UFh3n7s19ry/E8q31dPkPUJSmvxtJN/DC52oAACgPo6IdDqdl5fXpUuXzIcuXbrk6+trjZNNCub2XOF/LVeqcUEVAACwD9cno2JjYzdu3BgWFkZEGzZsiIuLs/iZTl9VktKktZE175/QhieLr4Rwb/XleQwLAcDxbdu27cyZM1Z6cVmWJUly8Qlzb2/vcePGWenFhaSkpI8//nj27NnR0dHV1dWFhYUHDx784osvLHsacwoujeQ7aJSVMfxffpK/PCVvzVeWRvKdmiEMAcCxLVq0qKysrGPHjmoX4pyKiopycnKsGITjx4/38PDo27fv/v37165d2759+/fff79169YWPMfpq8qENOnbSL67F6uoICL64EHewFNGgRKeLL4Tyk/rgoYRAHBs06ZNs94ntYs7ePDg5MmTrff6wvjx481fderU6eWXX7bGOdw19G0U363FDSO/v/Xjt+Ur/u70+A5pywVlUTjvrbPGyQEAAO7EFkMxPwO7KQXN4vxZtxZs72ghpCX1XG1KOe+6ewoBAEAt6s9Jajia15dfFSO88KM0fZd0DfdiAwAAG1I/CM0GtmEHHhKIqP9aMfsKhoYAAGAj9hKERNRMQ4vC+XdCuRFbxHkHpJsvTgoAAGAFdhSEZuOCuENjNfsKlfBkMacMYQgAANZld0FIRG31lDJUeKozF54sLj6Ji9AAgNM6WaqU1tzlMVdr6EQpRgVWZI9BSESMaFoXLm2EsOikPDZVulKldkEAAFaQW04jt4p3yMKrNZS4Vcwtv+0D4N7ZaRCadWvBfhwl9GxJPVebks9jaAgAzmZ4O/Z6b374ZrG4up6fXq2hkVvFl0O4Ee2aeAWujz/+uLCw8J5KvGcnT5787rvv6v3RmTNnvvrqK9uWUw+7DkL6fXPFD7HCiz/KT6RjcwUAOJthAeytvnzilpuz0JyCfwnhRrdv+gf1kiVLrly5ci/lPfLIIwcOHLiXVzhz5syGDRvq/dGvv/76v//9715e3CLsPQjNwnzYgYcEvUA9fxD3XMJcOQA4lVuz0CIpeJPKysqysrKbDpaVld35vnvHjx8vL2/czKyiKEVFRaJ4feAyYsSIb775pvan5eXl1dX1jX/V4xhBSL9vrnj/AW5sqjhnn2TCRCkAOJG6WWjxFFQU5dVXX+3YsWOvXr0GDx5cUFBARHv37u3atWu/fv18fX0//PBDIlq5cuWECRNqn+Xt7f3666/n5ORMmzatX79+K1asaMi5jh071r179/DwcF9f388++4yI1qxZM2bMmCNHjvTq1WvmzJl9+vTx8fF57733LPJPswgHuyf4w0FceFtuSoY4OEVcGsEHN8edKwDATu24qBgbuZoTH8ANShaJKKkj03BsY14jJsD0AkX51v+RuGLFih07duTm5ur1+unTp7/66qtLliyZMGHC22+/PWnSpLNnz4aGhvbv37+mpqbCfGMEIiIqKSl56aWXkpOTP/nkk4iIiIbUoCjK+PHj33rrrUcffTQnJyc0NDQqKspkMlVWVoqieOTIkdmzZy9cuDAnJ6dnz57Tpk1r+L/OqhwsCImojZ6ShwpLTsoDk8U3evOzejjMoBYAXMo3p+VLxsYt5Ygy/XZNUYj2XKJ9hVKjnuujY1G+fL0/Wrt27ZQpUwwGAxE9//zzkZGRx44dMxqNTzzxBBF16NBh3LhxGzZs6NGjRwPPVVlZeerUKSLSaDTdu3evPX7y5Mni4uJHHnmEiDp16jRs2LDNmzf7+fmZf9qyZcvHHnvM/CN/f//c3NxG/QOtx/GCkH7fXDG4LXs8Xdp8Qf5yiOBrULsmAIAb/d+Q+mPpdswzot9E8lqevXNAShkqtNRappJLly55eXmZv/by8iotLS0qKvLy8mLs+giyVatWRUVFdZ+iKMod1g7PnTv35z//mYi8vb3XrFlTe/zy5cstWrSofdmWLVuWlJTUBmGzZs1qH6nVamtq7raD0lYceDjVtQXbO0p40IeFrjWtP4c1QwBwYHXXBW/XR9pkiqIcPXrU/PWRI0fuu+++4ODg8+fPl5aWmg8eOnQoODjYy8urNg5Pnz5tDkKNRiNJN49Nu3btmpmZmZmZWTcFiUiW5by8vKtXr5q/PXz4cKdOnSzzb7AmhxwR1hI4mteXHxbATUyXVv2qfDqI99CoXRMAQCPd2h0zLIAR8YlbREuNC5ctWxYcHOzj4/Piiy++8MILQUFBI0aMmDhx4gsvvLB79+7s7OylS5cS0YkTJz755JOgoKDFixebB3YhISGfffZZQUFB//79g4OD73oiNze3qVOnzpgxY/v27fn5+ePGjUtOTrbAP8CaHHhEWGuAD8t+SNAL1Gu1uBubKwDAodyuR9Sy48K33377+PHjy5Yte/PNN5999lkiWr58eUxMzJdffnnt2rX9+/e3atWqVatWmzdvPnjw4Lp16+bPnz979my9Xr9gwYKoqKijR49eunSpISfq0KHD1KlTlyxZUl5evnv3br1e37Vr10ceeaRNmzbPPPNM7cOmTJkSEBBggX+YJTj2iLCWeXPF6rPyw6nik525/xfKa5wh4gHAyd15p4QFx4XNmzefP39+3SMajeaFF1646WFhYWFhYWHmr2t3ODz33HONOtfQoUOHDh1a+21ISEhISAgRzZ49u/bgiy++SESnT59u1CtbiVPFxdgO3OGxmp9LlPBk8dRVDA0BwN6t+lV+peed9gsOC2Bze/Orz6INwoqcZERYq42e1scLS07Kg5LFN3rzz/fgsNMQAOzWlPvvPhpJDGRETfwk++CDDwICAjIyMpr29MaKioo6dOhQwx/fs2fPN954w3r1NJCzBSH9vrliSFv2eLq06YL8f9hcAQCuKi4uTu0S7sTHx8fHx0ftKpxrarSuLi3Yj6OEAT4sdK1pHTZXAADAbThtENLvmytWxwov75WfSJcqTGoXBAAA9seZg9BsQO2dK1aLuwrQQQMAADdw/iAkIk8NLQrn//0gN3477lwBAAA3cIkgNHuoA3d4rOZ4CQ1KFn/B5goAACAip+wavQMfPa2L55eclMOxuQIALOrChQs///yz2lU4J2vfp8K1gpB+31wR688mpkubLshfDuH9DEhDALgnEydOnDt37pIlS6z0+oqi1N7SwTV17drVei/uckFo1tGT7Rwh/Puo3H+t9J9B3BgL3QYaAFzT6NGjR48ebaUXF0XRZDLp9XorvT64bgAIHM3uxa2O5V/9CZsrAABcl+sGodmDPuzQQ4KXlnquFjOxuQIAwPW4ehASkUGgBWH8p4P4x3ZIc/ZJNdhcAQDgShCE1w0LYNljhBOlFI7NFQAArgRB+AcfPa2L42d248KTxQXHZIQhAIArQBDe7Ilg7qfRwspf5WGbxN8qkYYAAE4OQViPIE+WPkKI9uP6r5XW4H6YAABODUFYP/PmijVx/Jx98hPpUjk2VwAAOCkE4Z080PqPzRUZ2FwBAOCMEIR3oRdoQRj/+SD+sR3SrCxsrgAAcDYIwgYZGsAOPiScq6B+a8WjxRgaAgA4DwRhQ7XW0do4/uUQLmajOP+wjN0VAADOAUHYOE8Ec3tHCyl58rDNYv41hCEAgMNDEDZakCfbMVyI8eP6rBGXn8GaIQCAY0MQNoV5c8WGocLbB+Sk7VJJtdoFAQBAUyEIm65/a3boIcHXQH3XijsvYpoUAMAhIQjvSe3mij+lY3MFAIBDQhBagHlzxfkK6rdWPILNFQAADgVBaBmtdbQmjn85hIvF5goAAIeCILQk850rNuTJQzeLF7C5AgDAESAILayDJ9sxQoj14/quEZflYM0QAMDeIQgtj2fXN1f8v4PYXAEAYO8QhNbSvzU7+JDga6A+a8R0bK4AALBXCEIrMm+uWDyYn5guzcqSqiW1CwIAgFsgCK0u3p8deEjIu0b91oprzslfnrp54dAk01vZkoj1RAAANQhqF+ASWutodSz/zWl5eqbUWseMJprR/fqfICaZJqRJUb5MwN8kAABqwKev7Zg3V7TS0dsH5b8fkqhOCv65O/5DAACoAyNCm+rgydJHCB8dk+fuk36+wlcqcrQfhxQEAFARgtDWOEYvhXD9W7O4TaZAD1oVixQEAFATPoVVYJLpw2Py6z3ksholZJVYZlK7IAAAF3bbIJQk6dixYz/88ENJSYktC3J6teuCf+kun0niq2TqulI8X4GNhgAA6rhtEK5YseKNN9548cUXf/31V1sW5Nxu6o7R83RyvNBKS71WS9lXkIUAACq4bRA+/vjja9euDQ4OtmU1Tu+Xq8rQgBt6RN042jdGiPKjhE1i8nnsJQQAsDU0y9hUDy/Ww4vddFDL0+pYYe9lZWyqdLYXoYkUAMCWBCL66quvVqxYUffon/70pyeeeEKlklzUgz5s10h+xBYpp0z5cADP3RyXAABgFQIRJSUlJSYm1j1qMBhUqselBXmy3SOFh1PFh1OlZVG8AcN1AADr45YvX24wGLxvZDAYKisrc3NzjUZjfn7+xYsX1a7TVXhpaWuC0EZPA9fj1r4AALbAPf/881VVVbf+4OTJk3PmzAkICFi6dOnixYttX5nLEjj6PJx/qjM3cL10qAhZCABgXSwoKGjhwoXDhw9v2vN79+7922+/abXaugcnT578yiuv1Pv4a9euGQwGxrACRpWVlTqdjuNu2xqzJo97ab+w6EEx3s/Ju0mNRqNGoxEEzAVTdXU1x3EajUbtQtRnMplkWb7ps8U1iaJoMpn0er3ahahPluWqqqpGLd4ZDIY7fMyaCR06dDh//nyTy/L09Fy4cOGAAQPqHmzWrJmHh8ftnuLu7o4gJCKO4+4chBO7UqdWysOp7M0+/LNdnbmVlOd5BKGZRqNBEJohCGshCGvJsiwIgsW7WAS9Xl9ZWdnk53Mc16ZNm8DAQAvWBLXCfFhmojBii3SyFK2kAABWwV25cqVt27ZqlwG3dV8ztmeUcKhISUqTjKLa1QAAOB3uyJEj/fr1U7sMuJOWWtqaIOh4it4oXjaqXQ0AgHPh4uLiOnfurHYZcBdanpZG8kMDWNh68WQpWkkBACyGW758udo1QIMwonl9+Tf7cJEbxPSLyEIAAMvg7tDeCXboyc7csihhQpq4NMfJ91QAANiGMzflO6sYP5Y2XPhrtjzvgISBIQDAPUIQOqTuXuyn0UJqvvJomlQlqV0NAIAjQxA6Km8dpQ4XiCh2o3ilnmvkAQBAgyAIHZiOpxXRfKw/C1svnrqKWVIAgKZAEDo2cyvp7F5cRIqYWYAsBABoNAShM5h6P/d1pDBuu7jiDFpJAQAaB0HoJOL9WWqC8No+ed4BNM8AADQCgtB5hLRke0bxKeeVKRmSCSNDAICGQRA6FT8DSx8hFFbR8C3i1Rq1qwEAcAQIQmfjoaE1cXyX5iw8WTxXgfYZAIC7QBA6IZ7RJwP5qfdzg5Kl/VeQhQAAd4IgdFqzenALB3LDN4vrzmHBEADgtgS1CwArGtOeC/Jko7ZKJ0tpdi/80QMAUA98ODq5Xi3ZnpH8/3Ll6bskESNDAIBbIAidn787y0gU8q4piVvFMpPa1QAA2BkEoUvw1ND6OCHIkw1OFvOuoX0GAOAPCEJXIXD02SB+cmdu4HrpYBGyEADgOgSha5nVg/twADd0k7ghD1kIAECEIHRB44K49fHC05nif46jeQYAAEHokgb4sF0jhYXH5VlZkoyRIQC4NgShi+royfaMFI4UK+O3S5Wi2tUAAKgHQei6vLS0JUFwFyh6o3jJqHY1AAAqQRC6NDeOvo7khwWwgevFE6WYJAUAV4QgdHWMaF5f/u/9uegN4o6LyEIAcDkIQiAimtCRWxUrPJYmLjqJVlIAcC0IQrhuUBuWOVL44ChaSQHAtSAI4Q+dmrE9o4SDRcojaZIRraQA4BoQhHCDVlraliAIHMVsFAur1K4GAMD6EIRwMy1Py6L4+AAWtl785SomSQHAySEIoR7mVtK5vbmIFDGjAFkIAM4MQQi3Nbkz922kMH67uCwHraQA4LQQhHAnsf5s+3DhjWx53gFJ7VoAAKwCQQh30cOL7RkpbDivPJUhmTAyBACngyCEu/M1UHqiUFJNCZvF0hq1qwEAsCgEITSIu0A/xPL9WrNB68Wz5WifAQDngSCEhuIZ/aM/P6sHF7Ze3H0JWQgATgJBCI0zrQv3VYTw0Dbxf2ewYAgAzgBBCI02NIClDhdm70MrKQA4AwQhNEXPlixrFL/+nPJ0JlpJAcCxIQihifwMLCNRuFipJG4Ry0xqVwMA0FQIQmg6Dw2tixeCm7PwZPF8BdpnAMAhIQjhnvCMFg7kp3TmBiVL2VeQhQDgeBCEYAGzenALwriEzWLyeSwYAoCDEdQuAJzE2A6cv4GNTZXO9qI/d8cfWADgMPCBBRbzoA/bNZL/7IQ8K0uSMUsKAA4CQQiWFOTJdo8UjhYrD6dKlaLa1QAANACCECzMS0tbE4Q2ehq4XrxwDQNDALB3CEKwPIGjz8P5pzpzA9dLh4qQhQBg1xCEYC2zenD/HsDFbxI35SELAcB+IQjBisYHcevihSmZ4ucnsK0CAOwUghCsK8yHZSYKHx1DKykA2CkEIVjdfc3YnlHCoSIlKU0yopUUAOwMghBsoaWWtiYIOp6iN4qXjWpXAwBQB4IQbETL09JIfmgAC1svnizFJCkA2AsEIdgOI5rXl3+zDxe5QUy/iCwEALuAIARbe7IztyxKGLdd7LdWrLmxmfSykeI2iVhHBABbQhCCCmL82M4RwrkK6rHKZJSuHyysopFbxdd68XpcCh4AbAhBCOro7sVOjBMUYvd9Z7paQ4VVlLhFfK8/H+3H1C4NAFwL/vYG1Xjr6PBYYVCy2G2d0M5dWTBQiPRFCgKArd1pRHju3Lns7Ozy8nKbVQOuxiDQpqECxyknSmlWljR3n7SrQJHQRgMANnTbEeHzzz+fk5PTrl27tLS0r776atCgQbYsC1xEYRWN3iYuHSSfNwpLz5C7wObsk44UKxG+bGQglxjI/AwYIwKAdd02COfPn6/X64nou++++/TTTxGEYHG164JhXiaNhgSeW31WTh0ulJso/aKcfE55bZ8U5Mli/VliO25gG8YhEwHACm4bhOYUJKLc3Nz27dvbqh5wFZeNNHKrOP8BPtKXGY1ERFPu50SFkrZLq2L58UHc+CCSFD7rkpKSJ7/wo3S+Qon05RID2chAzkurdvUA4EQEItqzZ8/GjRvrHn3wwQdHjhxJRFlZWStWrMjMzFSnOnBeOoE+DuMf9LlhlDe9C9epmcL/foxnFN6WhbflqT+dLVe25isp55VZWab7mjFzIvb1ZhglAsA9EojIz8/vpplP8xAwOzv7mWeeSU5Obt68uTrVgfNqpqGbUtAs5jbbJzp4smld2LQuVCXxuwqU1N/kx9Olkmol3p8b2Z4NDeCaaaxcMQA4KSE9PT0yMrJDhw43/eDw4cNTpkxZvXo15kXBruh4ivVnsf78P/pTbrmSmq98c1p+aqfUuxUbGciNbM+6tcAoEQAaQUhKSsrLy9Nqb151ee+99yoqKiZMmEBE/fr1+/zzz9UoD+BOOl4fJnKVIu25pCSfl4dtkjUcxfqzxEAW789pebVLBAC7xzp37vzee++NHTu2ac/v16+fTqfz9vauezAxMfGxxx6r9/EVFRXu7u4MKztElZWVOp2O43BxHzIajRqNRhAsc3mHE1fZpnwurYCyi1hoKyXBTxnVTm7nbpHXtrrq6mqO4zQazPOSyWSSZfnWv9FdkCiKJpOptoHRlcmyXFVVZTAYGv6UhnzMCt26dcvJyWlyWW5ubnFxcd27d697MDg4WKfT1ft4URR1Oh2CkIhkWUYQmimKYsEg7KOjPm1oLlFxNaVdVFLzlahtip5nI9rRyEBucFtys+O3nDGGIDTjeR5BaCaKIs/zt/tQdSmyLBNRo96KhnzGCs2aNSsuLm5yWRqNJjIyMiIiouE1cRyHIKTf3woEIVntrfDWU1JHSupIskIHi5TUfOXtg7Kd79bHb0Ut85uAt4LwW3Eja7wVQklJSc+ePS37ogB2hWMU6s1CvdnsXlxhFXbrA8ANhMOHDz/zzDNqlwFgI611VLtb/1CRknweu/UBXJ3AGIuPj1e7DABb468PE/l5fbFbH8ClCcnJyZZqUgBwULfu1p+YLhVjtz6AaxBCQkLUrgHAXmC3PoALwlgQoH533q0f58/psFsfwCkgCAHuwiBcHyYuCKPcciX5nPLxz/LEdKl/a5bYjnuoAwv0wDARwIEhCAEaoaMnm9WDzerBFVfT9t/k1Hxl/npJz1/vrxniy+x5tz4A1AtBCNAULbXXt2F8pvDm3frzDki1u/VHtGP+7hgmAjgGBCHAPal3t/7snyQ/AxvZHrv1ARwAghDAYurdrX+uQom62259o0j6W/5frPcgAFgcFjQALM+8W39eX37/GGH/GCHWn6WcVzp+Z+q3Vpx3QMq+oih1HnyxkiI2iBcrb3iFs+XK4BSxtMa2dQO4JPzBCWBd7T3uslvf10AfDeBHbhXXx/OteCKi8xXKw9ulReF8Cze1qwdwAQhCABu58279Gd24UVulVZEkE41PlxaF8/28sbQIYAsIQgAV1O7WrzDRtnx50wXlk+MKI+q3nvNwoxXRSEEA20EQAqjJQ0MPdeAe6kBEtClPnrRTKqyi+I3igDYszp+L92e9WuHa3wDWhSAEsAvnK5Q3suW10bKocC/8xKbezx0sUqZmyuZbRMX6s4R2rB32JgJYAYIQQH3nK5SHUqVF4XyIp8RxtHAg/3yWtD6e/0d/VmCkzILrdxJu4cZi/VmsPxsWwHnihhgAFoIgBFDZJSON2SZ9MYTv04pVVxMRDWzDFoTxo7dKqcOFtvrrexPl3y9hs/ik/HSm1L81i/XjYv0Z7psIcI8QhAAqa62jH2L5IM8b4mxQG/ZdDN+8zvaJupewMd8QI/U3efouOe+aEtGWi/Vnw9uxAMydAjQeghBAZRyjm1LQrGN9B81qb4hB/anASFsvyCnn/5g7xV2iABoFQQjg2Nrq6Ylg7ong69d1S83/4y5RmDsFaAgEIYCT4OvMnV4TKeuWuVPcEwOgXghCACfkXmfu9NdyZVu+kpqvzNkneWHuFOAWCEIAJxd0/So2mDsFqB+CEMBV1Dt3+kS6fKUac6fg0hCEAK6odu609grgtXOn5lsnhrdlmDsFF4EgBHB1HW+ZO51/RBq3XcHcKbgIBCEAXFd37rTCRD9evj53WlStDGnLmbts/AzIRHA2CEIAqIeHpp6509k/SX4GNrI9i/XjBrdlWsydglNAEALAXWDuFJwbghAAGurWudPk8/K47XKVqAxuy8X6s5GBnK9B7SoBGglBCABNUTt3uiDsj7nT1/aZgjxZrD/mTsGRIAgB4F7dde401BtTp2C/EIQAYDE3zZ3uuCinnFdq504TA1liINdSq3aVADfi1C4AAJyTh4ZGBnKLwvlfJwi7Rwmx/izlvNLpe1O/teKcfVJqvlIj3/yU0hqalSVJyg0HFaJXf5IuVtqscHA5CEIAsLqOnmxaF+77GL7wcc2icN7Ljc0/Ivl8a4rbJM4/LGdfuR59LdyovQd7bIck/p6RCtHMPRIRoQcHrAdTowBgO3XnTouqKe03OTVfeThVrpGVWD9uZHv2ZGeOSP5TuvTVIGJEM/dI7gK9/wC6bsCKEIQAoI5WWhofxI0PIqLrfacrc5Xpu0wdPVlzNxq2jYI9WUs9UhCsDkEIAOqr7Tutkfndl5StF+RPjivZV9ipJKQgWB3WCAHAjrhxFOnLykz0zP0U3ZaCvxezC5W7Pw3gHiAIAcCOKL+vC4ii26wAABdvSURBVL4XSiuj5HFBbFCKmHwOWQhWhCAEAHuh3NId83UEP70Ll7RDXPjzLZstACwEa4QAYC/KaqhrCzaz2w1/oC8I49vo6eOf5dNlyocDeA7XqAFLw4gQAOxFcze6KQXN5vbm944WjhQr47dLlaLt6wInhyAEAAfgpaUtCYKHhqI2iJeMalcDzgVBCACOwY2jryP4x+7j+q8VDxWhfQYsBkEIAI5kVg/u3wO4oZvFTXnIQrAMBCEAOJjxQdzaOGFKpvjZCbSSggUgCAHA8YT5sMxE4eOf5VlZkoyRIdwbBCEAOKT7mrHdI4UjxUpSGlpJ4Z4gCAHAUbXU0pYEwcBT9Ea0kkLTIQgBwIG5cfR1JP9Qe27gevF4KSZJoSkQhADg2BjR7F7ce/25qA1oJYWmQBACgDNI6sitiRWmZIqfo5UUGglBCABOYmAblpkofHQMraTQOAhCAHAe9zVje0YJh4uVpDTJiFZSaBgEIQA4lZZa2jJM0PMUvVG8jFZSaAAEIQA4Gy1P30TyQwNY2HrxBFpJ4W4QhADghBjRvL783/tzMRvFHReRhXAnCEIAcFoTOnKrYoTH0tBKCneCIAQAZzawDcsciVZSuBMEIQA4uU7N2J5RwqEiZQJaSaE+CEIAcH4ttbQ1QdChlRTqc5cgVBRMJQCAM0ArKdzObYPwtdde69u3b79+/SIiIi5evGjLmgAArMHcSvpWXy5qA1pJ4Q+3DcK33377wIED2dnZMTExn376qS1rAgCwnknB3A+xwuM7pK9Po5UUiO4QhG5ubkRUXFx86tSp+++/34YlAQBY16A2bGci//dDaCUFIiKBiA4cOJCXl1f3aJ8+fQIDA9999921a9cS0UcffaROdQAA1mFuJR27TXwkTfo6gtcLahcE6uGI6Pz584dvVFRURESvv/76vn37pk+fPnv2bLXrBACwsFZa2poguPEUg1ZS1yaUlZWNGTNmzJgxt3tEq1atKioqbFkTAIBtaHlaGsm/fUAamCymxPNdWjC1KwIVCM8+++yyZctu/UFCQkJAQIAkSVlZWd9++63tKwMAsAFzK2kHDzlyg/i/aCHSF1nocpibm1teXp6Pj89NPzAajcePHxcE4f7779fpdLd7/oABA0JCQjp06FD34AMPPBAZGVnv4ysqKtzd3RnDrxpVVlbqdDqOwzUNyGg0ajQaQcAqDVVXV3Mcp9Fo1C5EfSaTSZZlrVZrszPuuUyPptO7ofT4fTY7Z4OIomgymfR6vdqFqE+W5aqqKoPB0PCnCIJw18QRWrdufeLEiVuDUK/Xh4aG3vUciqKUl5eXlpbWPVhZWSnL9fcly7IsyzKCkH5/K9Suwi7Iv1O7EPWZ3wS8FaTGb8UAb0odSmPSuOwryr/6K/bzIYX/QWpZ6a0QWrZseenSpSY/X6fTPfvssxEREQ18vMlk0mq1CEIikiRJq9ViREhEsixjRFgLI0IzjuNsPCIkom5ayhpND20Tn9jF7KeVlOd5juNs/FbYJ1mWFUWx+FvBlZeXt2jRwrIvCgDgoFppaVuCoOEoZqNYWKV2NWAT3MWLF29a4QMAcGVanr6N4uMDWNh68SSuSuoCuJ49e3bu3FntMgAA7Ii5lfT13lzkBnEnrkrq7LjFixerXQMAgD16qjO3LEqYkCYuzUGjijPjevfurXYNAAB2KsaPZSYKfzsoz9knYWDorNCyCABwJ8HN2Z5Rwu5LyiNpUpWkdjVgBQhCAIC7aKWl1ARBYGgldU4IQgCAuzO3ksb5s7D14i9XMUvqVBCEAAANYm4lndubi0hBK6lTQRACADTC5N9bSb9FK6mzQBACADROjB/bPlx4K1uedwCtpM4AQQgA0GjdvdhPo4XtvymPopXU8SEIAQCawltH2xIEjlEsWkkdHIIQAKCJdDwti+Jj/dlAtJI6MgQhAEDTmVtJX+vNRaaIGQXIQoeEIAQAuFeTO3NLI4Xx29FK6pAQhAAAFhDrz9KGC2+ildQBIQgBACyjuxfLGiVszFMmZ0g1GBk6DgQhAIDFtNVTRqJQJVH0BvEKWkkdBIIQAMCSdDwtj+Jj/VnYevEUWkkdAYIQAMDCzK2kc3pxESliJlpJ7R6CEADAKqbcz30TKYzbLi5DK6l9QxACAFhLnD/bPlx4I1uedwDXYbNfCEIAACvq8Xsr6VNoJbVXCEIAAOsyt5IaRUrYLJbWqF0N3AJBCABgdTqeVkTzg9uyB9aJp9FKamcQhAAAtmBuJX21JzcEraR2BkEIAGA7U39vJV1+BguG9gJBCABgU+ZW0tf3o5XUXiAIAQBsrYcX2zNS2HBemZwhmTAyVBuCEABABb4GSk8UiqvRSqo+BCEAgDrcBVoTx4ejlVRtCEIAANWYW0lfCeGGpIi70EqqEgQhAIDKnu7CfR0pPLxdXIFWUjUgCAEA1Bfvz1IThNf2oZVUBQhCAAC7ENKSZY0SNpxXpqCV1LYQhAAA9sLcSlpUTQmbxd2X6lky/OWqUo0Ro6UhCAEA7Ii7QD/E8r1bsaGbxM9P3BB6h4qUR9Okwir01FiYoHYBAABwA57Rvx7k23uyWVlSQSW90YuI6FCRMjlDWhXLB7gztQt0NghCAAB79OduXHsPNj5VvGKkiR3Zs3ulVbF8R0+koOUhCAEA7NSoQJY2QojeIH7/K7d3DB+EFLQOrBECANgvd4Hu86QKkf1lr4S1QStBEAIA2CnzuuC6OHZ0tJT2mxKRIorYVmEFCEIAAHtU2x0T5EF+euVUkuZ4iTJgvWgU1a7M6SAIAQDsjlGk53ZLP9TpjvHR0S9JmrxrStRGscykbnXOBkEIAGB39AJlJAo3dce00tKFRzVhPix6g3jZqFZpTghBCABgj4T6Pp41HH04gB8fxIWtF8+UoXvGMhCEAAAOZnYv7qUQLmKDdKQYWWgBCEIAAMczoxv3/gNc3CZxT32XJIVGQRACADikx+7jlkcJY1PFzReQhfcEQQgA4Khi/Nj6eGFyhvhdLjYYNh2CEADAgT3QmqUOF17ZK390DFnYRAhCAADH1q0Fy0jkPz0hz9mHexU2BYIQAMDhdfBkmYnC1gvKzD2SjBXDRkIQAgA4gzZ6ykgUfrmqTEyXTJglbQwEIQCAk/DQUEq8UC3T8C1iBS7D1mAIQgAA56Hl6btovr0Hi9koFlWrXY2DQBACADgVntGSwXyELxuSLOZfw4Lh3SEIAQCcDSN6/wF+WhducIp0+iqy8C4EtQsAAACrmNWDa6GlqI1ScjzfpxW7+xNcFUaEAABOa1Iw95+B3PDNYmYBxoW3hSAEAHBmo9tzK6KFh1PF1WexqaJ+CEIAACcX6cs2DRNm7pG++AVZWA+sEQIAOL9Qb5Y2XBi2WSqupld6Ygh0A7wdAAAuoUsLljVKWJojz9knYcGwrrsE4T//+c+3337bNqUAAIBV+Rpo5wgh46Ly1E5JxCzp7+4UhBkZGRs3bty5c6fNqgEAAKvy0tK24UKBURm3XarCzSqI6A5BWFlZ+c4778ybN8+GxQAAgNW5C7Q+XtDylLBZLMMlSc3NMoWFhaIo1j3aunXrN998c9asWR4eHioVBgAA1uLG0fIofuYeKXqDuGmY0FqndkGqEoho7ty5BQUFdY8+99xzJ06cmDFjxrFjx4xG48WLF319fVWqEAAALI9n9Nkgfv5heUiKuGUYH+jhupeeEYhoyZIlNx3NyMjw8PCYM2dOSUlJbm7u4sWL//rXv6pRHgAAWNHsXpy7hgYmS5uG8iEtXTQLhQ0bNowYMeKmo0OGDBkyZAgRZWdnv/LKK0hBAABnNbMb5+VGcZvENXFCmI8rZqEwf/78W4OwVseOHe+cgqIoHjhwQJJu6D0KDAzs1KmTxWoEAABr+lMnroWWjd4qLo0Uhga4XBYynuerq6t5nm/a80NDQyVJat68ed2DDz/88OTJk+t9/LVr1wwGA2Mu90bfqrKyUqfTcRyuaUBGo1Gj0QgCrnNE1dXVHMdpNBq1C1GfyWSSZVmr1apdiPpEUTSZTHq93ton2l/EPbpLmN9HHBtop3sMZVmuqqoyGAwNf4rBYLjrx6wgSVJhYWHbtm2bVpaHh8c777wTERHR8Ke4u7sjCImI4zgEoRnP8whCM41GgyA0QxDWslkQRnpQajNl+BbOyHHTu9jjR5Msy4IgNCoIG4IjIp3OtTtnAQCAiIi6e7H0Efy/j8pz9rnQZnvO3d29RYsWapcBAAB2IciTZSYKWy4of94jya5xTVJu9OjRatcAAAB2pI2edowQDhcrT+yUTHa6XGhJ3Jtvvql2DQAAYF9auNG2BMEo0thU0Sje/fEOjevSpYvaNQAAgN3R8vR9DO+jY1EbxaJqtauxJnvsCwIAAHvAM/rvEH5wGxaRIuZfc9oFQwQhAADcFiP654P8xE7c4BTp9FXnzEIEIQAA3MXsXtxbfbiojdKhIifMQgQhAADc3ZOduYUDuYTN4q4CZ8tCBCEAADTImPbc8mhhbKq45qxTbapAEAIAQENF+bKNw4TndktfnnKeLMTVHQEAoBH6ebMdI4Shm6WSavpLiDOMppzh3wAAALbUpQXLGsV/dcpJLkmKIAQAgEbzM7CdicLOi8pTGZLo4LOkCEIAAGiKllpKHS78dk0Zv12qcuSRIYIQAACayF2g5KGCG0/DN4vlJrWraSoEIQAANJ0bR8uj+ODmLHqDWFildjVNgiAEAIB7wjP6PJwfEcgiUsQ8B7wkKYIQAADuFSOa15d/tis3JEX6xdEuSYp9hAAAYBl/7s55aSkiRVwbJwzwYWqX01AYEQIAgMU83on772Bh1FZxa77DjAsRhAAAYEmJgWxVrDAxXVz5q2NsMEQQAgCAhQ1py7YPF/7yo7z4pANkIYIQAAAsr4cX25nI/+uoPO+AvW+2RxACAIBVBHmyjERh7VllVpYk2/GKIYIQAACspa2e0hOFA0XKpJ2SyV5nSRGEAABgRS3caMswobhaeThVMopqV1MfBCEAAFiXQaB1cYK3joZtFq/WqF3NLRCEAABgdQJHXwzh+7dm4cnib5X2tWCIIAQAAFtgRP96kH+8Ezc4Wcops6MsRBACAIDtzO7FvdyTi0iRDhfbSxYiCAEAwKae7cotHMgN2yTuvmQXWYggBAAAW3uoA7csSng4VdyUp34WIggBAEAF0X4sOV6YnCH+3ymVNxjiNkwAAKCO/q1Z2ghh2GappJpeClFtYIYRIQAAqKZrC5aRyC8+Kc/Zp9olSRGEAACgpvYebM8oYedF5bnd6lySFEEIAAAqa6mlbQlCTpkybrtUbfORIYIQAADU56Gh5HiBZzR8i1husumpEYQAAGAXtDz9L5q/rxmL2SheqbLdeRGEAABgL3hGi8L54e3YkBQx75qNFgwRhAAAYEcY0by+/PQu3JAU6ZertshC7CMEAAC7M6sH11JHESniujjhQR9m1XNhRAgAAPZoYiduyWB+5FZx4fF6Lj3zy1WlzEI9NQhCAACwUyMDuf9FCX/5UXpm1w2bKn4uUR5NkwqNlpk4RRACAID9ivZnP44Wvjmt/GnH9Sw8Wao8ni4ti+Lva2aZKVMEIQAA2LU+rdjBsfy6c8rYVPlUGXt0h7Q8iu/awmILh2iWAQAAe3d/c3Z0LB+yWsq+wqeN4IObW7J9BiNCAABwAJUSdfSk+zxpx0UL76nAiBAAAOydeV1wRRTzdxOfzOJFWX6um8UGchgRAgCAXTtZqtSuC7px9H0MvzVf+bS+PRVNgyAEAAD7Jcr0zG7p++g/umPcOPoumt98QTlYZJk5UkyNAgCA/RI42j5c4G9sjtHytCaO5y3UMWPrEWFeXl5FRYWNT2qf8vPzy8rK1K7CLhQUFJSUlKhdhV0oLCy8cuWK2lXYhStXrhQWFqpdhV0oKSkpKChQuwo11QZeWVlZfn7+TQfvna2D8Pnnn9+1a5eNT2qf5s6du2nTJrWrsAvvvvvuypUr1a7CLixYsODLL79Uuwq78OWXXy5YsEDtKuzCypUr3333XbWrsAubNm2aO3euxV8Wa4QAAODSEIQAAODS7rVZRhTFw4cPN/zxJSUlR48eNRgM93heJ1BYWHj8+PGdO3eqXYj6CgoKTp8+jbeCiPLy8kpLS/FWENHZs2fLy8vxVhDR6dOnCwoK8FYQ0fHjxwsLCxv1VvTr18/d3f3Oj2GKck/tpytWrFi4cKEgNDRQKyoqdDpdwx/vxK5du+bm5qbRaNQuRH1Go5HneTc3N7ULUV9VVRVjTKvVql2I+qqrqxVF0el0aheivpqaGkmS9Hq92oWoz2Qy1dTU3DXY6vriiy86dep058fcaxACAAA4NKwRAgCAS0MQAgCAS0MQAgCAS0MQAgCAS7Np92ZJSUlubm5AQECbNm1seV57I8vy/v37f/31V19f3/DwcI5z3T9HysrK9u/ff/ny5VatWkVERKBxlIhyc3OvXbsWEhKidiGquXLlyrlz52q/7dq1q4tvuDpx4sThw4dbtGgxaNAgT09PtctRx9mzZ4uKimq/5Xm+d+/elnpx2wXhqFGjtmzZIgjC3/72txdffNFm57VDAwcOrKys7NGjx5EjR3Q63Y4dO1z2l3vGjBkXL1708/P75ZdfioqKMjMzfX191S5KTXl5ef379/fw8KibBK5m9erVb775Zs+ePc3ffvbZZ3dtf3diL7/88rJly4YMGXL16tX9+/e/8cYbalekjlWrVm3ZssX89ZkzZ/R6/c8//2yxV1ds5dSpU9XV1QkJCR988IHNTmqfjh8/bv7CZDJ17979k08+UbceeyDL8uDBg//xj3+oXYjKEhISXnrppcDAQLULUdOiRYvGjRundhV2Yc2aNQEBAZcvX1a7EPsyZMiQ+fPnW/AFbTcpFxwcjIkvs65du5q/EATBz8+vqqpK3XrshCiKXl5ealehpq+++srHxyc+Pl7tQtRXVlaWlpZ27NgxWbbYzVcd0bJly6ZNm1ZTU5OZmVlaWqp2OXYhNzc3Kytr4sSJFnxN112dsgdZWVn79u0bP3682oWoKSUlZezYsd26devdu/dTTz2ldjmqKSgo+Pvf//7Pf/5T7ULsQl5e3vz584cNGzZo0KC6K0Ou5syZMxkZGYmJie+//35wcPC2bdvUrkh9//3vf4cPH27ZNRQEoWpycnKSkpI+/fTT9u3bq12Lmnr27Pn0008/+eSTq1ev3rt3r9rlqGbGjBnz5s1r3bq12oWo76mnnjp+/PiWLVvOnDnj4eHx17/+Ve2KVFNZWVlZWbl///7k5OR33313xowZalekMlEUv/7668mTJ1v2ZXHNT3WcO3cuNjb2rbfeevTRR9WuRWWBgYGBgYEJCQkVFRUffvhheHi42hWp4OTJk5s3b/b29t65c+eFCxeKi4unT5/+3nvvtWzZUu3SVFB7AV6tVjthwoQvvvhC3XpU5OvrGxYWxvM8EcXExEyfPt1oNLryRUc3b94sSVJCQoJlXxZBqIILFy7ExMS8/PLLTz/9tNq12JHS0tJGXUvXmfj4+Hz44Yfmr/V6/Y8//hgaGoo1dSI6fPiwv7+/2lWoJjIy8uTJk+avc3NzW7Zs6copSERffvnlpEmTLH6vAttddHvlypXZ2dkrV67s0KFD//79H3nkEQvuAnEsffv2LSsrGzdunPnbAQMGjBkzRt2S1BIfHz9w4EBvb+9Dhw59//33O3bsCA0NVbsolW3ZsmXatGmuvH1iypQpbdq08fX1zc7O/uGHH9LT0132t+LSpUt9+/Z9/PHHg4KC3n///ZkzZ7700ktqF6Way5cvt2vX7tChQ7X9hpZiuyBMT08/depU7bfR0dEuuzfo22+/raysrP22W7durjkfSETbtm3LysoqLS0NDAxMSkry8/NTuyL1XbhwISMj47HHHlO7ENXs2rUrLS2ttLS0Xbt2SUlJrjwiJKL8/Pyvv/66rKwsLi4uJiZG7XLUlJub++OPP1rjfw3chgkAAFwaukYBAMClIQgBAMClIQgBAMCl/X+vMeG4OgQHIAAAAABJRU5ErkJggg==",
"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.5.3"
},
"kernelspec": {
"name": "julia-1.5",
"display_name": "Julia 1.5.3",
"language": "julia"
}
},
"nbformat": 4
}