{
"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 = [ElementPsp(el.symbol, psp=load_psp(\"hgh/lda/si-q4.hgh\"))\n",
" for el in load_atoms(silicon)]\n",
"positions = load_positions(silicon)\n",
"lattice = load_lattice(silicon);\n",
"\n",
"model = model_LDA(lattice, atoms, positions)\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\n",
"--- --------------- --------- --------- ---- ----\n",
" 1 -7.846712147581 -0.70 0.80 4.2\n",
" 2 -7.852281470674 -2.25 -1.54 0.80 1.0\n",
" 3 -7.852621342771 -3.47 -2.53 0.80 3.0\n",
" 4 -7.852621950103 -6.22 -3.35 0.80 2.2\n",
" 5 -7.852621959119 -8.04 -4.63 0.80 1.2\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+gvaeTAAAgAElEQVR4nO3deUATZ/4/8GeOAAZQQAQEPNCKB4iIoIDgUUAEhaIIrdvWVu3q9tuD7andrdZttau77XZpu22xv6213bZWQbk8EFBUBBVFvABBBQ/wQhC5yczk90daShUlkJBJMu/XX8k4mfkQYt48M8/Mh1IqlQQAAECqaLELAAAAEBOCEAAAJA1BCAAAkoYgBAAASUMQAgCApCEIAQBA0hCEAAAgaQhCAACQNAQhAABImo6CUBCEs2fPnj9/Xje7AwAAUBOlm1usPfPMM4SQ+vp6Nze39evXd7lOe3v7+fPnx48fr+Y2eZ5nGEZrJUIPCYJA0ziiIA68+eLCl4+IBEGgKIqiKC1uUxf/lwoLC69cufK///0vJSUlJSWlurq6y9Wqq6sjIyPV32xzc7OWCoTewPsvoubmZtwlWET48Iuora1NEATtblNHQejn50cIoWna29v71KlTOtgpAACAOnQRhPX19ebm5qrHFhYWd+/e1cFOAQAA1NGHQdjY2Kh6YGdnd+fOHdXj27dv29vb991OAQAAekSjINy/f/8TTzwxdOjQ8PDwzsuPHTs2atSo4cOHDxs2LCcnZ+bMmXv37m1tba2trT158uSUKVM0qxkAAEBrNApClmXnz5//zDPPNDQ0dCxUKpXPPPPMG2+8UVNT849//OMPf/iDvb39n//85+nTp8+ePfujjz7qOEzaCw0K0sp3sfx2a683CQAAkqaFyye+/vrr77777tChQ6qneXl5ERERN2/eZFmWEOLi4vL555/PmTOn2+1UVlZ6eXl5enp2LNm8ebO1tXXndXZW0ZsuMv+bqjBjSGNjo4WFBSHk24vMgVv0Jj+Fhj8I9EhTU5Mmf9OAJpqamuRyuXZnkIP6Or58QPdaWlpMTEzUv3xFLpd3e60Rq3FV97t48aKrq6sqBQkhY8aMuXjxopqvlcvl77zzTsdTW1tbU1PTzivEjCJNRPlMHp0UxPTrx8vl8m/KlLuuC4lBjBkj09aPAOoQBEEul4tdhUSp3nwEoVh4nseHXywURfUoCNW54lb7QVhfX9/5I2JpaVlXV6fma2UyWUhIyKPXWTKaUBS1YJ+weQqdWE52XFZuD5GZ4dpWnaNpGtd0i0X15iMIxYIPv4joX2lxm9oPwkGDBtXX13c8raurs7Oz0+4uFrvSSkKC95namAkZs1lTpCAAAPSW9v+oGTduXElJSUtLCyGE5/lTp065ublpfS+8QFo5kn9TuaKAv9uu9c0DAIBUaDQirKmpKSoqKikpuXv3blZWlp2dnYeHx/jx4z09Pd9555233347ISHBzs4uMDBQW+Wq/Pe8kHxZOBzatqXa/OPT/I8XhL96Mi+No1kcqwAAXcnKytq4caPYVUhIaGjo0qVL+2LLGgXhhQsXNmzYQAhxcHDYsGFDYGCgh4cHIWTLli2vvfZaYGDguHHjkpOTtXsm47/nhe2VQlIwq2gmr4yjLViy6byw55rwZYnw0RRm7lCcNQEAXTh27BjHcQsXLhS7EEnIy8vLzMzUxyD09fXNzMx8cPmQIUMSExM12fLDJF8Wki8L24NZU4aorpZQnS/ce035uT/92hE+/hz51xRmvA3iEAD63NixY2NiYsSuQhIEQdixY0cfbVz7k2X6VLAjHeZM3zc7Zokr/cQwMtCUnJzHflMmhOzmZjnRH01h7PqJVCUAABgOAzurZiEjXc4RHWhKCCEsTZaNoUsWyBzNybhExZpCvq2r29AAAAB0MLAgVIe1KVnvw+RHssV1xHUb9125ljtXAQCAMTHCIFQZNYDaGsRsmsZ8clbwS+WO3EIXUwAA6ILRBqHK447UiSj2xbF0dBYfm81faUQcAoCYfroodHuYKrNK+clZHMrSHSMPQkIITZFFo+iyWHacNfHawa0s4Btxd24AEMn84XRihXLzw7Mws0q5+gT/7GPG/+WsP6TyXpuzZI0XUzSfrW4iYxO5jaWCgMEhAOicKUO2BTFJD8lCVQqmzWJtzXqz8WHDhnV0QdeQIAhZWVk8r+mEw7a2Nltb27a2tkesM2HCBPV7M/QFqQShirM59d0MJimY2Vwu+KZyuTcQhgCgaw/LQg1TkBDS1tameWc9lfb29pCQkEcHmDpYll29enVHP6KH7UsQxDwUbGDXEWrF5EFUbgSbWCE8e4B3tyaf+jEulrgAHwB0R5WFMdk8IcJzo2iijRS8j1Kp/Pnnn4uKikaOHLlo0aKOlnZZWVm5ubl2dnYLFy5UNXzdv3//oEGD3N3dCSHXrl0rKCiYN29ecnIyIeSbb74xMTGJjIx0cHBQZ6fNzc0JCQk1NTV+fn5z585VLTQz++VHSk5O9vLyysnJKS8vDwoKmjFjhnZ+VI1Ja0TYgSIkxoU+F80G2NM+yVxcPn8PJw4BQIc6jwu1noKEkFdeeWXjxo2PPfZYSkpKZGSkauHatWtffvllOzu7U6dOeXp61tbWEkI2bdq0b98+1QrFxcV///vfe7fHlpaWwMDA0tLSoUOHvv322+vXryeEKBSK5cuXKxQKQsiGDRueeuqpsrKyAQMGREdH5+TkaP5jaoUUR4Qd5CxZMYF+dhT1t0LBdavinQnMy240g8EhAPTWnAyuvScH+QQlicsXOIF4D6IW7ud6tK9tQayVSdf/dPny5R9//PHy5cuWlpbPPfeci4tLbm7uhAkT1q1bV1hYOHbsWEJIRETEF1988e6773a5haioKELIkiVL1G9BvHXrVmtr64SEBEJIQECAr6/va6+9dt86YWFhq1atIoTU1tYmJyfryaBQ0kGo4iinEgKY5WPo14/yG0uFj32Z2c4IQwDojTfGMz2aiHeiRlndJJixxNeOCnbs2SE6+cO/v8+cOTN27FhLS0tCiEwmmzRp0tmzZy0sLMzNzVUpSAgJDAwsLCzs0R4JIbdu3XrjjTcIISzLbtq0qfM/nT592tfXV/XYzc2NoqjKysphw4Z1XsfT01P1wMnJKS8vr6d77yMIwl942VI5c9i0K8LLebyLJYn3Y8ZZIQ4BoGced+zB90ZmlTL5snAogjVnybwsbqyVUnW+UHNNTU2tra0dT1taWuRyuVwub21tVSqVqo5Azc3NqtEewzAds0NVrWQfwcLCIjo6mhDyYI/4pqamjseCILS3tz84mmQYfWykLtFzhA8TMZQuWcDOHUJPT+eW5/I1rd2/BACgFzqfF+zHkpQQ9mHXVPTOqVOniouLCSFXr149cuTI1KlTXVxcbGxstm/fTghpaWnZtm2b6sjk8OHDVUNDpVKZlJSkermpqWm/fv1UJxE7k8vlUVFRUVFRHecdO0tNTW1ubiaEJCUlDRkyxMnJSVs/Tp9CEN5PRpM4d7o0RmbGELckxYZTQo+O+AMAdOvB2TGPvr6wFzw8PBYvXhwZGenj4/PXv/515MiRMpls06ZNr776alhYmLu7+8SJE59++mlCyNKlS3NycqZPn+7t7d1x9QVFUS+//PKkSZO8vb2LiorU3OmoUaP8/f3nzp370ksvJSQkPDhq1E84NNq1gaYk3o95cSz95lF+U5nwgTcd42IYv1EA0HMPmyP64DUVmrC1tU1PTy8qKnJ2dnZ0dFQtDAoKunDhQklJyaBBg4YMGaJa6OzsXFpaWlxc7OzsbGtr23FM9R//+MfatWubmppU5xrVERAQ8OKLL5aVlY0dO7Z///6EEDMzs5s3b6quoNi7d2/HpRQvvPDC888/r+HPqC0IwkcZY0Wlh7JZVcrXjvAbS4WPpzAeaPkLABrYV61cfYJPD2UHmnbxr6YM+flxZl4Wx1LkaY3vsmZiYjJ58uT7Fvbr18/Ly+u+hebm5j4+PqrHMpms8xZMTB4yM/UhrK2tp0yZ0nmJnZ2d6kHnQDU1Ne24tFF0CMLuBTtRqpa/obu5ECf6n1MYe7T8BYBeGdmfPCwFVfqxZEcwW9mrDgFDhgxRzUYZO3bsk08+2esieyc0NNTKyqoXL3RycuqcvrqHIFSLquVvjAu94TTvlqh42Y1eOYEx08fZTwCg14ZZdH9UqR9LxvZq1npBQYHqgYeHh4eHRy+2oIl58+b17oVZWVnaraSncN6rBzq3/B29jfuuHDfuBgAweAjCHlO1/P12OvPJWcEfLX8BAAwcgrCXZg5Gy18AAGOAIOw9tPwFADACmCyjKVXL3xdG038pEMYmcqsm0i+MpmlcZAFg7FxcXN58882MjAyxC5GEurq6jhuZah2CUDtULX8LbitfO8J/XSr8y5cJdEAYAhizhQsXjh49WluNcKFbzs7OfbRlBKE2+QyiDkWwiRXCIrT8BZCAB69MB0OEc4Rapmr5W4yWvwAABgJB2Cf6sWTFBPp0NNvKE9etivizAo/DJwAAeglB2IdULX93hbI7LgseSdzuqwhDAAC9gyDsc6qWv+sn06/k8yG7uXN1iEMAAD2CINSRjpa/M3ai5S8AgB5BEOqOquXv+RiZGUPGJaLlLwCAXkAQ6pqNKYn3Yw7OZQ/dEDySuG0VCEMAADEhCMWhavn7uT/zwUkheBd3uhYnDgEAxIEgFFOwE1UYxcaOoEN3c4ty+JstYhcEACA9CEKRqVr+lsTIHM2JW6JiTSHfyotdEwCAlCAI9YKVCVr+AgCIA0GoRzq3/PVL5fLR8hcAoO8hCPWOquXv/42lF6DlLwBA30MQ6iO0/AUA0BkEof5Stfwtms9WN5GxidzGUgFnDgEAtA5BqO9ULX+3BzPflQtTUrhDNxCGAADahCA0DKqWv29PoBcd4CP2cpcaEIcAANqBIDQYnVv+TkbLXwAALUEQGhhVy98z0bJWnoxCy18AAI0hCA3SYDlJCGB2h7I7Lgvj0fIXAEADCEIDpmr5uwEtfwEANIAgNHj3tfy9jZa/AAA9gSA0Bp1b/rqh5S8AQE8gCI2HquXvoQg296YwHi1/AQDUgyA0NqMHUGmz2P+g5S8AgHoQhMYJLX8BANSEIDRaaPkLAKAOBKGRU7X8PfIEWv4CAHQNQSgJj/X/peXvv9HyFwDg9xCEEjJzMHW8U8vfy2j5CwCAIJSazi1/J+3gVhbwDbhzNwBIG4JQijq3/B2Hlr8AIG0IQulCy18AAIIghI6Wv8+h5S8ASBKCEH5p+Xsumg2wp6ekcHH5fH272DUBAOgKghB+oWr5e3q+rJUnrtvQ8hcApAJBCL+javm7ZzabfFkYn8TtQstfADB2CELowsSB1P457IbJdBxa/gKAsUMQwkNFDKWLF7Bzh9Azd3LLc/lbrSQqk79w73ehqCTk1Xw+uxpJCQCGCkEIj6Jq+VsaIzNjiHui4jFLsiCbL7n7S+wpCflzPs9SJMiRErdOAIBeQxBC9zpa/p6/p6xvI+EZfMldpSoFGYr8y5cRu0AAgN5DEIK6VC1/vwxgGEKmpHBPHpQhBQHACCAIoWdmO1MlMewYK2rvdVpJCCeIXRAAgGYQhNAzSkLePMoH2FO7HldsKhPck7iTdzBTBgAMGIIQeqDzeUH/QcKRSLZNIMG7uZUFvAJDQwAwTDoKwqysrFWrVr3//vu62R30hQdnx4yxonaHMoP7UcduKQPSuI7ZpAAABkRHQVhSUjJw4MDExETd7A76AkVIqDN93+yYMVbUlseZb6czL42jp6dzG07hxmwAYGB0FISvvPLK/PnzdbMv6DvhQ7q4XtDdmhpqQS0aRRdEsZlVQkAad74eYQgABgPnCEFrhllQmeHsYlc6MI3bcArNfgHAMLDa3VxISEhtbW3nJTt37nRwcNDuXkBvUYQsG0OHOFFLD/KpV4RN0xjXAbjpDADoNXWDsKWlpaio6Nq1a35+fs7Ozh3Lb9++vWXLlubm5nnz5rm6umZmZvZNnWBIXCyp7Dns16VCQBr3xnjmLQ+aRhoCgL5S99Do6NGjlyxZsmTJkvz8/I6FNTU1EydOPHHiRF1d3eTJk0+cOPGwl2dlZf3444+1tbUbN24sKyvTtGrQe6qh4ZEn2N3XhOk7uftu1Q0AoD/UHRGWlpbK5fLx48d3Xvj11197enp+++23hBAzM7MNGzZs3br1YVuwsbFZvXr1o/fC8/ylS5c6ng4bNoxhcAcvAzbCktoXzv6/84JfKveuJ/OqO0aGAKB31A1CuVz+4MKsrKwFCxaoHoeHh8fHxz/s5cHBwerspb6+PigoqOPpvn37Bg0a9LCVm5qaKArfq6JpampSc80/OBM/K+pPR9nUSvKfyQrnLj5K0DNNTU1KpRKff7Hgy0dELS0tJiYm6o+R5HI5TXdz7FOjyTLV1dUdE2EcHBzu3r3b1NRkbm7e6w3a2NhUVFSoubJSqbSwsOj1vkBz6r//4y3IgQjy8Rlh2l56rTfzxzEYGmrK3Nwc38ViwZePiBiG6VEQqkOjyycYhhGEX+6sJQgCRVE4kgkPw9JkxQR6/xz261IhbA93rQlnDQFAL2gUhIMHD75+/brq8fXr162trc3MzLRRFRgtN2sqP5KdOZielMxtLMX9SQFAfBoF4ezZs1NSUlSPU1JSwsLCtFESGDnV0DA7nE0oFeZkcFUYGgKAqNQNwg8++CA2Nvbq1auffPJJbGxsSUkJIWTp0qWVlZXz58//v//7v6+//nrlypV9WSoYFXdr6mgkO82B9sLQEABERSmVav09XlhY2PmWMd7e3lZWVoSQe/fupaSkNDc3z50718nJSZNSKisrZ86cqf5kmYaGBktLS032CJpobGzUynyB07XK5w7wzuZkYwA7GBNK1dPY2IjJMiLCl4+IejprVB3qzhr18vLqcnn//v2fffZZ7dUDkuNhQx17gv3XGcE7mYv3oxe44P63AKBT+NIB8closmICvSOEWX1CiM3ma1rFLggApARBCPpi8iCqcB47oj/x2K7YUYmzhgCgIwhC0CNmDFnvw2wPZt8pEGKz+TttYhcEABKAIAS942tHnVQNDZO4lMsYGgJA30IQgj7qx5L1Psy2IObtY0JsNl+LoSEA9BkEIegvf3uq6JezhlzaFQwNAaBPIAhBr6mGhltmMq8fERbl8A0KsQsCAKODIAQDEOBAnZrPWpsSj+1cdjVuyQYA2oQgBMMgZ0m8H5MQwCw9yC/P5RsxNAQALUEQgiGZ5USdjmYJIR7buf3XMTQEAC1AEIKB6S8jCQHMF1OZ53L45bl8Eyd2QQBg4BCEYJBmO1NnVEPDJO4AhoYAoAEEIRiqASYkIYD5zJ95Jodfnss3Y2gIAL2CIATDFj7k16Hhdu7QDQwNAaDHEIRg8KxMSEIA829f5un9/MoCvo0XuyAAMCgIQjASc4dSJ+axF++RScnc8RoMDQFAXQhCMB6DzMi2IOY9Lzoig1tZwLfjpmwAoAYEIRibGBe6aL6srJ5M2sEVYmgIAN1BEIIRsu9Htgczq73ocAwNAaA7CEIwWqqhYeld4pPMnbyDoSEAdA1BCMbMoR9JDmHenUjP3sOtLOAVGBoCwAMQhGD8Ylzoonmy4joSkMaV3MXQEAB+B0EIkjBYTlJnMS+No6encxtOCTzSEAB+hSAECVk0ii6IYjOrhIA07nw9whAACEEQgtQMs6Ayw9nFrnRgGoaGAEAIghAkiCJk2Rj62BNsxjVhWjpXhqEhgLQhCEGihltS2XPY50bRAWnchlOCgDQEkCoEIUiXamh45Al29zVh+k7uwj2EIYAUIQhB6kZYUvvC2Wcfo/1TufizGBkCSA6CEIDQFFk2hs6PZJMqhdDd3NUmpCGAhCAIAX4xsj+1L5wNcqS9k7mNpRgaAkgFghDgNyxNVkyg94WzX5cKYXu4axgaAkgAghDgfm7WVH4kO3MwPSmZ21iK+5MCGDkEIUAXVEPD7HA2oVSYk8FVYWgIYLwQhAAP5W5NHY1kpznQXhgaAhgvBCHAo6iGhplh7FclQsRerroZQ0MAY4MgBOiehw119Ak2wJ72SeYTKzA0BDAqCEIAtchosmICvSOEWX1CiM3ma1rFLggAtARBCNADkwdRhfPYEf2Jx3bFjkoMDQGMAYIQoGfMGLLeh9kezL5TIMRm83faxC4IADSDIAToDV876qRqaJjEpVzG0BDAgCEIAXqpH0vW+zDbgpi3jwmx2XwthoYAhglBCKARf3uq6JezhlzaFQwNAQwPghBAU6qh4ZaZzOtHhEU5fINC7IIAoCcQhADaEeBAnZrPWpsSj+1cdjWuuwcwGAhCAK2RsyTej9kYwCw9yC/P5RsxNAQwBAhCAC0LcaJOR7OEEI/t3P7rGBoC6DsEIYD29ZeRhADmy6nMczn88ly+iRO7IAB4OAQhQF8JdabOqIaGSdwBDA0B9BWCEKAPDTAhCQHMZ/7MMzn88ly+GUNDAP2DIAToc+FDfh0abucO3cDQEEC/IAgBdMHKhCQEMP/2ZZ7ez68s4Nt4sQsCgF8hCAF0Z+5Q6sQ89uI9MimZO16DoSGAXkAQAujUIDOyLYh5z4uOyOAwNATQBwhCABHEuNBF82Vl9cQ7mSvE0BBAVAhCAHHY9yPbg5nVXnR4BreygG/H/boBRIIgBBCTamhYepf4JHMn72BoCCACBCGAyBz6keQQ5t2JdNgebmUBr8DQEEC3EIQAeiHGhT45T1ZcRwLSuJK7GBoC6A6CEEBfDJaT1FnMS+Po6enchlMCjzQE0AkEIYB+WTSKLohiM6uEgDTufD3CEKDPIQgB9M4wCyoznF3sSgemYWgI0OcQhAD6iCJk2Rj62BNsxjVhWjpXVq/klaSlq3t2o/0vgIYQhAD6a7gllT2HfW4UHZDGvXSYj8zk7utfsa9aGbmXw4gRQBMIQgC9phoaHnmCLa1XXmtSBu3kOoaAh24oVxzjf3ycpUStEMDQIQgBDMAIS2r/HPY1d+bcXaXXDkUjRx26oXz9CJ8Wyjr0E7s4AAPHil0AAKhFNTQMdKAiMrjRKbLhlkJWODvITOyyAAwfRoQAhmSsFfVlAGvGkJK7yrUn+bo2sQsCMHwIQgBDcuiG8i8FfN5sxZdT6e0VgutWxYZTAno5AWgCQQhgMDrOC9qbKZe40h/7MiP7U0dvKcdv57ZVCJg7CtA7ugvC69ev37p1S2e7AzAyx24r3zrG75r92+yY2BH0K250C6/8YirzYZHgn8rl3UQaAvSYLibLtLe3e3t7Ozo6NjY22tvb//zzzyyLSToAPeNuTaXPYm1/Pzvm6cfoQAdqqAV1Ior93wUhJpufak9tmEy7WOKSCgB16WJEyLJsTk7Onj17cnNz7969e+DAAR3sFMDIyFli29Uc0aEWFCGEpsiiUXR5LDvJlvJJ5uLy+fp2XVcIYKB0EYQ0TdvY2BBClErl3bt3bW1tdbBTAAmSs2TFBPp0NNvKE9dtivizAofuhgDd0fIhyh9++OHOnTudlyxcuHDQoEGqx2vXrvXy8powYYJ2dwoAnTnKqYQAJs6dfuso/1Wp8P4kOsYF0+IAHkrLQWhhYcFxv7sZIsMwqgfx8fGFhYVbt27V7h4BoEvjrKidoWxWlfL1o3xCifDRFMZzIE4cAnRB3SD85ptvsrOzy8rKVq5cGR0d3bH8+++//+CDD+rr66OiouLj45944okuX/7f//43IyNjx44dMplMC1UDgHqCnajCKPabMiE8gwt2pDdMZgbLxa4JQM+oe8Dk6tWrU6dObWtru337dsfCs2fPvvzyy99++21xcXFpaemHH37Y5Wtra2v/+Mc/qrbg7e2dlpamhcIBQD0sTZaNoctiZCP6k/FJipUFfAM6NwF0QimVPbjwaObMmU8++eSf/vQn1dM33nijrq7um2++IYRkZmY+99xz1dXVvS6lsrLS29s7KiqqY8m6desGDBjwsPUbGhosLS17vTvQUGNjo4WFhdhVSFRjY6O5uTlF9fhQ57Vmsv4MvaeKeme88PxjSgbHSnsFXz4iamlpMTEx6Tjp1i0TExOa7mbIp9E5wpKSktDQUNVjT0/P69ev19fXPyK6uiWTyby9vTueyuXyRxxKlclkONAqIrz/IlK9+b0IQpcBJCGAHK8hbx0jX5yn/uFDzXbuiwKNHD78IuI4TiaTqR+E6vw30SgIa2trO/4s6t+/PyGkpqZGkyA0MzPrGG52i2EY9d8L0Dq8/yJSvfm9CEKVKfbkYASTdkWIOyK4FJNPfBl3a4wNewAffhExv9LiNjWaVG1jY9PQ0KB6fO/ePULIwIEDtVAUAPS9iKF0yQI2xoWetZtbnsvfbBG7IACRaBSEo0aNOnfunOrxuXPnBg0aZGVlpY2qAEAXZDRZNoY+Gy2zNiVuiYo1hXwrGlmA9KgbhLdu3bp06VJLS0tNTc2lS5eampoIIc8//3xiYuLZs2ebm5vXr1+/ePHiviwVAPqEjSlZ78PkR7LFdcR1K/ddORpZgLSoG4Rr164NCQm5ffv2pk2bQkJCDh8+TAiZOHHiunXrZs2a5ejoaGNjs3r16r4sFQD60KgB1NYgZmsQk1AqTEnhDt5AGoJU9OzyiT5VWVk5c+bMiooKNdfHDGZx4fIJEfX68gl1KAlJrBBWHBPcrMm/fZmR/TGP5n748hFRTy+fUAfuQAgAv0MREuNCl8awwY60byq3PJe/3Sp2TQB9CUEIAF0woUmcO10aIzNjiFuiYsMpoQ3zaMBIIQgB4KEGmpJ4P+bgXPZEjXL8dm5bBebRgBFCEAJAN8ZYUVuDmC+nMh8WCf6pXN5NpCEYFQQhAKglyJE6EcW+OJaOyeZjs/mKBsQhGAkEIQCoi6bIolF0eSw7yZbySebi8vn6drFrAtAYghAAekbOkhUT6NPRbCtPXLcp4s8KnCB2TQAaQBACQG84yqmEAGb/HHZvlaCaRyN2RQC9hCAEgN4bZ0XtDGU/82M+OCkE7+KK7uDEIRgeBCEAaCrYiSqMYmNH0OEZ3KIc/nqz2AUB9ASCEAC0gKXJsjF0WYxsRH8yPkmxsoBvUIhdE4B6EIQAoDUWMrLGiymaz9a1kfLct3cAABffSURBVHGJ3MZSgcexUtB7CEIA0DJncyohgNkRwvxwQfBI4nZdRRiCXkMQAkCf8LalDsxl10+m4/L5kN3cmVrEIegpBCEA9KGIoXTxAjbGhQ7dwy3P5W+2iF0QwAMQhADQt2Q0WTaGPhctszYlbomKNYV8KxpZgD5BEAKALlibkvU+TH4kW1xHXLdyG0sFdLIAPYEgBADdGTWA2hrEbA1iNpcLvqncwRsIQxAfghAAdM3XjsqNYN/yoJ8/wEfs5S7cQxyCmBCEACACipAYF7o0hg12pP1SueW5/O1WsWsCqUIQAoBoTGgS506XxsjMGOKWqNhwSmjDPBrQOQQhAIhsoCmJ92MORbAnapSqRhY4VAq6hCAEAL0wegC1NYj5cirzYZHgn8rl3UQago4gCAFAjwQ5Uiei2NfH00/n8LHZfEUD4hD6HIIQAPQLTZEYF/pcNDvJlvJJ5uLy+fp2sWsCo4YgBAB9JGfJign06Wi2lSeu2xTxZwVOELsmMFIIQgDQX45yKiGA2T+H3VslqObRiF0RGCEEIQDou3FW1M5Q9jM/5oOTQvAurugOThyCNiEIAcAwBDtRhVFs7Ag6PINblMNfbxa7IDAWCEIAMBgsTZaNoctiZCP6k/FJipUFfINC7JrA8CEIAcDAWMjIGi+maD5b10bGJXIbSwUex0pBAwhCADBIzuZUQgCzI4T54YLgkcTtuoowhF5CEAKAAfO2pQ7MZddPpuPy+ZDd3JlaxCH0GIIQAAxexFC6eAEb40KH7uGW5/I3W8QuCAwKghAAjIGMJsvG0OeiZdamxC1RsaaQb0UjC1APghAAjIe1KVnvw+RHssV1xHUrt7FUQCcL6BaCEACMzagB1NYgZlsQs7lcmJLCHbyBMIRHQRACgHGaYkflRrBvT6CfP8BH7OUu3EMcQtcQhABgtChCYlzo0hg22JH2S+WW5/K3W8WuCfQPghAAjJwJTeLc6dIYmRlD3BIVG04JbZhHA50gCAFAEgaakng/5lAEe6JGqWpkgUOloIIgBAAJGT2A2hrEfDmV+bBI8Evl8m4iDQFBCADSE+RInYhi3xhPP53Dx2bzFQ2IQ0lDEAKAFNEUiXGhz0Wzk2wpn2QuLp+vbxe7JhAJghAApEvOkhUT6NPRbCtPXLcp4s8KnCB2TaBzCEIAkDpHOZUQwOyfw+6tEtyTuG0VCENpQRACABBCyDgramco+7k/88FJIXgXV3QHJw6lAkEIAPCbYCeqMIqNHUGHZ3Cx2fyVRsSh8UMQAgD8DkuTZWPoshjZOGvitYNbWcA3KMiLh/nSu/eH4sZS4aeLOI5q8BCEAABdsJCRNV5M0Xy2ro2MS+Qc5eSp/fy5ut+y8L/nhZTLwrzh+BY1ePgVAgA8lLM5lRDA7AhhsqqUTQoSmflLFv73vLC9UkgKZs0YsUsEjbFiFwAAoO+8bakDc9mkCuH1o4J/GveHYeyVVqSg8cCIEABALdEudHks62dHfX2BDRpMIQWNBoIQAEBd35cLDEW+9W9/57jw0mH0sDASODQKAKCWjvOCiuaWYXOY0N38jRay5XFGhgGFgcMvEACge/fNjgm0pw9HsFnVwvR0DjcpNXQIQgCA7vVjyfbfz46ZOJA6HME6mlOB6RyuuzdoCEIAgO79YSRt+sDsGHdrKjGIWepK+6fxJ2qQhYYKQQgAoJE4d/ozPzpsD5d2BXeZMUiYLAMAoKl5w2lHOTU/i7/sSV4ehwGGgcEvDABAC6bYUbkRzBfFQlw+L+AoqUFBEAIAaIeLJXU4gj1Tq1yQzTdzYlcDakMQAgBojbUp2RPGWsrIzJ3czRaxqwH1IAgBALTJhCbfTmfmD6f9U7niBzo3gR5CEAIAaBlFyIoJ9N996Bnp3J5ryEJ9hyAEAOgTsSPo5BB28QEuoRSXVeg1BCEAQF/xt6cORbD/OoOppHoNQQgA0Ice60/lRbKFd5RP7eNb0a9CLyEIAQD61kBTkhXGsjQJ2sXdbhW7GngAghAAoM+ZMuSHmUyIE+WfypXV4yCpfkEQAgDoAkXIGi9m5QR6Wjp36AayUI/o4l6jPM9/+eWX5eXllpaWixYtcnV11cFOAQD00NLR9FALakE294kv84eRGIroBV38GgRBkMlksbGxbm5uwcHB9fX1OtgpAIB+CnGissPZvxQIawoxeUYv6CIIZTLZ8uXLp06dunDhQltb21u3bulgpwAAesvdmsqPZNOvKJce5BW4yFBsOhqYKxSK5cuXz5o1KywsbNSoUbrZKQCA3hosJwfmsrdbSXgGV98udjXSps1zhO3t7f7+/vctPH78OCGEYZhly5ZVVFS89957L730kqOjoxb3CwBgiMxZsiOE+XM+H5DGpYcywywosSuSqB4E4aVLly5duuTj4zNgwICOhfX19RkZGTRNh4aGWlpaqmLvQTRNT5o0adKkSbt27crLy1uwYIGmhQMAGD6GIp/5M/FnhalpfHII422LLBSBWodGFQrFwIEDvb29Z8+eXVpa2rG8qqrK3d39559//v777ydMmPCwk38VFRXp6enl5eW7du06cOCAj4+PdmoHADAKce705/50+B4u9TJOGIpArSBkGObYsWO1tbXm5uadl8fHx0+bNi0pKSklJWXixIlffPFFly+Xy+V5eXmrVq3KzMzcvn37sGHDtFA4AIARiRpGp4eyLx4WPj2HLNQ1tQ6N0jQ9cuTIB5enpaV9+OGHqscLFiz46KOP1qxZ8+Bq9vb2Has9Wk1NDUX9dmSgrKzMwcHhYSs3Njaqs03oI01NTUolLgoWR1NTkyAInf+zgC710ZfPWDOy93FqwUGTkhr+7xM5Br/errS0tJiYmDAMo+b6crm825U1mixTVVXl5OSkeuzk5FRVVaXJ1gghtra2DQ0N6q9vaWmp4R6h1yiKsrCwELsKiaIoytzcHEEooj768nGzJMfmkdhsbmGebGsQaynri50YNpZlexSE6tDo8gme5zuqYRiG4zhtlAQAIF2WMpI2ix1mQQWmcdeacNBFFzQKwsGDB3dMkLl58yYuigAA0BxLk68CmMWutH8qX3QHWdjnNArCGTNmZGRkqB7v3bt3xowZWqgIAAAIiXOnP/alZ+3mdl9FFvYtdc8Rrlq16tatW62trevXr7ezs/vwww8HDhz42muv+fv7W1pa8jz/888/Hzt2rE9rBQCQlBgX2tmcis7iVk9k/jQWd+juK+oGobu7e319/aRJk1RPTUxMCCFubm7Hjh3bsmWLTCY7fvx4lzNLAQCg1/zsqENz2TkZfMld5Se+DI0JUn2A0p8Z8JWVlTNnzqyoqFBz/YaGBswaFVFjYyNmjYqlsbERs0ZFpPsvn9o2Mi+Ts+tHfTed6aeL7nn6q6eXT6gDY20AAH1nY0r2hrEmNAnaxd1uFbsao4MgBAAwAKYM+d9MZpYz5ZfKna/XlyN5xgFBCABgGChC1ngxf/Gkp6dzB28gC7UGQQgAYEiWuNL/m8HGZHM/XMBdSbUDQQgAYGCCnajscPbdE8KaQl7sWowBghAAwPC4W1N5EezOK8olB3kFRoaaQRACABikwXKSM5e900bC9nB328WuxpAhCAEADJU5S7YHM27WVEAad7kR02d6CUEIAGDAGIrE+zGvutEBaXzBbWRhbyAIAQAM3rIx9NeBTOReLuUyThj2GIIQAMAYzHam9oSxr+RhKmmPIQgBAIzEBBsqP5JJuaxcnstzGBmqDUEIAGA8nMypg3PZq03KuXu5ewqxqzEQCEIAAKNiKSOpIayLJRWYxl1twvSZ7iEIAQCMDUuTL6cyS1xp/1T+5B1kYTcQhAAAxinOnf7Elw7dze28iix8FAQhAIDRWuBCp85i/3iI+6IYk2ceCkEIAGDMfO2o3Aj203NCXD4vYGTYFQQhAICRG2FJ5UWyp2qVsfv4Fk7savQPghAAwPjZmJKM2Ww/hjy+i7vVInY1egZBCAAgCaYM+W4GE+pM+aVypXdxkPQ3CEIAAKmgCFnjxayaSM/YyeVcRxb+AkEIACAtz7vSP8xkn9zH/e8CppISgiAEAJCgIEdqXzi76oSwppDHwBBBCAAgRW7WVH4ku+uqcvEBvl3aI0MEIQCARDn0Iwfnsi08CdvD3W0XuxrxIAgBAKTLjCFbHmcCHajJKVx5vUSPkiIIAQAkTTWV9K3x9LR07vBNKWYhghAAAMgfx9DfTmfnZXJbLkruhCGCEAAACCEk1JnKCmdXFAhrCnmxa9EpBCEAAPzCw4bKj2RSLyv/eIjnJDMyRBACAMBvHOXUwbns9WblnAzunkLsanQCQQgAAL9jISMps9jHBlABadzVJuOfPoMgBACA+zEU+Y8/s9SV9k/lC2uMPAsRhAAA0LU4d/rfvvTsPVz6FWPOQgQhAAA8VLQLnTqLXZbL/afYaCfPIAgBAOBRfO2o3Aj282IhLp8XjHFkiCAEAIBujLCk8iLY07XKmGy+mRO7Gm1DEAIAQPesTUlGGGvOksd3cTdbxK5GqxCEAACgFhOabJ7BzHam/FO5krvGc5AUQQgAAOpS3aF7tRc9cye3/7qRZCGCEAAAeua5UXRSMPvMfv67cmOYSoogBACAHptqTx2Yy6wrElYW8IY+MEQQAgBAbzzWn8qLZA/fVD61j2815H4VCEIAAOilgaYkM4ylKRK0i6tpFbua3kIQAgBA75kx5MeZTIgT5ZfKldUb5FFSBCEAAGhENZV0xQR6ejqXe8PwshBBCAAAWvDCaHrzDDY6m/vpooFNJUUQAgCAdsxyorLC2HcKhDWFhjR5BkEIAABaM96Gyotk0q4oXzjEKwxkZIggBAAAbXKUUwfmsDdblHMyuPp2satRA4IQAAC0zEJGkkNY1wFUQBp3pVHfp88gCAEAQPsYinzuz7wwmvZP40/U6HUWIggBAKCvxLnTn/nRYXu4tCv6e8KQFbsAAAAwZvOG007m1LxM/rIneXmcPo6+9LEmAAAwJpMHUbkRzBfFQlw+L+jfUVIEIQAA9DkXS+pwBHumVhmdxTdzYlfzewhCAADQBWtTsieM7W9CZu7kbraIXU0nCEIAANARE5psns78YSTtn8oV39WXg6QIQgAA0Kk4d/rvPvSMdG7PNb3IQgQhAADoWuwIOjmEXXyASygV/7IKBCEAAIjA3546FMH+64z4U0kRhAAAII7H+lN5kezJO8qn9vEt4k0lRRACAIBoBpqSzDCWpUnQLu52qzg1IAgBAEBMpgz5YSYzy5nyT+XO14twkBRBCAAAIqMIWePFrJxAT0/nDt7QdRYiCAEAQC8sHU1/P4ONyeZ+vKjTqaQIQgAA0BchTlR2OPuXAmFNIa+znSIIAQBAj7hbU/mRbPoV5dKDvEInI0MEIQAA6JfBcnJgLnu7lYRncPXtfb47BCEAAOgdc5bsCGHGWlEBadzlxr6dPqO7ILx48aKjo+OWLVt0tkcAADBcDEU+9WNeGE1PSeX/dPj+U4atPHkpj2/TxplEHXWoVyqVb7/99pQpU9ra2nSzRwAAMAJx7vRQS/L0fv5Wi3J78C+Z1S6Q2Gw+fAhlymhhFzoaEX711VchISFDhgzRze4AAMBozBtG58xh91Ypg3dxhJB2gSzI4sOHUH8aq50I0+aIsKamJjc3t/OSgQMHBgYGXrlyZceOHXv27Pnzn/+sxd0BAIBETB5ElS5gJ+7g/HczTubKuUNpbaUg0W4QNjY2nj17tvOSIUOGBAYGrly5MiAgICkpqby8XBCE6dOnDx8+XIv7BQAAo+dsTp2JlrklKZzkSi2mIFEzCO/cubNly5aTJ082NTX99NNPHcvb29vfe++9vXv32trarlmzxs/P7913333w5ZGRkZcvX7506dK9e/du377d0NCgtfIBAEAa2gWy7BD/nodQ2sC8fYz/x2RtnB4khKgZhJcvX87Ly7Ozs9u8eXPnIFyzZs3Bgwc3b9589OjR8PDw8vJyW1vbB1/+1FNPqR5UVVVNmjRp/PjxWikdAAAkouO84HPDBZkJG3eUaDELKaVS3eszzpw54+XlpVAoVE8VCsXgwYNTUlKmTp1KCAkNDZ09e/Zrr732iC3U1NSYmppaWlp2+a+VlZW+vr6vvvpqx5IXX3zRwsLiYVtraGh42KZABxobGx/x24E+1djYaG5uTlGU2IVIFL58dKxdIE/uJ2HOZNlo0tLSYmJiQjNM3FFizpK/T+rmtSzLdvs/pffnCK9fv37nzh0fHx/VUx8fn9OnTz/6JV2OFzsTBKGurq7jKc/zgvDQG+wIgvCIf4W+hvdfRKo3H0EoFnz4dexaI1kwjCwcQQThtw//Jz4kvphqUSg1v4Ki90F469YtuVxuYmKiemptbX3y5EkNqzE3N//nP/+p5soVFRX9+/fvKAB0qaWlpa6ubuDAgWIXIlHXrl0zMzMzNzcXuxApam9vr6mp6fbPetCi0aZk9K9fNrdu3bK0tLSysiKErJione33fuLNgAEDWltbef6Xy/obGxtVlelMTExMaWmpLvcIHXJzc5cvXy52FdL14osvHjx4UOwqJKq8vHz+/PliVyFdq1atSkxM1O42ex+Ejo6OLMtevHhR9bS8vBwXRQAAgMHpfRCam5tHRUV9+umnhJCKior09PSnn35ae4UBAADoglrnCGtqalxdXQkhlpaWNjY2gwcPPnfuHCHkn//8Z1RU1JAhQxobG//617+OGzeub4sFAADQth5cPvEwN27csLKyMjMz03A7d+7cmTt37q1bt9Rcv66uztLSkmV1dN9w6EyhUDQ3Nw8YMEDsQiSqvr6+X79+mCkmCp7n7927Z21tLXYhEtXQ0CCTydRPnPT09LFjxz56HS0EoRbhvjMAAKBFzs7O3f7JqF9BCAAAoGPoUA8AAJKGIAQAAElDEAIAgKQhCAEAQNIM8tqDmzdvHj9+vKqqKigoaOTIkWKXIy3Nzc179uw5ceIETdMhISHTpk0TuyJpOXjw4MGDB2trax0dHZ999ll7e3uxK5KisrKynJyciIiIwYMHi12LhBw5cqRza4clS5Zo6/I5gxwRTps27cMPP1yxYkVBQYHYtUjOF1988dlnn5mbm5uYmMyfPz8+Pl7siqRl27ZtCoVixIgRJ0+eHD9+fHV1tdgVSY5CoXj22Wfj4uLKysrErkVakpKSvv3220u/0mIDEIO8fEIQBJqmPT09V65c2dH1F3SjtbW141LW77//ft26dbj1uVg8PT1ff/31RYsWiV2ItLz//vuCIHz55Zdbt26dPn262OVIyFtvvWViYrJu3Tqtb9kgR4Q0bZBlG4fON3RobW1Fb16xVFZWVldXu7m5iV2ItJSUlCQlJa1cuVLsQiTq1KlTGzZs+Omnn1paWrS4WSQK9FJNTc3f/va3FStWiF2I5Lz//vvOzs6urq7vvvvupEnd9ecG7eF5fvHixf/5z380v6Mk9MLgwYOdnZ3v3bv32WefeXh41NbWamvLBjlZBkR37969OXPmLFiwICYmRuxaJOfNN99cvnz50aNHly5d6uHhMWPGDLErkoqPP/7Y29s7ICBA7EIk6vXXX1c9UCqV06ZN+/TTT9esWaOVLWNECD3W2NgYHh7u4+PzySefiF2LFMnlcnt7+8jIyOjo6KSkJLHLkZAff/xx37593t7e3t7etbW1y5Yt++abb8QuSoooivL397906ZK2NogRIfRMc3NzZGTk6NGjP/30U4qixC5HWnie5zjO1NSUEMJx3MmTJxcsWCB2URLy448/dpyamjVr1ptvvhkWFiZuSZLS0tLSr18/QkhbW9vevXsXLlyorS0b5KzRV155JT8/v7i42MHBwcbG5quvvvL29ha7KKlYu3bt6tWrPT09VVOWTExM8vLyxC5KKurr6x977DF/f/8BAwYcPnzYwcEhIyMD85VEYW9vj1mjOjZy5MgxY8ZYW1vn5uYOHTp0z549crlcK1s2yCAsLy+/d+9ex1NXV1dLS0sR65GU6urq69evdzylKMrLy0vEeqSmurr6+PHjLS0tI0eOxN9/Ijp16tSIESPwzaNLVVVVJ06caG5uVn34tXhEyiCDEAAAQFswWQYAACQNQQgAAJKGIAQAAElDEAIAgKQhCAEAQNIQhAAAIGkIQgAAkDQEIQAASBqCEAAAJA1BCAAAkoYgBAAASfv/PxKK0xhNb2IAAAAASUVORK5CYII=",
"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.2"
},
"kernelspec": {
"name": "julia-1.7",
"display_name": "Julia 1.7.2",
"language": "julia"
}
},
"nbformat": 4
}