Eigenvalues of the dielectric matrix

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

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.235436700040                   -0.50    8.0
  2   -7.250415738519       -1.82       -1.40    1.0
  3   -7.251332417431       -3.04       -2.23    4.0
  4   -7.251337475228       -5.30       -2.91    3.0
  5   -7.251338104798       -6.20       -3.15    4.0
  6   -7.251338739529       -6.20       -3.73    3.0
  7   -7.251338775385       -7.45       -3.96    2.0
  8   -7.251338798238       -7.64       -4.85    1.0
  9   -7.251338798588       -9.46       -5.03    3.0
 10   -7.251338798699       -9.95       -5.75    2.0
 11   -7.251338798704      -11.29       -6.40    3.0
 12   -7.251338798704      -12.94       -6.67    3.0
 13   -7.251338798705      -13.43       -7.02    2.0
 14   -7.251338798705      -13.97       -7.68    2.0
 15   -7.251338798705   +    -Inf       -8.29    3.0

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

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

eigsolve(eps_fun, randn(size(scfres.ρ)), 5, :LM; eager=true, verbosity=3);
[ Info: Arnoldi iteration step 1: normres = 0.07535729274171925
[ Info: Arnoldi iteration step 2: normres = 0.38348799463753536
[ Info: Arnoldi iteration step 3: normres = 1.0279336807217834
[ Info: Arnoldi iteration step 4: normres = 0.2903088388762029
[ Info: Arnoldi iteration step 5: normres = 0.6331088426501572
[ Info: Arnoldi schursolve in iter 1, krylovdim = 5: 0 values converged, normres = (6.65e-02, 9.69e-02, 5.89e-01, 2.01e-01, 3.58e-03)
[ Info: Arnoldi iteration step 6: normres = 0.19660171567421703
[ Info: Arnoldi schursolve in iter 1, krylovdim = 6: 0 values converged, normres = (9.52e-03, 9.11e-02, 1.03e-01, 1.03e-01, 5.69e-02)
[ Info: Arnoldi iteration step 7: normres = 0.11050195290378062
[ Info: Arnoldi schursolve in iter 1, krylovdim = 7: 0 values converged, normres = (4.98e-04, 8.62e-03, 1.00e-02, 7.52e-02, 6.99e-02)
[ Info: Arnoldi iteration step 8: normres = 0.10184636776541693
[ Info: Arnoldi schursolve in iter 1, krylovdim = 8: 0 values converged, normres = (2.24e-05, 6.51e-04, 8.43e-04, 3.01e-02, 4.70e-02)
[ Info: Arnoldi iteration step 9: normres = 0.08977356885036525
[ Info: Arnoldi schursolve in iter 1, krylovdim = 9: 0 values converged, normres = (8.62e-07, 4.10e-05, 5.89e-05, 8.64e-03, 2.58e-02)
[ Info: Arnoldi iteration step 10: normres = 0.0733526665909979
[ Info: Arnoldi schursolve in iter 1, krylovdim = 10: 0 values converged, normres = (2.80e-08, 2.23e-06, 3.57e-06, 2.63e-03, 2.11e-02)
[ Info: Arnoldi iteration step 11: normres = 0.06028794117864753
[ Info: Arnoldi schursolve in iter 1, krylovdim = 11: 0 values converged, normres = (7.12e-10, 9.17e-08, 1.63e-07, 4.38e-04, 6.80e-03)
[ Info: Arnoldi iteration step 12: normres = 0.08337114299872637
[ Info: Arnoldi schursolve in iter 1, krylovdim = 12: 0 values converged, normres = (2.55e-11, 5.38e-09, 1.06e-08, 1.12e-04, 3.26e-03)
[ Info: Arnoldi iteration step 13: normres = 0.08356642389521317
[ Info: Arnoldi schursolve in iter 1, krylovdim = 13: 1 values converged, normres = (9.14e-13, 3.17e-10, 6.91e-10, 2.95e-05, 1.80e-03)
[ Info: Arnoldi iteration step 14: normres = 0.4813972029903742
[ Info: Arnoldi schursolve in iter 1, krylovdim = 14: 1 values converged, normres = (3.62e-13, 4.95e-10, 1.98e-09, 4.76e-01, 7.88e-04)
[ Info: Arnoldi iteration step 15: normres = 0.03808466922525542
[ Info: Arnoldi schursolve in iter 1, krylovdim = 15: 1 values converged, normres = (6.57e-15, 2.73e-11, 1.28e-02, 5.44e-03, 1.50e-06)
[ Info: Arnoldi iteration step 16: normres = 0.6574648457470813
[ Info: Arnoldi schursolve in iter 1, krylovdim = 16: 1 values converged, normres = (3.30e-15, 4.47e-11, 2.62e-02, 2.21e-02, 6.56e-01)
[ Info: Arnoldi iteration step 17: normres = 0.04087377213330988
[ Info: Arnoldi schursolve in iter 1, krylovdim = 17: 1 values converged, normres = (8.04e-17, 4.63e-03, 2.02e-02, 5.72e-05, 1.96e-03)
[ Info: Arnoldi iteration step 18: normres = 0.025578705571474993
[ Info: Arnoldi schursolve in iter 1, krylovdim = 18: 1 values converged, normres = (8.56e-19, 1.74e-09, 3.54e-04, 1.46e-07, 3.67e-05)
[ Info: Arnoldi iteration step 19: normres = 0.24732182361432078
[ Info: Arnoldi schursolve in iter 1, krylovdim = 19: 1 values converged, normres = (9.40e-20, 3.57e-05, 5.44e-05, 1.69e-06, 7.34e-06)
[ Info: Arnoldi iteration step 20: normres = 0.03802503237314808
[ Info: Arnoldi schursolve in iter 1, krylovdim = 20: 1 values converged, normres = (1.66e-21, 8.42e-07, 1.83e-06, 6.10e-08, 2.56e-07)
[ Info: Arnoldi iteration step 21: normres = 0.029537514741369317
[ Info: Arnoldi schursolve in iter 1, krylovdim = 21: 1 values converged, normres = (2.03e-23, 4.99e-10, 3.94e-08, 5.51e-09, 1.45e-09)
[ Info: Arnoldi iteration step 22: normres = 0.05048572624272137
[ Info: Arnoldi schursolve in iter 1, krylovdim = 22: 1 values converged, normres = (4.23e-25, 7.43e-11, 1.31e-09, 6.11e-13, 2.08e-10)
[ Info: Arnoldi iteration step 23: normres = 0.43762140957050844
[ Info: Arnoldi schursolve in iter 1, krylovdim = 23: 1 values converged, normres = (1.55e-25, 1.05e-10, 2.01e-09, 1.35e-12, 6.27e-10)
[ Info: Arnoldi iteration step 24: normres = 0.017259277123691037
[ Info: Arnoldi schursolve in iter 1, krylovdim = 24: 1 values converged, normres = (1.25e-27, 2.56e-12, 4.59e-11, 4.45e-09, 2.16e-09)
[ Info: Arnoldi iteration step 25: normres = 0.11824216803456712
[ Info: Arnoldi schursolve in iter 1, krylovdim = 25: 2 values converged, normres = (6.24e-29, 2.09e-13, 3.72e-12, 2.59e-04, 9.23e-05)
[ Info: Arnoldi iteration step 26: normres = 0.025587943699987167
[ Info: Arnoldi schursolve in iter 1, krylovdim = 26: 3 values converged, normres = (7.04e-31, 3.94e-15, 7.02e-14, 4.98e-07, 6.18e-07)
[ Info: Arnoldi iteration step 27: normres = 0.018591990761149514
[ Info: Arnoldi schursolve in iter 1, krylovdim = 27: 3 values converged, normres = (5.37e-33, 4.78e-17, 8.53e-16, 6.78e-09, 4.19e-08)
[ Info: Arnoldi iteration step 28: normres = 0.06947518389816214
[ Info: Arnoldi schursolve in iter 1, krylovdim = 28: 3 values converged, normres = (1.55e-34, 2.21e-18, 3.94e-17, 2.23e-09, 4.86e-09)
[ Info: Arnoldi iteration step 29: normres = 0.09871929720426023
[ Info: Arnoldi schursolve in iter 1, krylovdim = 29: 3 values converged, normres = (6.72e-36, 1.60e-19, 2.85e-18, 3.95e-10, 2.75e-10)
[ Info: Arnoldi iteration step 30: normres = 0.09607205166889601
[ Info: Arnoldi schursolve in iter 1, krylovdim = 30: 3 values converged, normres = (3.15e-37, 1.34e-20, 2.39e-19, 3.71e-11, 2.76e-11)
[ Info: Arnoldi schursolve in iter 2, krylovdim = 19: 3 values converged, normres = (3.15e-37, 1.34e-20, 2.39e-19, 3.71e-11, 2.76e-11)
[ Info: Arnoldi iteration step 20: normres = 0.06054022007587351
[ Info: Arnoldi schursolve in iter 2, krylovdim = 20: 3 values converged, normres = (7.97e-39, 5.48e-22, 9.78e-21, 1.69e-12, 1.23e-12)
[ Info: Arnoldi iteration step 21: normres = 0.055301973383325057
┌ Info: Arnoldi eigsolve finished after 2 iterations:
│ *  6 eigenvalues converged
│ *  norm of residuals = (1.9373515906125374e-40, 2.217481810035541e-23, 3.9540648182649653e-22, 6.883645113096308e-14, 6.883645113096308e-14, 5.2097446049594056e-14)
└ *  number of operations = 32