Comparison of DFT solvers
We compare four different approaches for solving the DFT minimisation problem, namely a density-based SCF, a potential-based SCF, direct minimisation and Newton.
First we setup our problem
using DFTK
using LinearAlgebra
a = 10.26 # Silicon lattice constant in Bohr
lattice = a / 2 * [[0 1 1.];
[1 0 1.];
[1 1 0.]]
Si = ElementPsp(:Si; psp=load_psp("hgh/lda/Si-q4"))
atoms = [Si, Si]
positions = [ones(3)/8, -ones(3)/8]
model = model_LDA(lattice, atoms, positions)
basis = PlaneWaveBasis(model; Ecut=5, kgrid=[3, 3, 3])
# Convergence we desire in the density
tol = 1e-61.0e-6Density-based self-consistent field
scfres_scf = self_consistent_field(basis; tol);n Energy log10(ΔE) log10(Δρ) Diag Δtime
--- --------------- --------- --------- ---- ------
1 -7.847793415581 -0.70 4.8 39.1ms
2 -7.852411199040 -2.34 -1.53 1.0 17.1ms
3 -7.852610879115 -3.70 -2.54 1.0 16.8ms
4 -7.852645449757 -4.46 -2.77 2.2 21.5ms
5 -7.852646066618 -6.21 -2.87 1.0 16.9ms
6 -7.852646668063 -6.22 -3.84 1.0 16.9ms
7 -7.852646684876 -7.77 -4.54 1.2 17.9ms
8 -7.852646686696 -8.74 -5.09 1.8 20.0ms
9 -7.852646686727 -10.51 -5.83 1.2 18.3ms
10 -7.852646686728 -11.78 -5.66 2.0 20.8ms
11 -7.852646686730 -11.87 -6.15 1.0 17.7msPotential-based SCF
scfres_scfv = DFTK.scf_potential_mixing(basis; tol);n Energy log10(ΔE) log10(Δρ) α Diag Δtime
--- --------------- --------- --------- ---- ---- ------
1 -7.847731724603 -0.70 4.5 29.7ms
2 -7.852560384301 -2.32 -1.62 0.80 2.0 20.2ms
3 -7.852640324683 -4.10 -2.71 0.80 1.0 16.0ms
4 -7.852646562329 -5.20 -3.35 0.80 2.0 19.8ms
5 -7.852646680911 -6.93 -4.39 0.80 1.2 17.1ms
6 -7.852646686602 -8.24 -4.84 0.80 2.0 20.6ms
7 -7.852646686723 -9.92 -5.48 0.80 1.2 16.9ms
8 -7.852646686730 -11.17 -6.55 0.80 1.8 18.5msDirect minimization
Note: Unlike the other algorithms, tolerance for this one is in the energy, thus we square the density tolerance value to be roughly equivalent.
scfres_dm = direct_minimization(basis; tol=tol^2);n Energy log10(ΔE) log10(Δρ) Δtime
--- --------------- --------- --------- ------
1 +1.352591858994 -1.01 44.5ms
2 -1.293552376780 0.42 -0.63 24.6ms
3 -3.620594405474 0.37 -0.44 32.5ms
4 -4.667042118362 0.02 -0.56 41.8ms
5 -6.477287009564 0.26 -0.68 33.8ms
6 -7.228194857748 -0.12 -1.10 33.1ms
7 -7.468689403252 -0.62 -1.48 25.1ms
8 -7.648392209908 -0.75 -1.96 24.9ms
9 -7.692041625438 -1.36 -2.19 24.9ms
10 -7.725418022551 -1.48 -2.28 24.9ms
11 -7.749622607035 -1.62 -2.22 24.9ms
12 -7.781044032048 -1.50 -2.13 25.0ms
13 -7.815746693196 -1.46 -1.90 25.2ms
14 -7.832841364537 -1.77 -2.36 25.0ms
15 -7.845681530423 -1.89 -2.38 25.0ms
16 -7.846110910922 -3.37 -2.86 33.1ms
17 -7.846189294969 -4.11 -2.36 49.3ms
18 -7.846264180461 -4.13 -2.27 41.7ms
19 -7.846380213303 -3.94 -2.51 42.2ms
20 -7.847050852519 -3.17 -2.80 25.3ms
21 -7.849588583978 -2.60 -3.46 25.1ms
22 -7.852512311485 -2.53 -3.39 25.0ms
23 -7.852580897595 -4.16 -3.85 24.9ms
24 -7.852622859810 -4.38 -4.09 25.0ms
25 -7.852638726549 -4.80 -4.34 24.8ms
26 -7.852642982334 -5.37 -4.48 24.9ms
27 -7.852644939676 -5.71 -5.21 24.9ms
28 -7.852646174423 -5.91 -4.71 24.9ms
29 -7.852646492833 -6.50 -4.91 24.8ms
30 -7.852646597775 -6.98 -5.47 24.7ms
31 -7.852646653981 -7.25 -6.07 24.8ms
32 -7.852646676008 -7.66 -5.75 24.8ms
33 -7.852646684324 -8.08 -5.86 24.8ms
34 -7.852646686102 -8.75 -6.54 24.8ms
35 -7.852646686519 -9.38 -6.69 24.8ms
36 -7.852646686661 -9.85 -6.66 24.8ms
37 -7.852646686699 -10.41 -7.25 24.8ms
38 -7.852646686717 -10.74 -7.32 24.8ms
39 -7.852646686726 -11.05 -7.37 24.8ms
40 -7.852646686729 -11.60 -7.56 24.7ms
41 -7.852646686730 -12.16 -8.03 24.7ms
42 -7.852646686730 -12.69 -8.25 24.8ms
43 -7.852646686730 -13.18 -8.22 24.8ms
44 -7.852646686730 -13.60 -8.56 24.9ms
45 -7.852646686730 -14.15 -9.04 24.7ms
46 -7.852646686730 + -Inf -9.52 24.7ms
47 -7.852646686730 -14.57 -9.46 24.8ms
48 -7.852646686730 -15.05 -9.68 49.0ms
49 -7.852646686730 + -Inf -9.65 41.1ms
50 -7.852646686730 + -Inf -9.73 25.0ms
┌ Warning: DM not converged.
└ @ DFTK ~/work/DFTK.jl/DFTK.jl/src/scf/scf_callbacks.jl:60Newton algorithm
Start not too far from the solution to ensure convergence: We run first a very crude SCF to get close and then switch to Newton.
scfres_start = self_consistent_field(basis; tol=0.5);n Energy log10(ΔE) log10(Δρ) Diag Δtime
--- --------------- --------- --------- ---- ------
1 -7.847760454769 -0.70 4.8 30.3msRemove the virtual orbitals (which Newton cannot treat yet)
ψ = DFTK.select_occupied_orbitals(basis, scfres_start.ψ, scfres_start.occupation).ψ
scfres_newton = newton(basis, ψ; tol);n Energy log10(ΔE) log10(Δρ) Δtime
--- --------------- --------- --------- ------
1 -7.852645844855 -1.63 390ms
2 -7.852646686730 -6.07 -3.69 274ms
3 -7.852646686730 -13.22 -7.20 137msComparison of results
println("|ρ_newton - ρ_scf| = ", norm(scfres_newton.ρ - scfres_scf.ρ))
println("|ρ_newton - ρ_scfv| = ", norm(scfres_newton.ρ - scfres_scfv.ρ))
println("|ρ_newton - ρ_dm| = ", norm(scfres_newton.ρ - scfres_dm.ρ))|ρ_newton - ρ_scf| = 7.251449160836048e-7
|ρ_newton - ρ_scfv| = 2.6679838064898556e-7
|ρ_newton - ρ_dm| = 6.165357765969544e-10