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-6
1.0e-6
Density-based self-consistent field
scfres_scf = self_consistent_field(basis; tol);
n Energy log10(ΔE) log10(Δρ) Diag Δtime
--- --------------- --------- --------- ---- ------
1 -7.847754733182 -0.70 4.5 29.1ms
2 -7.852400871343 -2.33 -1.53 1.0 17.3ms
3 -7.852609952983 -3.68 -2.54 1.2 17.7ms
4 -7.852645443883 -4.45 -2.77 2.2 21.6ms
5 -7.852646077301 -6.20 -2.87 1.2 18.0ms
6 -7.852646666286 -6.23 -3.83 1.0 17.1ms
7 -7.852646684960 -7.73 -4.65 1.5 19.1ms
8 -7.852646686694 -8.76 -5.06 1.8 19.9ms
9 -7.852646686726 -10.50 -5.72 1.0 17.8ms
10 -7.852646686728 -11.59 -5.66 1.5 20.1ms
11 -7.852646686730 -11.87 -6.24 1.0 34.6ms
Potential-based SCF
scfres_scfv = DFTK.scf_potential_mixing(basis; tol);
n Energy log10(ΔE) log10(Δρ) α Diag Δtime
--- --------------- --------- --------- ---- ---- ------
1 -7.847780929274 -0.70 4.8 31.1ms
2 -7.852560337567 -2.32 -1.62 0.80 2.0 20.5ms
3 -7.852641769653 -4.09 -2.68 0.80 1.0 16.4ms
4 -7.852646527308 -5.32 -3.41 0.80 1.8 19.5ms
5 -7.852646681635 -6.81 -4.40 0.80 1.5 18.6ms
6 -7.852646686624 -8.30 -4.84 0.80 2.5 22.5ms
7 -7.852646686726 -9.99 -5.52 0.80 1.2 17.3ms
8 -7.852646686730 -11.47 -6.56 0.80 1.5 19.5ms
Direct 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.362678295933 -0.97 45.8ms
2 -1.203990905250 0.41 -0.62 25.5ms
3 -3.543812684930 0.37 -0.43 33.7ms
4 -4.836743009589 0.11 -0.57 33.6ms
5 -6.687730082863 0.27 -0.67 33.2ms
6 -6.758871660422 -1.15 -1.57 25.2ms
7 -7.460227505159 -0.15 -1.41 25.0ms
8 -7.552217784959 -1.04 -1.85 25.0ms
9 -7.710400228843 -0.80 -1.85 25.0ms
10 -7.751415011380 -1.39 -2.18 24.9ms
11 -7.778758673453 -1.56 -2.35 24.8ms
12 -7.802674537301 -1.62 -2.36 24.7ms
13 -7.830896915392 -1.55 -2.18 24.7ms
14 -7.843221345422 -1.91 -2.53 24.8ms
15 -7.848336473864 -2.29 -2.26 24.9ms
16 -7.851691769535 -2.47 -2.58 24.8ms
17 -7.852282190995 -3.23 -3.37 24.8ms
18 -7.852450603248 -3.77 -3.80 24.8ms
19 -7.852568707405 -3.93 -3.38 24.8ms
20 -7.852622392684 -4.27 -3.63 25.2ms
21 -7.852635062012 -4.90 -3.66 25.3ms
22 -7.852641809902 -5.17 -4.13 25.2ms
23 -7.852645429596 -5.44 -4.58 25.4ms
24 -7.852646348164 -6.04 -4.68 25.3ms
25 -7.852646585798 -6.62 -4.84 25.0ms
26 -7.852646657574 -7.14 -5.31 25.1ms
27 -7.852646677416 -7.70 -5.27 24.9ms
28 -7.852646683442 -8.22 -5.94 24.8ms
29 -7.852646685523 -8.68 -6.06 25.1ms
30 -7.852646686263 -9.13 -6.18 25.1ms
31 -7.852646686590 -9.49 -6.54 25.0ms
32 -7.852646686686 -10.02 -6.60 24.8ms
33 -7.852646686716 -10.52 -6.91 24.9ms
34 -7.852646686726 -11.02 -7.29 25.3ms
35 -7.852646686728 -11.56 -7.46 25.1ms
36 -7.852646686729 -12.03 -7.65 25.1ms
37 -7.852646686730 -12.51 -8.02 25.1ms
38 -7.852646686730 -13.05 -7.97 25.0ms
39 -7.852646686730 -13.43 -8.25 24.9ms
40 -7.852646686730 -13.80 -9.02 25.7ms
41 -7.852646686730 -13.91 -8.80 34.3ms
42 -7.852646686730 -15.05 -9.18 30.3ms
43 -7.852646686730 -15.05 -9.37 25.4ms
44 -7.852646686730 -14.75 -9.50 25.2ms
45 -7.852646686730 -15.05 -9.62 25.0ms
46 -7.852646686730 -15.05 -10.20 41.0ms
47 -7.852646686730 + -Inf -10.06 32.9ms
48 -7.852646686730 -14.75 -10.08 65.1ms
49 -7.852646686730 + -Inf -10.03 73.7ms
50 -7.852646686730 + -Inf -10.01 73.4ms
┌ Warning: DM not converged.
└ @ DFTK ~/work/DFTK.jl/DFTK.jl/src/scf/scf_callbacks.jl:60
Newton 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.847764830578 -0.70 4.8 43.1ms
Remove 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.852645831532 -1.63 376ms
2 -7.852646686730 -6.07 -3.68 274ms
3 -7.852646686730 -13.21 -7.19 121ms
Comparison 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| = 5.819742373427392e-7
|ρ_newton - ρ_scfv| = 4.557405065127045e-7
|ρ_newton - ρ_dm| = 7.082803768398802e-10