-
Notifications
You must be signed in to change notification settings - Fork 89
/
index.html
158 lines (155 loc) · 37 KB
/
index.html
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
<!DOCTYPE html>
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>Monitoring self-consistent field calculations · DFTK.jl</title><link rel="canonical" href="https://juliamolsim.github.io/DFTK.jl/stable/examples/scf_callbacks/"/><link href="https://fonts.googleapis.com/css?family=Lato|Roboto+Mono" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/fontawesome.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/solid.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/font-awesome/5.15.0/css/brands.min.css" rel="stylesheet" type="text/css"/><link href="https://cdnjs.cloudflare.com/ajax/libs/KaTeX/0.11.1/katex.min.css" rel="stylesheet" type="text/css"/><script>documenterBaseURL="../.."</script><script src="https://cdnjs.cloudflare.com/ajax/libs/require.js/2.3.6/require.min.js" data-main="../../assets/documenter.js"></script><script src="../../siteinfo.js"></script><script src="../../../versions.js"></script><link class="docs-theme-link" rel="stylesheet" type="text/css" href="../../assets/themes/documenter-dark.css" data-theme-name="documenter-dark" data-theme-primary-dark/><link class="docs-theme-link" rel="stylesheet" type="text/css" href="../../assets/themes/documenter-light.css" data-theme-name="documenter-light" data-theme-primary/><script src="../../assets/themeswap.js"></script><link href="../../assets/favicon.ico" rel="icon" type="image/x-icon"/></head><body><div id="documenter"><nav class="docs-sidebar"><a class="docs-logo" href="../../"><img src="../../assets/logo.png" alt="DFTK.jl logo"/></a><div class="docs-package-name"><span class="docs-autofit">DFTK.jl</span></div><form class="docs-search" action="../../search/"><input class="docs-search-query" id="documenter-search-query" name="q" type="text" placeholder="Search docs"/></form><ul class="docs-menu"><li><a class="tocitem" href="../../">Home</a></li><li><span class="tocitem">Getting started</span><ul><li><a class="tocitem" href="../../guide/installation/">Installation</a></li><li><a class="tocitem" href="../../guide/tutorial/">Tutorial</a></li><li><a class="tocitem" href="../../guide/input_output/">Input and output formats</a></li><li><a class="tocitem" href="../../guide/parallelization/">Timings and parallelization</a></li><li><a class="tocitem" href="../../guide/density_functional_theory/">Density-functional theory</a></li></ul></li><li><span class="tocitem">Examples</span><ul><li><a class="tocitem" href="../metallic_systems/">Temperature and metallic systems</a></li><li><a class="tocitem" href="../pymatgen/">Creating supercells with pymatgen</a></li><li><a class="tocitem" href="../ase/">Creating slabs with ASE</a></li><li><a class="tocitem" href="../collinear_magnetism/">Collinear spin and magnetic systems</a></li><li><a class="tocitem" href="../geometry_optimization/">Geometry optimization</a></li><li class="is-active"><a class="tocitem" href>Monitoring self-consistent field calculations</a></li><li><a class="tocitem" href="../scf_checkpoints/">Saving SCF results on disk and SCF checkpoints</a></li><li><a class="tocitem" href="../polarizability/">Polarizability by linear response</a></li><li><a class="tocitem" href="../gross_pitaevskii/">Gross-Pitaevskii equation in one dimension</a></li><li><a class="tocitem" href="../gross_pitaevskii_2D/">Gross-Pitaevskii equation with magnetism</a></li><li><a class="tocitem" href="../cohen_bergstresser/">Cohen-Bergstresser model</a></li><li><a class="tocitem" href="../arbitrary_floattype/">Arbitrary floating-point types</a></li><li><a class="tocitem" href="../custom_solvers/">Custom solvers</a></li><li><a class="tocitem" href="../custom_potential/">Custom potential</a></li></ul></li><li><span class="tocitem">Advanced topics</span><ul><li><a class="tocitem" href="../../advanced/conventions/">Notation and conventions</a></li><li><a class="tocitem" href="../../advanced/data_structures/">Data structures</a></li><li><a class="tocitem" href="../../advanced/useful_formulas/">Useful formulas</a></li><li><a class="tocitem" href="../../advanced/symmetries/">Crystal symmetries</a></li></ul></li><li><a class="tocitem" href="../../api/">API reference</a></li><li><a class="tocitem" href="../../publications/">Publications</a></li></ul><div class="docs-version-selector field has-addons"><div class="control"><span class="docs-label button is-static is-size-7">Version</span></div><div class="docs-selector control is-expanded"><div class="select is-fullwidth is-size-7"><select id="documenter-version-selector"></select></div></div></div></nav><div class="docs-main"><header class="docs-navbar"><nav class="breadcrumb"><ul class="is-hidden-mobile"><li><a class="is-disabled">Examples</a></li><li class="is-active"><a href>Monitoring self-consistent field calculations</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href>Monitoring self-consistent field calculations</a></li></ul></nav><div class="docs-right"><a class="docs-edit-link" href="https://github.com/JuliaMolSim/DFTK.jl/blob/master/examples/scf_callbacks.jl" title="Edit on GitHub"><span class="docs-icon fab"></span><span class="docs-label is-hidden-touch">Edit on GitHub</span></a><a class="docs-settings-button fas fa-cog" id="documenter-settings-button" href="#" title="Settings"></a><a class="docs-sidebar-button fa fa-bars is-hidden-desktop" id="documenter-sidebar-button" href="#"></a></div></header><article class="content" id="documenter-page"><h1 id="Monitoring-self-consistent-field-calculations"><a class="docs-heading-anchor" href="#Monitoring-self-consistent-field-calculations">Monitoring self-consistent field calculations</a><a id="Monitoring-self-consistent-field-calculations-1"></a><a class="docs-heading-anchor-permalink" href="#Monitoring-self-consistent-field-calculations" title="Permalink"></a></h1><p><a href="https://mybinder.org/v2/gh/JuliaMolSim/DFTK.jl/gh-pages?filepath=v0.3.0/examples/scf_callbacks.ipynb"><img src="https://mybinder.org/badge_logo.svg" alt/></a> <a href="https://nbviewer.jupyter.org/github/JuliaMolSim/DFTK.jl/blob/gh-pages/v0.3.0/examples/scf_callbacks.ipynb"><img src="https://img.shields.io/badge/show-nbviewer-579ACA.svg" alt/></a></p><p>The <code>self_consistent_field</code> function takes as the <code>callback</code> keyword argument one function to be called after each iteration. This function gets passed the complete internal state of the SCF solver and can thus be used both to monitor and debug the iterations as well as to quickly patch it with additional functionality.</p><p>This example discusses a few aspects of the <code>callback</code> function taking again our favourite silicon example.</p><p>We setup silicon in an LDA model using the ASE interface to build the silicon lattice, see <a href="../ase/#Creating-slabs-with-ASE">Creating slabs with ASE</a> for more details.</p><pre><code class="language-julia">using DFTK
using PyCall
silicon = pyimport("ase.build").bulk("Si")
atoms = load_atoms(silicon)
atoms = [ElementPsp(el.symbol, psp=load_psp(el.symbol, functional="lda")) => position
for (el, position) in atoms]
lattice = load_lattice(silicon);
model = model_LDA(lattice, atoms)
kgrid = [3, 3, 3] # k-point grid
Ecut = 5 # kinetic energy cutoff in Hartree
basis = PlaneWaveBasis(model, Ecut; kgrid=kgrid);</code></pre><p>DFTK already defines a few callback functions for standard tasks. One example is the usual convergence table, which is defined in the callback <code>ScfDefaultCallback</code>. Another example is <code>ScfPlotTrace</code>, which records the total energy at each iteration and uses it to plot the convergence of the SCF graphically once it is converged. For details and other callbacks see <a href="https://dftk.org/blob/master/src/scf/scf_callbacks.jl"><code>src/scf/scf_callbacks.jl</code></a>.</p><div class="admonition is-info"><header class="admonition-header">Callbacks are not exported</header><div class="admonition-body"><p>Callbacks are not exported from the DFTK namespace as of now, so you will need to use them, e.g., as <code>DFTK.ScfDefaultCallback</code> and <code>DFTK.ScfPlotTrace</code>.</p></div></div><p>In this example we define a custom callback, which plots the change in density at each SCF iteration after the SCF has finished. For this we first define the empty plot canvas and an empty container for all the density differences:</p><pre><code class="language-julia">using Plots
p = plot(yaxis=:log)
density_differences = Float64[];</code></pre><p>The callback function itself gets passed a named tuple similar to the one returned by <code>self_consistent_field</code>, which contains the input and output density of the SCF step as <code>ρin</code> and <code>ρout</code>. Since the callback gets called both during the SCF iterations as well as after convergence just before <code>self_consistent_field</code> finishes we can both collect the data and initiate the plotting in one function.</p><pre><code class="language-julia">using LinearAlgebra
function plot_callback(info)
if info.stage == :finalize
plot!(p, density_differences, label="|ρout - ρin|", markershape=:x)
else
push!(density_differences, norm(info.ρout - info.ρin))
end
info
end
callback = DFTK.ScfDefaultCallback() ∘ plot_callback;</code></pre><p>Notice that for constructing the <code>callback</code> function we chained the <code>plot_callback</code> (which does the plotting) with the <code>ScfDefaultCallback</code>, such that when using the <code>plot_callback</code> function with <code>self_consistent_field</code> we still get the usual convergence table printed. We run the SCF with this callback ...</p><pre><code class="language-julia">scfres = self_consistent_field(basis, tol=1e-8, callback=callback);</code></pre><pre class="documenter-example-output">n Energy Eₙ-Eₙ₋₁ ρout-ρin Diag
--- --------------- --------- -------- ----
1 -7.844829555675 NaN 1.98e-01 4.0
2 -7.850328788370 -5.50e-03 2.93e-02 1.0
3 -7.850646851250 -3.18e-04 3.00e-03 3.0
4 -7.850647501898 -6.51e-07 4.63e-04 2.8
5 -7.850647510750 -8.85e-09 1.61e-05 1.0</pre><p>... and show the plot</p><pre><code class="language-julia">p</code></pre><?xml version="1.0" encoding="utf-8"?>
<svg xmlns="http://www.w3.org/2000/svg" xmlns:xlink="http://www.w3.org/1999/xlink" width="600" height="400" viewBox="0 0 2400 1600">
<defs>
<clipPath id="clip880">
<rect x="0" y="0" width="2400" height="1600"/>
</clipPath>
</defs>
<path clip-path="url(#clip880)" d="
M0 1600 L2400 1600 L2400 0 L0 0 Z
" fill="#ffffff" fill-rule="evenodd" fill-opacity="1"/>
<defs>
<clipPath id="clip881">
<rect x="480" y="0" width="1681" height="1600"/>
</clipPath>
</defs>
<path clip-path="url(#clip880)" d="
M189.496 1486.45 L2352.76 1486.45 L2352.76 47.2441 L189.496 47.2441 Z
" fill="#ffffff" fill-rule="evenodd" fill-opacity="1"/>
<defs>
<clipPath id="clip882">
<rect x="189" y="47" width="2164" height="1440"/>
</clipPath>
</defs>
<polyline clip-path="url(#clip882)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:2; stroke-opacity:0.1; fill:none" points="
250.72,1486.45 250.72,47.2441
"/>
<polyline clip-path="url(#clip882)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:2; stroke-opacity:0.1; fill:none" points="
760.923,1486.45 760.923,47.2441
"/>
<polyline clip-path="url(#clip882)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:2; stroke-opacity:0.1; fill:none" points="
1271.13,1486.45 1271.13,47.2441
"/>
<polyline clip-path="url(#clip882)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:2; stroke-opacity:0.1; fill:none" points="
1781.33,1486.45 1781.33,47.2441
"/>
<polyline clip-path="url(#clip882)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:2; stroke-opacity:0.1; fill:none" points="
2291.53,1486.45 2291.53,47.2441
"/>
<polyline clip-path="url(#clip880)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
189.496,1486.45 2352.76,1486.45
"/>
<polyline clip-path="url(#clip880)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
250.72,1486.45 250.72,1469.18
"/>
<polyline clip-path="url(#clip880)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
760.923,1486.45 760.923,1469.18
"/>
<polyline clip-path="url(#clip880)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
1271.13,1486.45 1271.13,1469.18
"/>
<polyline clip-path="url(#clip880)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
1781.33,1486.45 1781.33,1469.18
"/>
<polyline clip-path="url(#clip880)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
2291.53,1486.45 2291.53,1469.18
"/>
<path clip-path="url(#clip880)" d="M241.102 1543.18 L248.741 1543.18 L248.741 1516.82 L240.431 1518.49 L240.431 1514.23 L248.695 1512.56 L253.371 1512.56 L253.371 1543.18 L261.01 1543.18 L261.01 1547.12 L241.102 1547.12 L241.102 1543.18 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M755.576 1543.18 L771.895 1543.18 L771.895 1547.12 L749.951 1547.12 L749.951 1543.18 Q752.613 1540.43 757.196 1535.8 Q761.803 1531.15 762.983 1529.81 Q765.229 1527.28 766.108 1525.55 Q767.011 1523.79 767.011 1522.1 Q767.011 1519.34 765.067 1517.61 Q763.145 1515.87 760.043 1515.87 Q757.844 1515.87 755.391 1516.63 Q752.96 1517.4 750.182 1518.95 L750.182 1514.23 Q753.006 1513.09 755.46 1512.51 Q757.914 1511.93 759.951 1511.93 Q765.321 1511.93 768.516 1514.62 Q771.71 1517.31 771.71 1521.8 Q771.71 1523.93 770.9 1525.85 Q770.113 1527.74 768.006 1530.34 Q767.428 1531.01 764.326 1534.23 Q761.224 1537.42 755.576 1543.18 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M1275.37 1528.49 Q1278.73 1529.2 1280.61 1531.47 Q1282.5 1533.74 1282.5 1537.07 Q1282.5 1542.19 1278.98 1544.99 Q1275.47 1547.79 1268.98 1547.79 Q1266.81 1547.79 1264.49 1547.35 Q1262.2 1546.93 1259.75 1546.08 L1259.75 1541.56 Q1261.69 1542.7 1264.01 1543.28 Q1266.32 1543.86 1268.85 1543.86 Q1273.24 1543.86 1275.54 1542.12 Q1277.85 1540.38 1277.85 1537.07 Q1277.85 1534.02 1275.7 1532.31 Q1273.57 1530.57 1269.75 1530.57 L1265.72 1530.57 L1265.72 1526.73 L1269.93 1526.73 Q1273.38 1526.73 1275.21 1525.36 Q1277.04 1523.97 1277.04 1521.38 Q1277.04 1518.72 1275.14 1517.31 Q1273.27 1515.87 1269.75 1515.87 Q1267.83 1515.87 1265.63 1516.29 Q1263.43 1516.7 1260.79 1517.58 L1260.79 1513.42 Q1263.45 1512.68 1265.77 1512.31 Q1268.11 1511.93 1270.17 1511.93 Q1275.49 1511.93 1278.59 1514.37 Q1281.69 1516.77 1281.69 1520.89 Q1281.69 1523.76 1280.05 1525.75 Q1278.41 1527.72 1275.37 1528.49 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M1784.34 1516.63 L1772.53 1535.08 L1784.34 1535.08 L1784.34 1516.63 M1783.11 1512.56 L1788.99 1512.56 L1788.99 1535.08 L1793.92 1535.08 L1793.92 1538.97 L1788.99 1538.97 L1788.99 1547.12 L1784.34 1547.12 L1784.34 1538.97 L1768.74 1538.97 L1768.74 1534.46 L1783.11 1512.56 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M2281.81 1512.56 L2300.17 1512.56 L2300.17 1516.5 L2286.09 1516.5 L2286.09 1524.97 Q2287.11 1524.62 2288.13 1524.46 Q2289.15 1524.27 2290.17 1524.27 Q2295.95 1524.27 2299.33 1527.44 Q2302.71 1530.62 2302.71 1536.03 Q2302.71 1541.61 2299.24 1544.71 Q2295.77 1547.79 2289.45 1547.79 Q2287.27 1547.79 2285 1547.42 Q2282.76 1547.05 2280.35 1546.31 L2280.35 1541.61 Q2282.43 1542.74 2284.66 1543.3 Q2286.88 1543.86 2289.36 1543.86 Q2293.36 1543.86 2295.7 1541.75 Q2298.04 1539.64 2298.04 1536.03 Q2298.04 1532.42 2295.7 1530.31 Q2293.36 1528.21 2289.36 1528.21 Q2287.48 1528.21 2285.61 1528.62 Q2283.75 1529.04 2281.81 1529.92 L2281.81 1512.56 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><polyline clip-path="url(#clip882)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:2; stroke-opacity:0.1; fill:none" points="
189.496,1403.6 2352.76,1403.6
"/>
<polyline clip-path="url(#clip882)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:2; stroke-opacity:0.1; fill:none" points="
189.496,1071.61 2352.76,1071.61
"/>
<polyline clip-path="url(#clip882)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:2; stroke-opacity:0.1; fill:none" points="
189.496,739.614 2352.76,739.614
"/>
<polyline clip-path="url(#clip882)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:2; stroke-opacity:0.1; fill:none" points="
189.496,407.619 2352.76,407.619
"/>
<polyline clip-path="url(#clip882)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:2; stroke-opacity:0.1; fill:none" points="
189.496,75.6241 2352.76,75.6241
"/>
<polyline clip-path="url(#clip880)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
189.496,1486.45 189.496,47.2441
"/>
<polyline clip-path="url(#clip880)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
189.496,1403.6 215.455,1403.6
"/>
<polyline clip-path="url(#clip880)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
189.496,1071.61 215.455,1071.61
"/>
<polyline clip-path="url(#clip880)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
189.496,739.614 215.455,739.614
"/>
<polyline clip-path="url(#clip880)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
189.496,407.619 215.455,407.619
"/>
<polyline clip-path="url(#clip880)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
189.496,75.6241 215.455,75.6241
"/>
<path clip-path="url(#clip880)" d="M51.6634 1423.4 L59.3023 1423.4 L59.3023 1397.03 L50.9921 1398.7 L50.9921 1394.44 L59.256 1392.77 L63.9319 1392.77 L63.9319 1423.4 L71.5707 1423.4 L71.5707 1427.33 L51.6634 1427.33 L51.6634 1423.4 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M91.0151 1395.85 Q87.404 1395.85 85.5753 1399.41 Q83.7697 1402.96 83.7697 1410.09 Q83.7697 1417.19 85.5753 1420.76 Q87.404 1424.3 91.0151 1424.3 Q94.6493 1424.3 96.4548 1420.76 Q98.2835 1417.19 98.2835 1410.09 Q98.2835 1402.96 96.4548 1399.41 Q94.6493 1395.85 91.0151 1395.85 M91.0151 1392.15 Q96.8252 1392.15 99.8808 1396.75 Q102.959 1401.34 102.959 1410.09 Q102.959 1418.81 99.8808 1423.42 Q96.8252 1428 91.0151 1428 Q85.2049 1428 82.1262 1423.42 Q79.0707 1418.81 79.0707 1410.09 Q79.0707 1401.34 82.1262 1396.75 Q85.2049 1392.15 91.0151 1392.15 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M102.959 1386.25 L127.071 1386.25 L127.071 1389.44 L102.959 1389.44 L102.959 1386.25 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M145.71 1375.15 L136.118 1390.14 L145.71 1390.14 L145.71 1375.15 M144.713 1371.84 L149.49 1371.84 L149.49 1390.14 L153.496 1390.14 L153.496 1393.3 L149.49 1393.3 L149.49 1399.92 L145.71 1399.92 L145.71 1393.3 L133.033 1393.3 L133.033 1389.63 L144.713 1371.84 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M52.585 1091.4 L60.2238 1091.4 L60.2238 1065.04 L51.9137 1066.7 L51.9137 1062.44 L60.1776 1060.78 L64.8535 1060.78 L64.8535 1091.4 L72.4923 1091.4 L72.4923 1095.34 L52.585 1095.34 L52.585 1091.4 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M91.9366 1063.85 Q88.3255 1063.85 86.4969 1067.42 Q84.6913 1070.96 84.6913 1078.09 Q84.6913 1085.2 86.4969 1088.76 Q88.3255 1092.3 91.9366 1092.3 Q95.5709 1092.3 97.3764 1088.76 Q99.2051 1085.2 99.2051 1078.09 Q99.2051 1070.96 97.3764 1067.42 Q95.5709 1063.85 91.9366 1063.85 M91.9366 1060.15 Q97.7468 1060.15 100.802 1064.76 Q103.881 1069.34 103.881 1078.09 Q103.881 1086.82 100.802 1091.42 Q97.7468 1096.01 91.9366 1096.01 Q86.1265 1096.01 83.0478 1091.42 Q79.9923 1086.82 79.9923 1078.09 Q79.9923 1069.34 83.0478 1064.76 Q86.1265 1060.15 91.9366 1060.15 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M103.881 1054.25 L127.993 1054.25 L127.993 1057.45 L103.881 1057.45 L103.881 1054.25 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M147.703 1052.79 Q150.43 1053.37 151.954 1055.21 Q153.496 1057.05 153.496 1059.76 Q153.496 1063.92 150.637 1066.2 Q147.778 1068.47 142.512 1068.47 Q140.744 1068.47 138.863 1068.11 Q137.002 1067.78 135.008 1067.08 L135.008 1063.41 Q136.588 1064.33 138.469 1064.8 Q140.349 1065.27 142.399 1065.27 Q145.973 1065.27 147.835 1063.86 Q149.716 1062.45 149.716 1059.76 Q149.716 1057.28 147.966 1055.89 Q146.236 1054.48 143.133 1054.48 L139.86 1054.48 L139.86 1051.36 L143.283 1051.36 Q146.086 1051.36 147.571 1050.25 Q149.057 1049.12 149.057 1047.01 Q149.057 1044.85 147.515 1043.7 Q145.992 1042.54 143.133 1042.54 Q141.572 1042.54 139.785 1042.87 Q137.998 1043.21 135.854 1043.93 L135.854 1040.54 Q138.017 1039.94 139.898 1039.64 Q141.797 1039.34 143.471 1039.34 Q147.797 1039.34 150.317 1041.31 Q152.838 1043.27 152.838 1046.62 Q152.838 1048.95 151.502 1050.57 Q150.167 1052.16 147.703 1052.79 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M53.3561 759.406 L60.995 759.406 L60.995 733.041 L52.6848 734.707 L52.6848 730.448 L60.9487 728.781 L65.6246 728.781 L65.6246 759.406 L73.2634 759.406 L73.2634 763.341 L53.3561 763.341 L53.3561 759.406 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M92.7078 731.86 Q89.0967 731.86 87.268 735.425 Q85.4624 738.967 85.4624 746.096 Q85.4624 753.203 87.268 756.767 Q89.0967 760.309 92.7078 760.309 Q96.342 760.309 98.1475 756.767 Q99.9762 753.203 99.9762 746.096 Q99.9762 738.967 98.1475 735.425 Q96.342 731.86 92.7078 731.86 M92.7078 728.156 Q98.5179 728.156 101.573 732.763 Q104.652 737.346 104.652 746.096 Q104.652 754.823 101.573 759.429 Q98.5179 764.013 92.7078 764.013 Q86.8976 764.013 83.8189 759.429 Q80.7634 754.823 80.7634 746.096 Q80.7634 737.346 83.8189 732.763 Q86.8976 728.156 92.7078 728.156 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M104.652 722.258 L128.764 722.258 L128.764 725.455 L104.652 725.455 L104.652 722.258 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M140.236 732.734 L153.496 732.734 L153.496 735.931 L135.666 735.931 L135.666 732.734 Q137.829 730.496 141.553 726.734 Q145.296 722.954 146.255 721.863 Q148.079 719.813 148.794 718.402 Q149.527 716.973 149.527 715.6 Q149.527 713.362 147.948 711.951 Q146.387 710.541 143.866 710.541 Q142.08 710.541 140.086 711.161 Q138.111 711.782 135.854 713.042 L135.854 709.205 Q138.149 708.284 140.142 707.813 Q142.136 707.343 143.791 707.343 Q148.155 707.343 150.75 709.525 Q153.345 711.707 153.345 715.355 Q153.345 717.086 152.687 718.647 Q152.048 720.189 150.336 722.295 Q149.866 722.841 147.346 725.455 Q144.826 728.051 140.236 732.734 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M53.0552 427.411 L60.694 427.411 L60.694 401.046 L52.3839 402.712 L52.3839 398.453 L60.6477 396.787 L65.3236 396.787 L65.3236 427.411 L72.9625 427.411 L72.9625 431.347 L53.0552 431.347 L53.0552 427.411 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M92.4068 399.865 Q88.7957 399.865 86.967 403.43 Q85.1615 406.972 85.1615 414.101 Q85.1615 421.208 86.967 424.772 Q88.7957 428.314 92.4068 428.314 Q96.0411 428.314 97.8466 424.772 Q99.6753 421.208 99.6753 414.101 Q99.6753 406.972 97.8466 403.43 Q96.0411 399.865 92.4068 399.865 M92.4068 396.162 Q98.217 396.162 101.273 400.768 Q104.351 405.351 104.351 414.101 Q104.351 422.828 101.273 427.435 Q98.217 432.018 92.4068 432.018 Q86.5967 432.018 83.518 427.435 Q80.4625 422.828 80.4625 414.101 Q80.4625 405.351 83.518 400.768 Q86.5967 396.162 92.4068 396.162 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M104.351 390.263 L128.463 390.263 L128.463 393.46 L104.351 393.46 L104.351 390.263 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M137.321 400.739 L143.528 400.739 L143.528 379.317 L136.776 380.671 L136.776 377.21 L143.49 375.856 L147.289 375.856 L147.289 400.739 L153.496 400.739 L153.496 403.936 L137.321 403.936 L137.321 400.739 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M82.7903 95.4165 L90.4291 95.4165 L90.4291 69.0509 L82.119 70.7176 L82.119 66.4583 L90.3828 64.7917 L95.0587 64.7917 L95.0587 95.4165 L102.698 95.4165 L102.698 99.3517 L82.7903 99.3517 L82.7903 95.4165 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M122.142 67.8703 Q118.531 67.8703 116.702 71.4351 Q114.897 74.9768 114.897 82.1064 Q114.897 89.2128 116.702 92.7776 Q118.531 96.3193 122.142 96.3193 Q125.776 96.3193 127.582 92.7776 Q129.41 89.2128 129.41 82.1064 Q129.41 74.9768 127.582 71.4351 Q125.776 67.8703 122.142 67.8703 M122.142 64.1667 Q127.952 64.1667 131.008 68.7731 Q134.086 73.3564 134.086 82.1064 Q134.086 90.8332 131.008 95.4396 Q127.952 100.023 122.142 100.023 Q116.332 100.023 113.253 95.4396 Q110.198 90.8332 110.198 82.1064 Q110.198 73.3564 113.253 68.7731 Q116.332 64.1667 122.142 64.1667 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M143.791 46.3627 Q140.857 46.3627 139.371 49.2591 Q137.904 52.1367 137.904 57.9295 Q137.904 63.7035 139.371 66.5999 Q140.857 69.4774 143.791 69.4774 Q146.744 69.4774 148.211 66.5999 Q149.697 63.7035 149.697 57.9295 Q149.697 52.1367 148.211 49.2591 Q146.744 46.3627 143.791 46.3627 M143.791 43.3534 Q148.512 43.3534 150.994 47.0962 Q153.496 50.8201 153.496 57.9295 Q153.496 65.02 150.994 68.7627 Q148.512 72.4867 143.791 72.4867 Q139.07 72.4867 136.569 68.7627 Q134.086 65.02 134.086 57.9295 Q134.086 50.8201 136.569 47.0962 Q139.07 43.3534 143.791 43.3534 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><polyline clip-path="url(#clip882)" style="stroke:#009af9; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
250.72,87.9763 760.923,363.232 1271.13,691.535 1781.33,961.205 2291.53,1445.72
"/>
<line clip-path="url(#clip882)" x1="250.72" y1="87.9763" x2="234.72" y2="71.9763" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="250.72" y1="87.9763" x2="234.72" y2="103.976" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="250.72" y1="87.9763" x2="266.72" y2="103.976" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="250.72" y1="87.9763" x2="266.72" y2="71.9763" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="760.923" y1="363.232" x2="744.923" y2="347.232" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="760.923" y1="363.232" x2="744.923" y2="379.232" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="760.923" y1="363.232" x2="776.923" y2="379.232" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="760.923" y1="363.232" x2="776.923" y2="347.232" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="1271.13" y1="691.535" x2="1255.13" y2="675.535" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="1271.13" y1="691.535" x2="1255.13" y2="707.535" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="1271.13" y1="691.535" x2="1287.13" y2="707.535" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="1271.13" y1="691.535" x2="1287.13" y2="675.535" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="1781.33" y1="961.205" x2="1765.33" y2="945.205" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="1781.33" y1="961.205" x2="1765.33" y2="977.205" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="1781.33" y1="961.205" x2="1797.33" y2="977.205" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="1781.33" y1="961.205" x2="1797.33" y2="945.205" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="2291.53" y1="1445.72" x2="2275.53" y2="1429.72" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="2291.53" y1="1445.72" x2="2275.53" y2="1461.72" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="2291.53" y1="1445.72" x2="2307.53" y2="1461.72" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip882)" x1="2291.53" y1="1445.72" x2="2307.53" y2="1429.72" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<path clip-path="url(#clip880)" d="
M1792.09 216.178 L2280.65 216.178 L2280.65 95.2176 L1792.09 95.2176 Z
" fill="#ffffff" fill-rule="evenodd" fill-opacity="1"/>
<polyline clip-path="url(#clip880)" style="stroke:#000000; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
1792.09,216.178 2280.65,216.178 2280.65,95.2176 1792.09,95.2176 1792.09,216.178
"/>
<polyline clip-path="url(#clip880)" style="stroke:#009af9; stroke-linecap:butt; stroke-linejoin:round; stroke-width:4; stroke-opacity:1; fill:none" points="
1816.13,155.698 1960.35,155.698
"/>
<line clip-path="url(#clip880)" x1="1888.24" y1="155.698" x2="1862.64" y2="130.098" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip880)" x1="1888.24" y1="155.698" x2="1862.64" y2="181.298" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip880)" x1="1888.24" y1="155.698" x2="1913.84" y2="181.298" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<line clip-path="url(#clip880)" x1="1888.24" y1="155.698" x2="1913.84" y2="130.098" style="stroke:#009af9; stroke-width:3.2; stroke-opacity:1"/>
<path clip-path="url(#clip880)" d="M1988.32 136.751 L1988.32 184.158 L1984.38 184.158 L1984.38 136.751 L1988.32 136.751 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M2000.72 151.681 Q2002.14 149.32 2005.63 147.237 Q2007 146.427 2011.21 146.427 Q2015.93 146.427 2018.87 150.177 Q2021.84 153.927 2021.84 160.038 Q2021.84 166.149 2018.87 169.899 Q2015.93 173.649 2011.21 173.649 Q2008.36 173.649 2006.3 172.538 Q2004.27 171.403 2002.92 169.089 L2002.92 182.839 L1998.64 182.839 L1998.64 160.269 Q1998.64 154.922 2000.72 151.681 M2017.41 160.038 Q2017.41 155.339 2015.47 152.677 Q2013.55 149.992 2010.17 149.992 Q2006.79 149.992 2004.84 152.677 Q2002.92 155.339 2002.92 160.038 Q2002.92 164.737 2004.84 167.422 Q2006.79 170.084 2010.17 170.084 Q2013.55 170.084 2015.47 167.422 Q2017.41 164.737 2017.41 160.038 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M2038.94 150.038 Q2035.52 150.038 2033.52 152.723 Q2031.53 155.385 2031.53 160.038 Q2031.53 164.691 2033.5 167.376 Q2035.49 170.038 2038.94 170.038 Q2042.34 170.038 2044.34 167.353 Q2046.33 164.667 2046.33 160.038 Q2046.33 155.431 2044.34 152.746 Q2042.34 150.038 2038.94 150.038 M2038.94 146.427 Q2044.5 146.427 2047.67 150.038 Q2050.84 153.649 2050.84 160.038 Q2050.84 166.404 2047.67 170.038 Q2044.5 173.649 2038.94 173.649 Q2033.36 173.649 2030.19 170.038 Q2027.04 166.404 2027.04 160.038 Q2027.04 153.649 2030.19 150.038 Q2033.36 146.427 2038.94 146.427 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M2057.46 162.746 L2057.46 147.052 L2061.72 147.052 L2061.72 162.584 Q2061.72 166.265 2063.15 168.116 Q2064.59 169.945 2067.46 169.945 Q2070.91 169.945 2072.9 167.746 Q2074.91 165.547 2074.91 161.751 L2074.91 147.052 L2079.17 147.052 L2079.17 172.978 L2074.91 172.978 L2074.91 168.996 Q2073.36 171.357 2071.3 172.515 Q2069.27 173.649 2066.56 173.649 Q2062.09 173.649 2059.77 170.871 Q2057.46 168.093 2057.46 162.746 M2068.18 146.427 L2068.18 146.427 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M2092.16 139.691 L2092.16 147.052 L2100.93 147.052 L2100.93 150.362 L2092.16 150.362 L2092.16 164.436 Q2092.16 167.607 2093.02 168.51 Q2093.89 169.413 2096.56 169.413 L2100.93 169.413 L2100.93 172.978 L2096.56 172.978 Q2091.63 172.978 2089.75 171.149 Q2087.88 169.297 2087.88 164.436 L2087.88 150.362 L2084.75 150.362 L2084.75 147.052 L2087.88 147.052 L2087.88 139.691 L2092.16 139.691 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M2119.45 158.093 L2131.93 158.093 L2131.93 161.89 L2119.45 161.89 L2119.45 158.093 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M2155.7 151.681 Q2157.11 149.32 2160.61 147.237 Q2161.97 146.427 2166.19 146.427 Q2170.91 146.427 2173.85 150.177 Q2176.81 153.927 2176.81 160.038 Q2176.81 166.149 2173.85 169.899 Q2170.91 173.649 2166.19 173.649 Q2163.34 173.649 2161.28 172.538 Q2159.24 171.403 2157.9 169.089 L2157.9 182.839 L2153.62 182.839 L2153.62 160.269 Q2153.62 154.922 2155.7 151.681 M2172.39 160.038 Q2172.39 155.339 2170.45 152.677 Q2168.52 149.992 2165.14 149.992 Q2161.76 149.992 2159.82 152.677 Q2157.9 155.339 2157.9 160.038 Q2157.9 164.737 2159.82 167.422 Q2161.76 170.084 2165.14 170.084 Q2168.52 170.084 2170.45 167.422 Q2172.39 164.737 2172.39 160.038 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M2183.87 147.052 L2188.13 147.052 L2188.13 172.978 L2183.87 172.978 L2183.87 147.052 M2183.87 136.959 L2188.13 136.959 L2188.13 142.353 L2183.87 142.353 L2183.87 136.959 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M2218.59 157.329 L2218.59 172.978 L2214.33 172.978 L2214.33 157.468 Q2214.33 153.788 2212.9 151.959 Q2211.46 150.13 2208.59 150.13 Q2205.14 150.13 2203.15 152.33 Q2201.16 154.529 2201.16 158.325 L2201.16 172.978 L2196.88 172.978 L2196.88 147.052 L2201.16 147.052 L2201.16 151.08 Q2202.69 148.742 2204.75 147.584 Q2206.83 146.427 2209.54 146.427 Q2214.01 146.427 2216.3 149.205 Q2218.59 151.959 2218.59 157.329 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /><path clip-path="url(#clip880)" d="M2232.57 136.751 L2232.57 184.158 L2228.64 184.158 L2228.64 136.751 L2232.57 136.751 Z" fill="#000000" fill-rule="evenodd" fill-opacity="1" /></svg>
<p>The <code>info</code> object passed to the callback contains not just the densities but also the complete Bloch wave (in <code>ψ</code>), the <code>occupation</code>, band <code>eigenvalues</code> and so on. See <a href="https://dftk.org/blob/master/src/scf/self_consistent_field.jl#L101"><code>src/scf/self_consistent_field.jl</code></a> for all currently available keys.</p><div class="admonition is-success"><header class="admonition-header">Debugging with callbacks</header><div class="admonition-body"><p>Very handy for debugging SCF algorithms is to employ callbacks with an <code>@infiltrate</code> from <a href="https://github.com/JuliaDebug/Infiltrator.jl">Infiltrator.jl</a> to interactively monitor what is happening each SCF step.</p></div></div></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="../geometry_optimization/">« Geometry optimization</a><a class="docs-footer-nextpage" href="../scf_checkpoints/">Saving SCF results on disk and SCF checkpoints »</a><div class="flexbox-break"></div><p class="footer-message">Powered by <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> and the <a href="https://julialang.org/">Julia Programming Language</a>.</p></nav></div><div class="modal" id="documenter-settings"><div class="modal-background"></div><div class="modal-card"><header class="modal-card-head"><p class="modal-card-title">Settings</p><button class="delete"></button></header><section class="modal-card-body"><p><label class="label">Theme</label><div class="select"><select id="documenter-themepicker"><option value="documenter-light">documenter-light</option><option value="documenter-dark">documenter-dark</option></select></div></p><hr/><p>This document was generated with <a href="https://github.com/JuliaDocs/Documenter.jl">Documenter.jl</a> on <span class="colophon-date" title="Thursday 3 June 2021 17:00">Thursday 3 June 2021</span>. Using Julia version 1.6.1.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>