{
"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.844713381043 NaN 1.97e-01 0.80 3.8\n",
" 2 -7.850287787564 -5.57e-03 2.93e-02 0.80 1.0\n",
" 3 -7.850646596278 -3.59e-04 3.04e-03 0.80 3.0\n",
" 4 -7.850647499895 -9.04e-07 4.64e-04 0.80 2.2\n",
" 5 -7.850647510771 -1.09e-08 1.88e-05 0.80 1.5\n",
" 6 -7.850647511582 -8.11e-10 5.78e-06 0.80 4.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+gvaeTAAAgAElEQVR4nO3deUATZ94H8GeOhBBAQAEVEERFBW/ACooKAgKeoOK5Uq3a2lrrutvDfXdrbde2alu7dmtb3VqrrbZWRVQE5a4XHhXFA60nonjigRzhmOP9I5Z6m5CQySTfz1/JmGR+TCnfPL955hlKFEUCAABgrWipCwAAAJASghAAAKwaghAAAKwaghAAAKwaghAAAKwaghAAAKwaghAAAKwaghAAAKwaghAAAKwaghAAAKyaiYIwNTU1IiIiKipq//79T3tNbW3tsWPHdP9MnueNUZoVwRHTlyAIUpcgMzhi+sIR05cgCEZfGZQywVqjt27d6tWr1/79+8vLy6Oiok6cOKFUKh9/WVFRUXh4+IULF3T82PLycgcHB6NWauFwxPRVUVFhb28vdRVyUllZqVarKYqSuhDZqKqqUqlUNI3mnK40Go1SqWQYxoifaYqjn56eHhUV1axZs9atW/v5+T1jUAgAAGBipgjCq1evtmjRQvvY3d39ypUrJtgpAACALkwRhHZ2djU1NdrHGo3Gzs7OBDsFAADQhUFBePv27W3bts2fP/+bb755cHt1dfWiRYvGjx8/f/78iooKf3//goIC7T8dPXrU39/fkJ0CAAAYEWvIm9esWbNu3TrtifHp06fXb58yZcq1a9dee+21VatWjRkzJiUlRaPRzJ07986dO/7+/m3atGnwHsvriIImqsfOkt6sJq6qBn8qAABYLyPMGv3f//63evXqXbt2aZ9eunSpXbt2ly9fdnV1raysbN68+f79+9u0aZORkaFUKiMjI1n2yelbVFTUu3fv//znP/VbhgwZolI9lG9bisXlp4QNEYyK+XMO5Le/i9lXxbVhmHb1HJg1qi/MGtUXZo3qC7NG9aXvrFFdjq1BI8InOnjwYIcOHVxdXQkhdnZ2QUFBeXl5nTp1GjZs2HPfW1NT8/PPP9c/7d27t5OT04MviHQh1zyYuHRhTZ86vlrDMMzq88yWy/SPfeqqqoz+o1gajUZj3DnHFk+j0eAvlF40Gg0hBEGou6qqKkEQ8GumO41Gw3Gc7n/K1Gr1cw+v8YPw2rVrzZo1q3/q4uJy9epVHd/bpEmTpKSkZ79mehdiYyMk7mNW9RLXlai3XRWSo1kVY9Pwiq2GKIoY3+gLR0wvFEVhRKgXmqYxItQLwzBGv47Q+EFoa2tbV1dX/7S6ulqtVht3F5Pb0yIhkdk2zVTC9hjWBoMcAABoKON/DfHw8CguLq5/Wlxc7OnpafS98ALRcGTvdfHjAqEGC4cBAEBDGX9E2L9/f41Gk5OTEx4efuTIkQsXLsTExBh3Fyt+F5IvCnuja364pP7sGP/jGeHLPkyMJ7oxAGA6mZmZy5cvN/BDeJ6naRrNZF1ER0dPmTKlMT7ZoCDMzs6eNm1aeXl5eXl527ZtY2Jili5damNjs3jx4oSEhODg4AMHDixcuNDR0dFY5RJCVvwuJBUJGyPZuiryty6Msw319Unhr3m8ryNZ2pvxssfvEwCYwoEDBziOGzdunNSFWIW9e/dmZGSYYxAGBwdnZGTUP61fMmbixIlRUVHHjx/v2LGjcfuiyReF5ItCUiRrwxDteUjt+cK0S2InZ9JjEzezE/1/3RklTjwDQOPz8/NLSEiQugqrIAjCpk2bGunDDQpCtVr9tKvjW7RoUb++qBFFutOxnvQjs2Neak8P9ybNbMjEdvTMPL7rRm5pHybCHUNDAAB4PpkNnewV5IlzRJvZEEJI2yZUajT7SS966i4+MZe/oTFxdQAAID8yC0JdDPWiT4xk2zQhXZPqlhwX+Ea/3yIAAMiYBQYhIUTNknkBzK9D2JRLwgubuf03EIYAAPBklhmEWh0cqfRYdlYnOi6DS8zlb9VIXRAAWL2fzgmrzwjPfk1Gifj58ee8BozIkoOQEEIRkuhLn0xQONsQ/w11y08JGBsCgIRGtKY3XBBXPT0LM0rEuYf4ie0s/I+zWbGKY+2kJEtCmO0x7MrTQv8U7vgdpCEASMOGIesjmI1PyUJtCm4dyLo06L5y3t7et27dMrREQgghgiBkZmbyvKELd9XU1Li4uNTfm/2JunXrdu7cOQN3ZAirCEKtHs2oPUPZqR3oqFRuVh5fXvf8twAAGN3TstDAFCSE1NTUGH5nPa3a2tqoqKhnB5guWJadO3fu0+6+V78vQZCyFWz8JdbMGU2RRF96iBf9fj7vv4H7MIhO9LWirwIAYCa0WZiQxRMivOhLE2Ok4CNEUVy3bt2RI0fatm2bmJhoY3P/Fj2ZmZm7d+92c3MbN26cs7MzISQnJ8fV1bVz586EkMuXLx88eDA+Pj45OZkQ8t133ymVymHDhul4XXhVVdWyZctKS0tDQkKGDBmi3Vh/W9nk5OSAgIDc3NwzZ85ERESEhYUZ50c1mDXGQFMbsiSEWRvOfHpMiEjlTt1FpxQATO3BcaHRU5AQMnPmzOXLl7dr127z5s31t4OdP3/+66+/7ubmVlBQ0L1799u3bxNCVq5cmZ2drX1BYWHhxx9/3LA9ajSavn37njp1ysvL6+23316wYAEhpK6u7pVXXtHekmjhwoVjx449ffq0o6PjyJEjc3NzDf8xjcK6RoQP6tuCyo9jlxYKoVu5qR3peQGMCrdzAgDDDN7B1erT5BNEMitP4AQS5EqNy+H02tf6CNZJ+eR/unjx4tq1ay9evOjg4PDiiy/6+Pjs3r27W7duH374YX5+vp+fHyFk6NChX3311b/+9a8nfkJcXBwh5KWXXtL9Pnq//PKLs7PzsmXLCCGhoaHBwcGzZ89+5DWxsbHvvvsuIeT27dvJyclmMii03iAkhLA0mdWZTmhDzTkgdN7IfRHCDGqFhdkAoOH+3oXRa276oVLxSqWgYkmwGxXprl+LTv30v9/Hjh3z8/NzcHAghCgUisDAwOPHj9vb29vZ2WlTkBDSt2/f/Px8vfZICLlx48bf//53QgjLsitXrnzwn44ePRocHKx93KlTJ4qiioqKvL29H3xN9+7dtQ88PDz27t2r794biVUHoZa7mlodxuRcFWfs4b8+Sb7szXjjFhYA0CAD9FnlOKNETL4o7BrK2rEkPpPzcxJfNNKshcrKyurq6vqnGo1GrVar1erq6mpRFLV3faqqqtKO9hiGqZ8dqtE8Z2lKe3v7kSNHEkJo+tFSKysr6x8LglBbW/v4aNK4d5Y3Fms8R/hE4S2pghFspDvdM5mbl8/r1dwAANDXg+cFbVmyOYp92jUVDVNQUFBYWEgIuXTp0r59+/r06ePj49O0adOkpCRCiEajWb9+vbYz2bp1a+3QUBTFjRs3at9uY2Nja2urPYn4ILVaHRcXFxcXV3/e8UFbtmypqqoihGzcuLFVq1YeHh7G+nEaFYLwTwqazOpM7xvO/nZT7LKRyyjBJBoAaBSPz4559vWFDdC1a9fJkycPGzasZ8+e//znP9u2batQKFauXPnGG2/ExsZ27ty5R48eEyZMIIRMmTIlNze3f//+QUFB9VdfUBT1+uuvBwYGBgUFHTlyRMed+vr69u7de8iQITNmzFi2bNnjo0bzhNboo9o4UCnR7NZi4eXdfE8X6r+9mea2UtcEABbkaXNEH7+mwhAuLi4pKSlHjhzx9PR0d3fXboyIiDh79uzJkyddXV1btWql3ejp6Xnq1KnCwkJPT08XF5f6nuqiRYvmz59fWVmpPdeoi9DQ0FdfffX06dN+fn5NmjQhhKhUquvXr2uvoEhPT6+/lGLq1KmTJk0y8Gc0FgThkw31oiPc6UVH+a5Jdf/XjXm9E83gvCEAGCz7ijj3EJ8SzTazecK/2jBk3QAmPpNjKTLB4FXWlErlCy+88MhGW1vbgICARzba2dn17NlT+1ihUDz4CUrlU2amPoWzs3OvXr0e3OLm5qZ98GCg2tjY1F/aKDl5jFslob2Fxe4hbOolISiZy8MtLADAYG2bkKeloJYtSzZFsgEuDfnq3apVK+1sFD8/vzFjxjS4yIaJjo4OCQlpwBs9PDweTF/Tw4jwOXwdqR2x7NZiYWw2378FtTiYMeIVrwBgbXSZlG7LEj+nhgThwYMHtQ+6du3atWvXBnyCIeLj4xv2xszMTONWoi+MCHUy1Is+OoLV3sJiyXEB97AAALAYCEJdOSrJkhAmPZZdd17otZn7rRRhCABgCRCE+unejNozjJ3ZiR6Wzs3K4+/hFhYAADKHINSb9ma/J0YqCCH+G7jn3mwaAADMGSbLNJCzDVkSwkz0FV/bw688LSztw/g36OQ2AMiUj4/Pm2++uWPHDkM+pH7BM3i2O3fu1C9kanQIQoMEuVB7h7JLC4V+W7kJ7egPgxh7KecAA4DpjBs3rkOHDgbeCLe6ulqpVMplBRZpeXp6NtInIwgNpb2Fxeg29DsHeL8N3OfB9Cgf/E4DWIXHr0zXV1VVlUqlQhBKC0ffOFqqyeow5scw5r1DwtB0rqgcc0oBAOQBQWhM/VtSR0awke70C5u5efl8DS91QQAA8DwIQiPT3sLicDx7/h7pksTtuIyhIQCAWUMQNgoPO2p1GPNZL/rVPfzQdO5SJeIQAMBMIQgb0VAv+sRINtCF6p7ELSwQOFxwCABgfhCEjcuWJfMCmP3D2ZyrQlAyt/c6hoYAAOYFQWgK7ZpQ22PYfwfR43P4xFz+ZrXUBQEAwB8QhKYz1Is+MYpt04R02YhbWAAAmAsEoUnZsWReAJMey66/ILywmTtwE2EIACAxBKEEujaldg1l3+hED0/nXtnNl9VKXRAAgBVDEEpDewuLwlEKFUM6beRWn0GjFABAGghCKWlvYbElillaKIRv407cQRoCAJgaglB6AS5U3jD2pfZ0ZCo3K4+vwM1+AQBMCEFoFmiKJPrSh0co7tQQP9zsFwDAhBCEZqSFLVkdxqwJZz45KkSmcr+XoVMKANDoEIRmp18LKj+eHepF90vh5uXz1biFBQBAY0IQmqOHbmGxkduOW1gAADQaBKH5cldTq8OYZaHM7H380HSuuAJxCABgfAhCczfAnToczwa6UD02cfPy+VpMowEAMCoEoQyoGDIvgDkwnD1wU+y6kcu6gqEhAIDRIAhlo20TKjWa/aQXPXUXn5jL39BIXRAAgEVAEMqM9ma/bZqQLkl1S44LPAaHAACGQRDKj5ol8wKYXUPYlEvCC5u5/TcQhgAADYcglKv2jlRGLDunGx2XwSXm8rdqpC4IAECeEITyluBDn0xQONsQ/w11y0/hFhYAAHpDEMqek5IsCWG2x7ArTwv9U7jjuIUFAIA+EIQWokczau8wdmoHOiqVm5XHl+MWFgAAukEQWg7tzX5PjFIQQvz/uIVFDU9qnrRaaVmtiasDADBTCEJL09SGLAlhfgpnPj0mRKRyq84IIzO5R7IwqUj4Sy4nUYEAAOYFQWiZQltQ+XHsMC/6/w7yAiHDM7j6u1gkFQlfnBDWhrOSFggAYC4QhBaLpcmszvSxkQoXG+rwLdI3havm76fg1oGsg0Lq+gAAzAOGBRaupZqsDmPSS8TEXL7DFlUHJyE9FikIAPAnjAitwkAPakkIraDFw7dErNkNAPAgBKFVSCoSvj4pHBlc+2YXekIONyGHr8JcGQAAQojJgnDevHmBgYGDBg0yze7gQfXnBe1ZcX4Qs6gXk3tVCErmjt7G0BAAwFRBOGzYsJUrV16+fNk0u4N6mSXifx+eHTPDj/4gkLFjSUQqt+Q4VmUDAGtnoskyAQEBxcXFptkXPKiXG7XlsTmiUzrQsa2o8joyPofPvSp+249pZiNRfQAAUjNyEObk5NTUPHQfhLCwMJVKZdy9gO6eNkHUXU0RQvYNYz88wnfZWLeyHxvtSZm0MgAA82DkINy3b19FRcWDW0JCQhCEZktBk3kBTN8W9Iu/8iNbU5/0YpSYPgUAVkbXINy2bdu+fftOnz49e/bs4ODg+u1paWlLlizRaDTjxo2bPn36P/7xj8apExpRhDt1OJ59aScXupVbG860a4KhIQBYEV2//69atUqj0ezatevSpUv1G48dOzZ27NgpU6bMnz9/wYIFa9eufdrbFy9e/NprrxUXF48ePTo3N9fAosHoXFVky0B2agc6dCv341lB6nIAAEyHEkU9pg126dJl7ty5CQkJ2qczZswQBOHrr78mhKxYseLbb7/Ny8t74huvX79eWVmpfezq6urg4PD4a4qKioKDg9944436La+++qq9vf3TiikvL3/i58DT6HLECu+Sv+wkHZqQr3sTJ6Vp6jJfFRUVz/gNhMdVVlaq1WqKQlNBV1VVVSqViqZxTkJXGo1GqVQyDKPj61mWfe4vpEHnCAsKCqZNm6Z93KtXr5kzZz7tlc2bN9flAwVBuH37dv1TjuN4/kn3ECKEEMLz/DP+FR6nyxHr4EB2RpN3j9A9t1Lf9xGCXa368gr8julLe8QQhLrTHjG9BiRWTt//KxmGadwgvHHjhpOTk/Zx06ZNNRpNWVmZo6Njgz/Qzs7u008/1fHFdXV1mIajFx2PmIqQL0NJRJEwdieZ7kfP7cHQ1vpnjeM4/I7phed5lUqFINSdIAgYEepFFEW9RoS6MOjoOzg4VFVVaR9XVlYyDIM+ksWIb00fjGNzr4qRqVxJJb6uAoDFMigIvb29z507p3189uxZT09P46Y0SKuVHZU9iB3uTffczG25iBk0AGCZDArC8ePH//DDD+Xl5aIofv311+PHjzdWWWAmaIrM6kwnR7F/2y8k5mKpbgCwQLoGYVRUFEVRx48fHz16NEVRO3fuJITEx8f36dOnbdu2rVu3Lisre/vttxuzVJDMC65UfjwrENIzmSvAUt0AYFl0nSyTkZHx+EaGYb777rvS0tLa2lp3d3ejFgbmpYmC/BjGrD4jRKZy/+rOvNHZaifQAIClMcJUJRcXF6SglUj0pQ8MZ9edF+Iz+Fs1z389AID5w5xd0I+PA7VzCNu9GemysW77ZbRJAUD2EISgN5Ym8wKYH8PYabv4WXl8LeaTAoCcIQihgQa4U/nx7IVyErqVO3sPQ0MAkCsEITScq4psHshol+pefgoDQwCQJQQhGIQi5OWOdPZgdmmhMDqLv1srdUEAAHpCEIIR+DtR+4ezLdWkxyZuz3W0SQFAThCEYBwqhiwJYT4PpkdmcvPyeR5pCAAygSAEY4rzpn+LY3OvilFYqhsAZAJBCEbmaUflDGaHe9M9NnHrzmMGDQCYOwQhGB9FyKzO9LZo9l+/CYm5fCWW6gYAM4YghMbS05U6FM+KWKobAMwbghAaURMF+SGMeT+QjknjlhwXEIYAYIYQhNDoEnzovGHsLxeEmDTumkbqagAAHoYgBFNo7UD9OpgNaU4FbMJS3QBgXhCEYCJYqhsAzBOCEExqgDt1OJ4tqiB9tnBnyjA0BADpIQjB1FxUZHMUM7MT3TcFS3UDgPQQhCCNRN/7S3UnZPF3cLN7AJAOghAko12q211Nemzidl9DmxQApIEgBClpl+peEkKPysJS3QAgDQQhSG+4N10wQrHvhhiZyl3GUt0AYFoIQjALzW1JWgwb500HbOJ+PocZNABgOghCMBf1S3W/ewhLdQOA6SAIwbz0dKXy41lCSFAyd+QW2qQA0OgQhGB2HBRkdRjzQSAdu51bWCBgrW4AaFQIQjBT2qW6txYLMduxVDcANCIEIZiv1g5U7mC2d3MqYFNd2iUMDAGgUSAIwaxpl+peE86+vBtLdQNAo0AQggyEt6SOjWSvaUjvLdxpLNUNAEaFIAR5cFKSdQOYNzrRfbZyS45jYAgARoMgBDlJ9KV3DmFXnhZGYaluADASBCHIjJ8TtW8464GlugHASBCEID/apbqXhTJjc3gs1Q0ABkIQglxFe1KH4tj9N8TQrdyFcoQhADQQghBkrLktSY1hx7ahe23GUt0A0EAIQpA37VLdqTHs3Hws1Q0ADYEgBEsQ5EIdimNtWRKUzB3GUt0AoA8EIVgIBwVZFsp8EEgPwlLdAKAPBCFYlAQfet/w+0t1X62SuhoAkAMEIVgab/s/l+pOxVLdAPA8CEKwQNqlun8awE7fzc/K42t4qQsCADOGIASLFdaSOjaSva4hvbdiqW4AeCoEIVgyRyX5eQAzC0t1A8DTIQjB8mmX6v7+jDAyk7+NpboB4GEIQrAKfk5U3jDW044EbOJ2YaluAHgAghCshXap7uV9mXE5/JyDfB0apQBACEEQgrUZ6EHlx7FHb4t9U7jzWKobABCEYIXcbMm2aHZcGzp4M/cTluoGsHoIQrBG2qW602LY9/KFxFy+ok7qggBAOghCsF6BLlRBPOtsg6W6AawaghCsmi1LloQw84PomO3cvHweS3UDWCEEIQAZ5UMfHM5mXRGjt3N7r4n5pY/moUjI+gtISQDLhCAEIIQQL3sqZxDbpzkVl8mNzuYP3vwz9URCXt/LH7wpUhLWBwCNBkEIcJ92qe4tA9k6XhySzu25LpI/UtCOJYteYKQuEAAaBYIQ4CHBbtTRkYpgNyoyjVtfTCMFASweK3UBAGbHUUk2R7EfFfDT9iqCXMVdQ/G/CYAlw4gQ4AlEQkoqSbw3f7JM7Liey72KiTIAFgtBCPCo+vOCK0O4YyNYiiLjc/ih6VxxBeIQwAKZLghv3rx5584dk+0OoGEemR3jZU9lD2Ka25KWtlSPTdy8fL4Wi7IBWBZTBGFtbW337t3Hjx8fGxs7duxYjuNMsFOAhskqEZ2VD82O8bKnkqOYe3Vk/3D24E2xy0Zux2UMDQEshylmAbAsm5mZ6eLiQggJDw/fuXPngAEDTLBfgAaI9KAiPR6dI+ptT/08gCGEbItmtxYLr+7hOzmTpb0ZL3tcWwgge6YYEdI0rU1BURQrKiqaNm1qgp0CNJKhXvSJkWygy/1OaQ0vdUEAYBgjjwh//vnn27dvP7hl9OjR2hQkhCxYsKBz587du3c37k4BTMyWJfMCmInt6Fl5fNck7osQJtoTQ0MAuTJyECqVShsbmwe30PT9QedXX32Vl5e3YcMG4+4RQCptm1Ap6JQCyJ+uQfj9999nZ2f//vvvb7/99siRI+u3//TTT//+97/Lysri4uIWL148YsSIJ779u+++27JlS3JyslKpNELVAGZjqBcd6U4vPMr32MTN7ET/oxtjg1VoAGRF13OEFy5cCAoKqqqqunnzZv3GwsLC6dOnf/PNN4cPHz569OjHH3/8xPfeuXNn2rRply9fDg0NDQoK2rp1qxEKBzAb2k7pgeHsoVKxSxK3HXNKAWSFEkU9/qcNDw8fM2bM9OnTtU/ffPPN0tLS77//nhCSnp4+efLkkpKSBpdSVFQUFBQUFxdXv+XDDz90dHR82uvLy8sdHBwavDsrhCOmr4qKCnt7e73esu0yefM3xt9J/Lyn4GXXSHWZr8rKSrVaTVFoEeuqqqpKpVLVn0KC59JoNEqlkmF0bbwolcrnHl6DzhGePHkyKipK+7h79+5XrlwpKyt7RnQ9l0KhCAwMrH+qUqme8dMyDKP7sQCCI6a/BhyxYd4kyoN8eoIKSWVm+JG3OxOr6pRqjxiCUHfaI4Yg1B3zBx1fr8tvo0FBeOvWrSZNmmgfa/OvtLTUkCBUqVSvvvqqji9WKBQKhaLB+7JCOGL6atgRUyjIB0FkUgdxVh4fsJV8EcLEWM2cUu0RQxDqTnvEEIS64zhOoVAY9zu9QUff2dm5vLxc+/jevXuEkGbNmhmhKAD5a+NAbR3IftaLnrGHH5rOXcQ6pQDmyqAgbNeuXWFhofZxYWGhi4uLk5OTMaoCsBBDveiTCWykO90zGVffA5gpXYPwxo0b58+f12g0paWl58+fr6ysJIRMmjRp/fr1J0+erK6uXrRo0aRJkxqxUgB5UtJkVmc6P549f4903silXcLQEMC86BqE8+fPj4qKunnz5sqVK6Oiovbs2UMICQwMfP/99/v37+/m5mZvb//ee+81ZqkAMuZpR60OY77uw/xtPz80nSsqRxwCmAv9Lp9oVEVFReHh4RcuXNDx9bgYQF84YvpqwOUTz1UrkK8LhQ+P8K/5W+DV97h8Ql+4fEJf+l4+oQscfQCTeqRTmopOKYDUEIQAEtB2Sr8JZd5EpxRAaghCAMlEuFMFI9hId/qFzdy8fL4ac0oBpIAgBJCS4oFOaRd0SgGkgCAEkN4jndIL6JQCmBCCEMBc1HdKe6FTCmBCCEIAM6LtlB7+Y07pNnRKARofghDA7HjYUavDmOWhzNvolAI0PgQhgJka4E4dQacUoPEhCAHM1yOd0pRiDA0BjA9BCGDutJ3S//Vl3jmATimA8SEIAeQhvOX9TmnwFnRKAYwJQQggG+iUAjQGBCGAzLir73dK5xzko9K4U3cRhwAGQRACyFJ4S+pwPDukFR26lZuVx1dyUhcEIFsIQgC50nZKj45k79QQv/Xc6jOC1BUByBKCEEDetJ3S1WHMoqMCOqUADYAgBLAEYQ93SivqpC4IQD4QhAAW4qFO6QZ0SgF0hSAEsCjaTukPYcwnR4XIVO4kOqUAz4MgBLBAYS2p/Hh2qBfdF51SgOdBEAJYJm2n9NhIBTqlAM+GIASwZC3VZHUY82MY88lRIQKdUoAnQRACWL7+LanD8ewwdEoBngRBCGAVWHRKAZ4CQQhgReo7pZ8eEyJSuUJ0SgEQhABWqH9LKj+OHeZF99vKzcrjy9EpBeuGIASwRuiUAtRDEAJYL22ndG048+kxYcA2dErBSiEIAaxdvxZUfhw73BudUrBSCEIAuN8pPT4KnVKwRghCALivhS1ZHcb8FM58dkwYsI07cQedUrAKCEIAeEjfFtShOHa4Nx22DZ1SsAoIQgB4FOaUglVBEALAkz3YKQ1HpxQsF4IQAJ5F2ykd15aOTEWnFCwTghAAnoOlycsd6eOjFIQQ/w3c6jMCxoZgSRCEAKCTZjZkSQiTFMl8WSgM2MYdR6cULAWCEN8/cSUAABkLSURBVAD00NOV2jeMndyejkrlZuXx9+rIrRpy9t4TQnH/DSQlyAOCEAD0Q1Mk0ffPTum3p/jRWfzpsodib/ExYfFxAS1UkAUEIQA0hLZTmhzFJBWJLE2GZ/D1S5V+flzYd0NcE8bQlLQ1AugEQQgADRfkQuUNY1/3p0urxbAU7tBt+vPjQt51cW04w+KvC8gEK3UBACBv2k5pjCc9fQ8Xkc4GuYq7hrBIQZAR/LYCgBG42ZLQ5nRXZ/G3UjGpCCvRgJxgRAgARqA9L5gdVbfximpCLi8SMqYNvmeDPCAIAcBQ9ecFazRkcntaINRfcnmaIgk+yEKQAQQhABhkz3Vx3437s2NqCCGETO1Aa3jxLzm8jz0V5IqZo2DuEIQAYJA+zakQt0evlJjpz7jbUsMyuPRYtrMzshDMGoIQAAz1xOsFR/rQIiHRaXzGIMbfCVkI5gtBCACNZZQPXcWRyFQuexDbEVkI5gpBCACNKNGXFkQSvZ3PGcy0cUAWgjlCEAJA45rUnuZFMmAbnzOY8UEWgvlBEAJAo5vSgRZEMiCVzx3MeNsjC8G84CofADCFaR3pWZ3oqDT+ShXuSQHmBUEIACby18709I50+Db+apXUpQA8AEEIAKbzty701A70gFTumkbqUgD+gCAEAJN6qyud4ENFp3G3aqQuBYAQgiAEANP7IJAZ4kVFpnK3kYVgBhCEACCBD4OYaE8qMpW7gywEqZkoCK9du7Zr164jR46IIiaMAQAhhCzoyUR6UIN3cOV1UpcC1s0U1xFWV1cnJiZ27ty5uLi4tLQ0MzOTZXH9IgCQhS8wr+/lY7ZzO2JYe4XU1YC1MkUgqVSq9PR07eOgoKAzZ874+fmZYL8AYOYoQr7szby2h4/dwW2PYe3wDRmkYKLWqCAImZmZ33zzTZMmTdq2bWuanQKA+aMI+aoP4+dExWVw1bzU1YBVMuYXMI7j5syZ88jGTz/9lBDC83xmZubFixdbtmxpxD0CgAWgCPm6DzPpV354Ord5IKtipC4IrIwRgrC2tpYQolQqaZoePHjwE1+jUCgWLFhACBk9enRaWtrw4cMN3y8AWAyGIt/3Zybm8vEZXHIUa4MsBBPSqTXK8/zLL78cGBjYtGnTQ4cO1W+vq6ubPHmyq6urq6vr9OnTRVEMfwwhpLy8vLq6mhBSXV1dVFTk5ubWSD8MAMgXQ5EfwhhHJTUuh68TpK4GrImu5wjbtWv3ySef1NXVcRxXv3HFihVHjhwpKSkpLi7evXv3mjVrnvje4uLi0NDQwMDAF154YcSIESEhIUYoHAAsjjYLeZGMy+E5ZCGYCqXXhX2Ojo7p6em9evXSPu3du/dLL700depUQsh///vf5OTkrKysBpdSVFQUEBDQvXv3+i2rVq1ydnZ+2usrKirs7e0bvDsrhCOmr8rKSjs7O6mrkJOqqipbW1uKMuhGS7UCmbCbVbPUiuA61tLX/KiqqlKpVDRt6T+n8Wg0GqVSyTC6ds/VavVzD69B5wjPnj1bfyGEn5/f4sWLDfk0Qoharf7HP/5R/9TFxcXGxuZpL+Z5Xq1WG7hHq4Ijpi9BEHDE9CKKolqtNjAI1YQkRZERWfysfNWKvjRt6bcvRBDqhaIovYJQl2NrUBCWlZXVf1+2t7e/c+eOIZ9GCFEoFFFRUTq+mKZp/PboBUdMXzhi+tIeMQODkBBiS5ONkfSQHdwre8T/9WUsOAvpP0hdiGw0xhEz6LNcXV3v3r2rfXznzh3MggEAY1GzJCWaPVcuvrKbx8KM0KgMCkJ/f//8/Hzt48OHD/v7+xujJAAAQghRs2RbNPt7mfjXPFxpD41I19bo3r17q6qqOI47ePBgeXl5nz59bG1tX3nlldmzZw8YMIDjuC+//HLFihWNWisAWBs7lmwdyA5M42bv4z8PxtWF0Ch0nTU6adKkkpKS+qc//vhj8+bNCSH/+c9/vvvuO5qmX3nllVdffdWQUoqKisLDwy9cuKDj68vLyx0cHAzZo7XBEdMX5tnqq7Ky0vDJMo+7W0uiUrn+LalPe1laFmLWqL70nTWqC/0un2hUCMLGhiOmLwShvhopCAkhd2tJRCo31IuaF2BRWYgg1FdjBCGOPgDIgJOSbI9hN1wQ5x/GlfZgZAhCAJAHVxXJHsT+dE74uABZCMaEIAQA2XCzJdmD2dVnhEVHkYVgNAhCAJCT5rYkI5ZZfkpYWogsBOPADaEBQGY87ajsQUz4Np6myKt++DYPhkIQAoD8eNlT2YOZsG08Q5GXOyILwSAIQgCQJW97Kj2GGZDK27JkYjtkITQcghAA5MrXkcoaxESk8gxFxrdFFkIDIQgBQMbaO1IZsUxkGk8TMhZZCA2CIAQAeevoRKVGM9FpnC1LhnsjC0FvCEIAkL2uTam0GDZ2O0dTZKgXshD0g98YALAE3ZtR26LZabv4bZfMZf1kkAsEIQBYiAAXalMU+9JOLvsKshD0gCAEAMsR4kYlRbLjcrjcq8hC0BWCEAAsSp/m1MYIdkw2t/MashB0giAEAEsT2oL6MYwdncXtv4EshOdDEAKABYryoFaHscMzuIM3kYXwHAhCALBMAz2olf3YoencoVJkITwLghAALFZsK+rrPszQdO74HWQhPBUuqAcASxbfmuZFEp3Gp8cynZwpqcsBc4QgBAALN8qHFkQyMI3PHMT4OSEL4VEIQgCwfKPb0BqexGznswcxbZsgC+EhCEIAsAov+tK8SAak8rmDGR8HZCH8CUEIANbipfa0IJLwbXzuYKY1shD+gCAEACsytQNdxZGoND53MONhhywEQhCEAGBt3uh0v0eaM5hxVyMLAdcRAoD1md2ZfrkjPWAbf00jdSlgBhCEAGCN/t6FHtuWGpjGlVZLXQpIDUEIAFZqXgAzzIuKTOVu1UhdCkgKQQgA1mt+EDOoFRWVyt1GFloxBCEAWLWPejJRHtTgHdy9OqlLAYkgCAHA2i14gQl0oWK3c+XIQquEIAQAa0cR8t/eTLemVOx2rgJZaH0QhAAAhCJkaR+mkzMVn8lpOKmrAdNCEAIAEEIIRcjXfZgWttTwDK6al7oaMCEEIQDAfTRFvu/PuKqo+AyuBlloNRCEAAB/YiiyOoxxsqHGZvN1gtTVgEkgCAEAHsJQZHV/RiBkbDbPIQutAIIQAOBRCpqsj2BqBXFcDrLQ8iEIAQCeQEmTjZFsFSdO2cULotTVQGNCEAIAPJmSJhsi2MuVyEILhyAEAHgqW5ZsHcgWlYsv7+YRhZYKQQgA8CxqlqREs2fKxFl5yELLhCAEAHgOO5ZsjWb33xD/tg9XF1ogBCEAwPM1UZD0WHbPdfHv+5GFlgZBCACgE0clSY9lf70qvncIWWhREIQAALpyUpK0GDapSPz3YVxdaDkQhAAAenBVkaxB7M/nhI+OIAstBIIQAEA/brYkezD741lhYQGy0BIgCAEA9NbclqTHMv/7Xfj0GLJQ9hCEAAAN4WlH5Qxmvi4UvipEFsobK3UBAABy1cqOyhnMhG3jGZq80hHjCrlCEAIANJyXPZUey4Rv42lCpiEL5QlBCABgkHZNqOxBzIBUXs2SCe2QhfKDIAQAMJSvI5U5iIlI5WmKjGuLLJQZBCEAgBF0cKTSopmoNI6hyOg2yEI5QRACABhHl6bUjlg2Jo1TMWSYN7JQNhCEAABG060plRrDxm7nGJoa3IqSuhzQCb6zAAAYU49mVHIUO+lXLu0Sbl8oDyYNwp07d169etWUewQAML1gN2pzFDtpJ5dzFVkoA6YLwvXr18fFxaWnp5tsjwAAUundnNoYwY7N5n5FFpo9EwXhrVu3vv322zFjxphmdwAAkgttQa0JY0dlcbuuIQvNmomC8K9//etHH32kUChMszsAAHMQ6UGtCWcTsrgDN5GF5suYs0aPHTv23nvvPbjFz8/vww8/TElJ8fDwCAwMXLVqlRF3BwBg/gZ6UN/3Z4elcynRbJAL5pGaIz2CkOO48vJyZ2fnRzaeOnXK1dW1efPmvr6+n3/++YP/qlKpCCFLlizRaDRRUVG///57dnZ2y5YtBw4caJTqAQDMX4wn9U0oM2g7tz2GDUAWmh+dgvDEiRNTpkwpKCjgOK6urq5++6lTpwYNGuTo6FhSUjJ58uSFCxd6e3s//vakpCSO4wgh77zzTo8ePcLCwoxUPACAPMR504JIhqRzGbFsJ2dkoXnR6Ryhs7PzvHnzUlNTH9n+zjvvjBkz5vDhwwUFBStWrMjPz3/i2x0cHJydnZ2dnbt169auXTulUmlo1QAAcjOiNb0khBmYxhfexflC80KJoq7/SY4dOxYQEFA/ItS2Sc+cOePj40MImTx5squr66JFixpcSlFRUbdu3Zo2bVq/JSsry9XV9Wmvr6iosLe3b/DurBCOmL4qKyvt7OykrkJOqqqqbG1tKQojnqf6qYh57yib1K+2s5NICNFoNDY2NjRNE0KqeaJipK7P7Gk0GqVSyTC6Him1Wv3cFzd8skxJSQkhpL4X2qZNmxMnTjT407QcHR2zsrLqn3p7ez/7B3BwcDBwj9YGR0wvFEXhq4NeaJpWq9UIwmd4uQu5WscPyLQ5GM92dqIYhlGpVDRNH7klztjL/zqYZbHe1zOxLKtXEOr0mQ1+Z2VlZf0XGUKIWq0uLy83sBqGYdq0aWPghwAAmLP3ApjrGtJzE3cwjm1jQwghBbfFl3by6yMYpKAkGh6Ebm5uVVVV1dXV2qmht27datGihfEKAwCwWF/1YWp40jOZ2x1D0VXilF3C+gimbROMpKXR8K8f7u7uLVu2zMvL0z7Ny8sLDAw0UlUAABZuRT8mwYcO3c6OzRY2RiIFpaTTiFCj0fzwww8lJSWCICxfvtzOzm7ChAkMw8yYMWP27NmLFy8+cODA8ePHk5OTG7tcAACL8bcu9I7L/MVKKj6Df6kDPb4t7aKSuiarpFMQchx36NAhQsjUqVMPHTrk5OQ0YcIEQsicOXPUavUnn3zi6uqam5vr6OjYuMUCAFgK7XnBrCjuRIXNwqPikVvi3EN1PV2plzvScd60AicLTUiPyycaW1FRUXh4+IULF3R8fXl5OeZA6gVHTF+44ERflZWVmDWqi4Lb4uRf+fURTEtWo1Kpki6SpYXCmnA2s0T44axwuFQc6UO/0pHGMjSP0/fyCV3gDvUAACZ16q44+Vd+YyTj40BVVRFCyCgfWhDJxFxuRwyb6EsXV4g/nRMTsngVQxJ96cntaTdbqYu2aAhCAACT8ranNg9kWtk9NNob3Ybu3ozSXj7hZU+90416qyu997r4w1nBb0Nd7+ZUoi9apo0FBxUAwKRsWfJICmq1d3xoI02R0BbUslDmwlhFgg+9/JTQck3dK7v5w7fM5XyWxcCIEADArDVRkERfur5lOjKTt2VIoi/9UgfaFbNMjQEjQgAAefCyp97pRp8dzS4LZc6Xix3X1w1N59ZfEOoEqSuTOQQhAICcPN4y9f65blYefwQt04ZCaxQAQJbqW6any8S154QRaJk2FEaEAADy1t6RmhfA1LdMO6BlqicEIQCAJfizZToGLVP9oDUKAGBRHJX3W6a/l4k/nRPiM3k1QxJ96SkdsJbpk2FECABgmTo4UvMCmHOj2SUhzIk7Yvs/WqYcWqYPw4gQAMCS0RSJ9KAiPZiyWmbzRWH5KeGvecIoH+qlDnS3pljLlBAEIQCAlXikZTo8nXdRkYnt6AntrL1litYoAIB10bZMz49hF/RkDpWiZYoRIQCAVULLtB6CEADAqtW3TE/dFX8+b40tU7RGAQCAEEI6Oj3aMh2dxW8ttvyWKUaEAADwp/qW6d1a5pfzwsICYfpuYZQPNaUD3dVCW6YYEQIAwBM4KcnLHendQ9msQYyzDRmazgclc0uOC7dqpK7M2BCEAADwLNqW6YU/Wqa+v1hayxStUQAAeL7HW6av7hFGtqamdqC7yLxlihEhAADoob5lmhHL2LJkYBon95YpghAAABrCz4la0JO5PF4h95YpWqMAANBwzJNapn9pR01pT/s6yqNlihEhAAAYwYMtU0JIv5T7LdPbZt8yRRACAIAxPdgy3XNd9Pn5fsuUN9c7BCMIAQDA+LQt018imKKxikgPamGB4P0zN+cgf6bM7PIQQQgAAI3I2eZ+yzQ9liGE9E3hgpK55aeE8jqpK/sDghAAAEzB34la0JMpGa9Y0JPJLBE915pLyxRBCAAApvN4y7T1z9ycg/zZe5LlIYIQAAAkUN8y3RHLEEJCt0rWMkUQAgCAlJ7YMs0sEetHiL+Xif8+/OhV+tU8mbGXr+GNUACCEAAApPdIy/Sdg7zXT/dbpu0dqWsa8e0Df4ZerUBGZ/FdnCkbxgi7RhACAIAZ0bZMD8Wx22MYQkifrVzPZK5rU+puLdFmYa1ARmXyg1pR0/2ME2FYYg0AAMxRJ2dqQU/mg0BmW7Hw/Rlx1zXBzZY6X0bXCOJgL9pYKUgwIgQAAHOmpEl8a3pzFHNqlOKl9vSOKzRNRCOmIEEQAgCALDjZkN3XxA97CB721IPnCw2HIAQAAHNXf15wmq/wRTBVyREjZiGCEAAAzNojs2MoQr7szRgxCxGEAABg1koqxTFtHpojqs3C5raUtV9HWFJSUltbK3UVcnLlypWaGrO/M5g5uXbtmkajkboKObl+/XpVVZXUVcjJzZs3KyoqpK7C3Pk4UBPa3U+r0tLS8vJyQghFyN+70NZ+HWFCQsKpU6ekrkJOJk6cePjwYamrkJNp06bt3btX6irkZObMmVlZWVJXISdvvfXWtm3bpK5CTt59990NGzYY9zNlHIQAAACGQxACAIBVQxACAIBVo0RR6lsi/uHWrVtDhgy5ceOGjq+/c+eOg4MDy2KVOF3hiOnr7t27dnZ2CoVC6kJko6yszNbWVqlUSl2IbNy7d8/GxsbGxkbqQmSjvLxcoVCoVCodX5+SkuLn5/fs15hREBJCbt68qZ0OBAAAYDhPT8/nfjMzryAEAAAwMZwjBAAAq4YgBAAAq4YgBAAAq4YgBAAAqybLmfTXr1//7bffSkpKIiIi2rZtK3U5MlBYWLh9+/aSkhIfH5/ExMQmTZpIXZG5y8nJ2b179507dzw9PSdOnOjq6ip1RbLx448/NmvWLDY2VupCzN3u3bsLCwvrn7788ssSFiMXd+/eXbVqVVFRkfZ/TDc3N6N8rCxHhP369fvoo4/eeeedgwcPSl2LPERHR58/f97Ly2vHjh09evQoKyuTuiJz98svv3Ac16ZNm4MHD3bt2vX69etSVyQPP/zww4wZM7744gupC5GBtWvXrlmz5vwfpC5HBq5evdqjR4/du3e3bt26pKTk9OnTxvpkWV4+IQgCTdPdu3efM2fO2LFjpS5HBqqrq7XXnwqC4Ovr+/HHH48ePVrqomTD39//3XffHTdunNSFmLtr165FRESMGDHit99+S0tLk7occ/faa695eHj885//lLoQ2XjxxRcVCsW3335r9E+W5YiQpmVZtoQeXIWhurra3t5ewmLk5dy5czdu3OjUqZPUhcjA66+//v777zs7O0tdiGwcOnRo0aJF69atw/3RdJGampqQkLB69eqlS5eePXvWiJ+MRLEuH330UfPmzQcOHCh1ITIwd+5cT09Pf3//Dz74oGvXrlKXY+5++eWX6urqUaNGSV2IbHh6erZo0eLu3bufffZZQEDAvXv3pK7IrJWXl5eWlr711lsnTpw4f/58z5499+zZY7RPF2WrW7duP/30k9RVyMnq1atbtmx56tQpqQuRh8rKymvXriUlJTVr1mzXrl1Sl2PWSktLfX19L126JIriZ599FhMTI3VFcsJxXFBQ0KJFi6QuxKxpb1/80UcfaZ/+61//Gjx4sLE+XJazRqEB1q1bN2fOnKysrA4dOkhdizyo1Wq1Wh0fH5+SkpKUlBQaGip1ReYrPT29tLQ0Li6OEHL9+vV79+7169dv586dUtclDwzDBAcHY77Ms9nZ2Tk7O/v7+2ufdurUKTk52VgfjtaoVUhKSpo9e/aOHTueuwo7EEJ4nq+trdU+rqurO3LkiJeXl7Qlmbno6OiMjIxly5YtW7YsISGhU6dOmDj6XBqNRvugqqoqMzMT56GfKy4ubt++fdrHeXl5Rjxispw1OnPmzLy8vMLCwhYtWjRt2vSbb74JCgqSuijzVVdXZ2dn5+Li4u7urt0yc+bMF198UdqqzFlpaam/v3/v3r0dHBz27Nnj6em5fft2tVotdV3ysHjx4oyMDMwafS4PD4+AgABHR8dff/21Y8eOKSkpuBnTs507d65///6hoaEcx+3fvz8rK6t9+/ZG+WRZBuGZM2cePLHcvn17BwcHCesxc6Io5ufnP7jFw8OjRYsWUtUjC5cvX87Pz9doNO3atQsMDJS6HDnRtkZ9fX2lLsTcXbp0KT8/v7q6Gr9jurt3715WVpZKperTp48RFwaRZRACAAAYC84RAgCAVUMQAgCAVUMQAgCAVUMQAgCAVUMQAgCAVUMQAgCAVUMQAgCAVUMQAgCAVUMQAgCAVUMQAgCAVUMQAgCAVft/4zQ18mH2NNMAAAAASUVORK5CYII=",
"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
}