Eigenvalues of the dielectric matrix

We compute a few eigenvalues of the dielectric matrix ($q=0$, $ω=0$) iteratively.

In [1]:
using DFTK
using Plots
using KrylovKit
using Printf

# Calculation parameters
kgrid = [1, 1, 1]
Ecut = 5

# Silicon lattice
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, Si]
positions = [ones(3)/8, -ones(3)/8]

# Compute the dielectric operator without symmetries
model = model_LDA(lattice, atoms, positions, symmetries=false)
basis = PlaneWaveBasis(model; Ecut, kgrid)
scfres = self_consistent_field(basis, tol=1e-14);
n     Energy            log10(ΔE)   log10(Δρ)   Diag
---   ---------------   ---------   ---------   ----
  1   -7.235314313664                   -0.50    7.0
  2   -7.250106077970       -1.83       -1.39    1.0
  3   -7.250778594233       -3.17       -1.75    4.0
  4   -7.250915885485       -3.86       -1.85    2.0
  5   -7.251331684373       -3.38       -2.69    2.0
  6   -7.251338277330       -5.18       -3.21    2.0
  7   -7.251338779684       -6.30       -3.88    2.0
  8   -7.251338794830       -7.82       -4.17    3.0
  9   -7.251338798525       -8.43       -4.73    2.0
 10   -7.251338798680       -9.81       -5.32    2.0
 11   -7.251338798701      -10.66       -5.82    3.0
 12   -7.251338798704      -11.53       -6.55    2.0
 13   -7.251338798705      -13.12       -6.95    3.0
 14   -7.251338798705      -13.97       -7.42    2.0
 15   -7.251338798705      -14.05       -8.10    2.0

Applying $ε^† ≔ (1- χ_0 K)$ …

In [2]:
function eps_fun(δρ)
    δV = apply_kernel(basis, δρ; ρ=scfres.ρ)
    χ0δV = apply_χ0(scfres, δV)
    δρ - χ0δV
end;

… eagerly diagonalizes the subspace matrix at each iteration

In [3]:
eigsolve(eps_fun, randn(size(scfres.ρ)), 5, :LM; eager=true, verbosity=3);
[ Info: Arnoldi iteration step 1: normres = 0.08851872482341214
[ Info: Arnoldi iteration step 2: normres = 0.4524401464036587
[ Info: Arnoldi iteration step 3: normres = 0.9171121180848835
[ Info: Arnoldi iteration step 4: normres = 0.23502035278006775
[ Info: Arnoldi iteration step 5: normres = 0.6175782816699701
[ Info: Arnoldi schursolve in iter 1, krylovdim = 5: 0 values converged, normres = (3.70e-02, 6.93e-02, 5.74e-01, 2.14e-01, 2.46e-03)
[ Info: Arnoldi iteration step 6: normres = 0.23542500804335315
[ Info: Arnoldi schursolve in iter 1, krylovdim = 6: 0 values converged, normres = (7.23e-03, 1.56e-01, 1.04e-01, 1.01e-01, 5.02e-02)
[ Info: Arnoldi iteration step 7: normres = 0.09183570814984225
[ Info: Arnoldi schursolve in iter 1, krylovdim = 7: 0 values converged, normres = (3.23e-04, 1.31e-02, 6.83e-03, 6.51e-02, 5.53e-02)
[ Info: Arnoldi iteration step 8: normres = 0.12839471373667652
[ Info: Arnoldi schursolve in iter 1, krylovdim = 8: 0 values converged, normres = (1.79e-05, 1.19e-03, 6.97e-04, 2.65e-02, 4.09e-02)
[ Info: Arnoldi iteration step 9: normres = 0.07768277815756745
[ Info: Arnoldi schursolve in iter 1, krylovdim = 9: 0 values converged, normres = (6.10e-07, 6.74e-05, 4.41e-05, 8.10e-03, 3.53e-02)
[ Info: Arnoldi iteration step 10: normres = 0.0816457446242494
[ Info: Arnoldi schursolve in iter 1, krylovdim = 10: 0 values converged, normres = (2.19e-08, 4.02e-06, 2.93e-06, 2.54e-03, 2.80e-02)
[ Info: Arnoldi iteration step 11: normres = 0.07240215304447777
[ Info: Arnoldi schursolve in iter 1, krylovdim = 11: 0 values converged, normres = (6.76e-10, 2.03e-07, 1.64e-07, 5.62e-04, 1.21e-02)
[ Info: Arnoldi iteration step 12: normres = 0.07331149865050755
[ Info: Arnoldi schursolve in iter 1, krylovdim = 12: 0 values converged, normres = (2.11e-11, 1.03e-08, 9.24e-09, 1.20e-04, 4.74e-03)
[ Info: Arnoldi iteration step 13: normres = 0.10445245277205002
[ Info: Arnoldi schursolve in iter 1, krylovdim = 13: 1 values converged, normres = (9.49e-13, 7.63e-10, 7.57e-10, 4.00e-05, 3.27e-03)
[ Info: Arnoldi iteration step 14: normres = 0.3673600834842453
[ Info: Arnoldi schursolve in iter 1, krylovdim = 14: 1 values converged, normres = (3.01e-13, 1.12e-09, 2.54e-09, 3.66e-01, 9.86e-04)
[ Info: Arnoldi iteration step 15: normres = 0.17036936933576477
[ Info: Arnoldi schursolve in iter 1, krylovdim = 15: 1 values converged, normres = (2.38e-14, 2.60e-10, 4.74e-02, 5.65e-04, 1.03e-05)
[ Info: Arnoldi iteration step 16: normres = 0.19592603284920262
[ Info: Arnoldi schursolve in iter 1, krylovdim = 16: 1 values converged, normres = (4.89e-15, 1.76e-09, 1.80e-01, 8.01e-03, 7.44e-02)
[ Info: Arnoldi iteration step 17: normres = 0.028900213563365175
[ Info: Arnoldi schursolve in iter 1, krylovdim = 17: 1 values converged, normres = (6.01e-17, 1.46e-03, 3.20e-03, 1.27e-06, 1.49e-03)
[ Info: Arnoldi iteration step 18: normres = 0.019985545474946024
[ Info: Arnoldi schursolve in iter 1, krylovdim = 18: 1 values converged, normres = (4.97e-19, 1.21e-05, 4.48e-05, 3.31e-08, 2.16e-05)
[ Info: Arnoldi iteration step 19: normres = 0.19165269645955182
[ Info: Arnoldi schursolve in iter 1, krylovdim = 19: 1 values converged, normres = (4.12e-20, 3.32e-08, 6.33e-06, 1.63e-06, 2.84e-06)
[ Info: Arnoldi iteration step 20: normres = 0.0402260087202923
[ Info: Arnoldi schursolve in iter 1, krylovdim = 20: 1 values converged, normres = (7.84e-22, 2.09e-07, 4.21e-08, 7.86e-08, 9.76e-08)
[ Info: Arnoldi iteration step 21: normres = 0.03569143430508035
[ Info: Arnoldi schursolve in iter 1, krylovdim = 21: 1 values converged, normres = (1.16e-23, 5.06e-09, 2.93e-10, 2.35e-09, 2.29e-09)
[ Info: Arnoldi iteration step 22: normres = 0.02537879670053612
[ Info: Arnoldi schursolve in iter 1, krylovdim = 22: 1 values converged, normres = (1.22e-25, 8.31e-11, 1.90e-11, 4.47e-11, 4.13e-11)
[ Info: Arnoldi iteration step 23: normres = 0.545405848584182
[ Info: Arnoldi schursolve in iter 1, krylovdim = 23: 1 values converged, normres = (3.12e-26, 3.70e-11, 8.32e-12, 2.25e-11, 2.07e-11)
[ Info: Arnoldi iteration step 24: normres = 0.0298484169173863
[ Info: Arnoldi schursolve in iter 1, krylovdim = 24: 1 values converged, normres = (7.80e-28, 6.30e-12, 1.42e-12, 7.30e-10, 6.99e-10)
[ Info: Arnoldi iteration step 25: normres = 0.053924432323077184
[ Info: Arnoldi schursolve in iter 1, krylovdim = 25: 3 values converged, normres = (1.73e-29, 2.24e-13, 5.06e-14, 9.52e-10, 8.83e-10)
[ Info: Arnoldi iteration step 26: normres = 0.06069937248059313
[ Info: Arnoldi schursolve in iter 1, krylovdim = 26: 3 values converged, normres = (4.73e-31, 1.04e-14, 2.35e-15, 4.67e-05, 3.91e-06)
[ Info: Arnoldi iteration step 27: normres = 0.034266516183009264
[ Info: Arnoldi schursolve in iter 1, krylovdim = 27: 3 values converged, normres = (6.70e-33, 2.35e-16, 5.31e-17, 2.09e-09, 1.35e-09)
[ Info: Arnoldi iteration step 28: normres = 0.08441934556941956
[ Info: Arnoldi schursolve in iter 1, krylovdim = 28: 3 values converged, normres = (2.48e-34, 1.44e-17, 3.26e-18, 8.09e-08, 2.26e-08)
[ Info: Arnoldi iteration step 29: normres = 0.01760559206301117
[ Info: Arnoldi schursolve in iter 1, krylovdim = 29: 3 values converged, normres = (1.82e-36, 1.71e-19, 3.86e-20, 3.40e-11, 1.14e-09)
[ Info: Arnoldi iteration step 30: normres = 0.16863310639846676
[ Info: Arnoldi schursolve in iter 1, krylovdim = 30: 3 values converged, normres = (1.30e-37, 1.98e-20, 4.48e-21, 4.31e-12, 1.49e-10)
[ Info: Arnoldi schursolve in iter 2, krylovdim = 19: 3 values converged, normres = (1.30e-37, 1.98e-20, 4.48e-21, 4.31e-12, 1.49e-10)
[ Info: Arnoldi iteration step 20: normres = 0.09527878615188579
[ Info: Arnoldi schursolve in iter 2, krylovdim = 20: 4 values converged, normres = (5.97e-39, 1.62e-21, 3.65e-22, 4.00e-13, 1.38e-11)
[ Info: Arnoldi iteration step 21: normres = 0.06726132883967972
┌ Info: Arnoldi eigsolve finished after 2 iterations:
│ *  6 eigenvalues converged
│ *  norm of residuals = (1.741272246909043e-40, 7.815673750792932e-23, 3.512897732544951e-23, 2.1513527375105488e-14, 7.377435491512547e-13, 7.594400714654044e-14)
└ *  number of operations = 32