/
index.html
68 lines (65 loc) · 11.1 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
<!DOCTYPE html>
<html lang="en"><head><meta charset="UTF-8"/><meta name="viewport" content="width=device-width, initial-scale=1.0"/><title>Custom solvers · DFTK.jl</title><link rel="canonical" href="https://juliamolsim.github.io/DFTK.jl/stable/examples/custom_solvers/"/><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"/><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/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><a class="tocitem" href="../scf_callbacks/">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 class="is-active"><a class="tocitem" href>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>Custom solvers</a></li></ul><ul class="is-hidden-tablet"><li class="is-active"><a href>Custom solvers</a></li></ul></nav><div class="docs-right"><a class="docs-edit-link" href="https://github.com/JuliaMolSim/DFTK.jl/blob/master/examples/custom_solvers.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="Custom-solvers"><a class="docs-heading-anchor" href="#Custom-solvers">Custom solvers</a><a id="Custom-solvers-1"></a><a class="docs-heading-anchor-permalink" href="#Custom-solvers" title="Permalink"></a></h1><p><a href="https://mybinder.org/v2/gh/JuliaMolSim/DFTK.jl/gh-pages?filepath=v0.2.4/examples/custom_solvers.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.2.4/examples/custom_solvers.ipynb"><img src="https://img.shields.io/badge/show-nbviewer-579ACA.svg" alt/></a></p><p>In this example, we show how to define custom solvers. Our system will again be silicon, because we are not very imaginative</p><pre><code class="language-julia">using DFTK, LinearAlgebra
a = 10.26
lattice = a / 2 * [[0 1 1.];
[1 0 1.];
[1 1 0.]]
Si = ElementPsp(:Si, psp=load_psp("hgh/lda/Si-q4"))
atoms = [Si => [ones(3)/8, -ones(3)/8]]
# We take very (very) crude parameters
model = model_LDA(lattice, atoms)
kgrid = [1, 1, 1]
Ecut = 5
basis = PlaneWaveBasis(model, Ecut; kgrid=kgrid);</code></pre><p>We define our custom fix-point solver: simply a damped fixed-point</p><pre><code class="language-julia">function my_fp_solver(f, x0, max_iter; tol)
mixing_factor = .7
x = x0
fx = f(x)
for n = 1:max_iter
inc = fx - x
if norm(inc) < tol
break
end
x = x + mixing_factor * inc
fx = f(x)
end
(fixpoint=x, converged=norm(fx-x) < tol)
end;</code></pre><p>Our eigenvalue solver just forms the dense matrix and diagonalizes it explicitly (this only works for very small systems)</p><pre><code class="language-julia">function my_eig_solver(A, X0; maxiter, tol, kwargs...)
n = size(X0, 2)
A = Array(A)
E = eigen(A)
λ = E.values[1:n]
X = E.vectors[:, 1:n]
(λ=λ, X=X, residual_norms=[], iterations=0, converged=true, n_matvec=0)
end;</code></pre><p>Finally we also define our custom mixing scheme. It will be a mixture of simple mixing (for the first 2 steps) and than default to Kerker mixing. In the mixing interface <code>δF</code> is <span>$(ρ_\text{out} - ρ_\text{in})$</span>, i.e. the difference in density between two subsequent SCF steps and the <code>mix</code> function returns <span>$δρ$</span>, which is added to <span>$ρ_\text{in}$</span> to yield <span>$ρ_\text{next}$</span>, the density for the next SCF step.</p><pre><code class="language-julia">struct MyMixing
α # Damping parameter
end
MyMixing() = MyMixing(0.7)
function DFTK.mix(mixing::MyMixing, basis, δF::RealFourierArray, δF_spin=nothing; n_iter, kwargs...)
if n_iter <= 2
# Just do simple mixing on total density and spin density (if it exists)
(mixing.α * δF, isnothing(δF_spin) ? nothing : mixing.α * δF_spin)
else
# Use the KerkerMixing from DFTK
DFTK.mix(KerkerMixing(α=mixing.α), basis, δF, δF_spin; kwargs...)
end
end</code></pre><p>That's it! Now we just run the SCF with these solvers</p><pre><code class="language-julia">scfres = self_consistent_field(basis;
tol=1e-8,
solver=my_fp_solver,
eigensolver=my_eig_solver,
mixing=MyMixing());</code></pre><pre class="documenter-example-output">n Energy Eₙ-Eₙ₋₁ ρout-ρin Diag
--- --------------- --------- -------- ----
1 -7.234302105583 NaN 3.17e-01 0.0
2 -7.246539783547 -1.22e-02 1.44e-01 0.0
3 -7.248795495550 -2.26e-03 6.55e-02 0.0
4 -7.249077010080 -2.82e-04 3.61e-02 0.0
5 -7.249170222206 -9.32e-05 2.00e-02 0.0
6 -7.249201046259 -3.08e-05 1.12e-02 0.0
7 -7.249211378481 -1.03e-05 6.34e-03 0.0
8 -7.249214920041 -3.54e-06 3.64e-03 0.0
9 -7.249216167364 -1.25e-06 2.12e-03 0.0
10 -7.249216619477 -4.52e-07 1.25e-03 0.0
11 -7.249216787981 -1.69e-07 7.50e-04 0.0
12 -7.249216852389 -6.44e-08 4.56e-04 0.0
13 -7.249216877549 -2.52e-08 2.81e-04 0.0
14 -7.249216887557 -1.00e-08 1.75e-04 0.0
15 -7.249216891595 -4.04e-09 1.10e-04 0.0</pre><p>Note that the default convergence criterion is on the difference of energy from one step to the other; when this gets below <code>tol</code>, the "driver" <code>self_consistent_field</code> artificially makes the fixpoint solver think it's converged by forcing <code>f(x) = x</code>. You can customize this with the <code>is_converged</code> keyword argument to <code>self_consistent_field</code>.</p></article><nav class="docs-footer"><a class="docs-footer-prevpage" href="../arbitrary_floattype/">« Arbitrary floating-point types</a><a class="docs-footer-nextpage" href="../custom_potential/">Custom potential »</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="Tuesday 29 December 2020 12:37">Tuesday 29 December 2020</span>. Using Julia version 1.5.3.</p></section><footer class="modal-card-foot"></footer></div></div></div></body></html>