Polarizability by linear response
We compute the polarizability of a Helium atom. The polarizability is defined as the change in dipole moment
\[μ = ∫ r ρ(r) dr\]
with respect to a small uniform electric field $E = -x$.
We compute this in two ways: first by finite differences (applying a finite electric field), then by linear response. Note that DFTK is not really adapted to isolated atoms because it uses periodic boundary conditions. Nevertheless we can simply embed the Helium atom in a large enough box (although this is computationally wasteful).
As in other tests, this is not fully converged, convergence parameters were simply selected for fast execution on CI,
using DFTK
using LinearAlgebra
using PseudoPotentialData
a = 10.
lattice = a * I(3) # cube of ``a`` bohrs
pseudopotentials = PseudoFamily("cp2k.nc.sr.lda.v0_1.largecore.gth")
# Helium at the center of the box
atoms = [ElementPsp(:He, pseudopotentials)]
positions = [[1/2, 1/2, 1/2]]
kgrid = [1, 1, 1] # no k-point sampling for an isolated system
Ecut = 30
tol = 1e-8
# dipole moment of a given density (assuming the current geometry)
function dipole(basis, ρ)
rr = [(r[1] - a/2) for r in r_vectors_cart(basis)]
sum(rr .* ρ) * basis.dvol
end;Using finite differences
We first compute the polarizability by finite differences. First compute the dipole moment at rest:
model = model_DFT(lattice, atoms, positions;
functionals=LDA(), symmetries=false)
basis = PlaneWaveBasis(model; Ecut, kgrid)
scfres = self_consistent_field(basis; tol)
μref = dipole(basis, scfres.ρ)-0.00013457376814935573Then in a small uniform field:
ε = .01
model_ε = model_DFT(lattice, atoms, positions;
functionals=LDA(),
extra_terms=[ExternalFromReal(r -> -ε * (r[1] - a/2))],
symmetries=false)
basis_ε = PlaneWaveBasis(model_ε; Ecut, kgrid)
res_ε = self_consistent_field(basis_ε; tol)
με = dipole(basis_ε, res_ε.ρ)0.017612220294307476polarizability = (με - μref) / ε
println("Reference dipole: $μref")
println("Displaced dipole: $με")
println("Polarizability : $polarizability")Reference dipole: -0.00013457376814935573
Displaced dipole: 0.017612220294307476
Polarizability : 1.774679406245683The result on more converged grids is very close to published results. For example DOI 10.1039/C8CP03569E quotes 1.65 with LSDA and 1.38 with CCSD(T).
Using linear response
Now we use linear response (also known as density-functional perturbation theory) to compute this analytically; we refer to standard textbooks for the formalism. In the following, $χ_0$ is the independent-particle polarizability, and $K$ the Hartree-exchange-correlation kernel. We denote with $δV_{\rm ext}$ an external perturbing potential (like in this case the uniform electric field).
# `δVext` is the potential from a uniform field interacting with the dielectric dipole
# of the density.
δVext = [-(r[1] - a/2) for r in r_vectors_cart(basis)]
δVext = cat(δVext; dims=4)54×54×54×1 Array{Float64, 4}:
[:, :, 1, 1] =
5.0 5.0 5.0 5.0 … 5.0 5.0 5.0
4.81481 4.81481 4.81481 4.81481 4.81481 4.81481 4.81481
4.62963 4.62963 4.62963 4.62963 4.62963 4.62963 4.62963
4.44444 4.44444 4.44444 4.44444 4.44444 4.44444 4.44444
4.25926 4.25926 4.25926 4.25926 4.25926 4.25926 4.25926
4.07407 4.07407 4.07407 4.07407 … 4.07407 4.07407 4.07407
3.88889 3.88889 3.88889 3.88889 3.88889 3.88889 3.88889
3.7037 3.7037 3.7037 3.7037 3.7037 3.7037 3.7037
3.51852 3.51852 3.51852 3.51852 3.51852 3.51852 3.51852
3.33333 3.33333 3.33333 3.33333 3.33333 3.33333 3.33333
⋮ ⋱
-3.33333 -3.33333 -3.33333 -3.33333 … -3.33333 -3.33333 -3.33333
-3.51852 -3.51852 -3.51852 -3.51852 -3.51852 -3.51852 -3.51852
-3.7037 -3.7037 -3.7037 -3.7037 -3.7037 -3.7037 -3.7037
-3.88889 -3.88889 -3.88889 -3.88889 -3.88889 -3.88889 -3.88889
-4.07407 -4.07407 -4.07407 -4.07407 -4.07407 -4.07407 -4.07407
-4.25926 -4.25926 -4.25926 -4.25926 … -4.25926 -4.25926 -4.25926
-4.44444 -4.44444 -4.44444 -4.44444 -4.44444 -4.44444 -4.44444
-4.62963 -4.62963 -4.62963 -4.62963 -4.62963 -4.62963 -4.62963
-4.81481 -4.81481 -4.81481 -4.81481 -4.81481 -4.81481 -4.81481
[:, :, 2, 1] =
5.0 5.0 5.0 5.0 … 5.0 5.0 5.0
4.81481 4.81481 4.81481 4.81481 4.81481 4.81481 4.81481
4.62963 4.62963 4.62963 4.62963 4.62963 4.62963 4.62963
4.44444 4.44444 4.44444 4.44444 4.44444 4.44444 4.44444
4.25926 4.25926 4.25926 4.25926 4.25926 4.25926 4.25926
4.07407 4.07407 4.07407 4.07407 … 4.07407 4.07407 4.07407
3.88889 3.88889 3.88889 3.88889 3.88889 3.88889 3.88889
3.7037 3.7037 3.7037 3.7037 3.7037 3.7037 3.7037
3.51852 3.51852 3.51852 3.51852 3.51852 3.51852 3.51852
3.33333 3.33333 3.33333 3.33333 3.33333 3.33333 3.33333
⋮ ⋱
-3.33333 -3.33333 -3.33333 -3.33333 … -3.33333 -3.33333 -3.33333
-3.51852 -3.51852 -3.51852 -3.51852 -3.51852 -3.51852 -3.51852
-3.7037 -3.7037 -3.7037 -3.7037 -3.7037 -3.7037 -3.7037
-3.88889 -3.88889 -3.88889 -3.88889 -3.88889 -3.88889 -3.88889
-4.07407 -4.07407 -4.07407 -4.07407 -4.07407 -4.07407 -4.07407
-4.25926 -4.25926 -4.25926 -4.25926 … -4.25926 -4.25926 -4.25926
-4.44444 -4.44444 -4.44444 -4.44444 -4.44444 -4.44444 -4.44444
-4.62963 -4.62963 -4.62963 -4.62963 -4.62963 -4.62963 -4.62963
-4.81481 -4.81481 -4.81481 -4.81481 -4.81481 -4.81481 -4.81481
[:, :, 3, 1] =
5.0 5.0 5.0 5.0 … 5.0 5.0 5.0
4.81481 4.81481 4.81481 4.81481 4.81481 4.81481 4.81481
4.62963 4.62963 4.62963 4.62963 4.62963 4.62963 4.62963
4.44444 4.44444 4.44444 4.44444 4.44444 4.44444 4.44444
4.25926 4.25926 4.25926 4.25926 4.25926 4.25926 4.25926
4.07407 4.07407 4.07407 4.07407 … 4.07407 4.07407 4.07407
3.88889 3.88889 3.88889 3.88889 3.88889 3.88889 3.88889
3.7037 3.7037 3.7037 3.7037 3.7037 3.7037 3.7037
3.51852 3.51852 3.51852 3.51852 3.51852 3.51852 3.51852
3.33333 3.33333 3.33333 3.33333 3.33333 3.33333 3.33333
⋮ ⋱
-3.33333 -3.33333 -3.33333 -3.33333 … -3.33333 -3.33333 -3.33333
-3.51852 -3.51852 -3.51852 -3.51852 -3.51852 -3.51852 -3.51852
-3.7037 -3.7037 -3.7037 -3.7037 -3.7037 -3.7037 -3.7037
-3.88889 -3.88889 -3.88889 -3.88889 -3.88889 -3.88889 -3.88889
-4.07407 -4.07407 -4.07407 -4.07407 -4.07407 -4.07407 -4.07407
-4.25926 -4.25926 -4.25926 -4.25926 … -4.25926 -4.25926 -4.25926
-4.44444 -4.44444 -4.44444 -4.44444 -4.44444 -4.44444 -4.44444
-4.62963 -4.62963 -4.62963 -4.62963 -4.62963 -4.62963 -4.62963
-4.81481 -4.81481 -4.81481 -4.81481 -4.81481 -4.81481 -4.81481
;;; …
[:, :, 52, 1] =
5.0 5.0 5.0 5.0 … 5.0 5.0 5.0
4.81481 4.81481 4.81481 4.81481 4.81481 4.81481 4.81481
4.62963 4.62963 4.62963 4.62963 4.62963 4.62963 4.62963
4.44444 4.44444 4.44444 4.44444 4.44444 4.44444 4.44444
4.25926 4.25926 4.25926 4.25926 4.25926 4.25926 4.25926
4.07407 4.07407 4.07407 4.07407 … 4.07407 4.07407 4.07407
3.88889 3.88889 3.88889 3.88889 3.88889 3.88889 3.88889
3.7037 3.7037 3.7037 3.7037 3.7037 3.7037 3.7037
3.51852 3.51852 3.51852 3.51852 3.51852 3.51852 3.51852
3.33333 3.33333 3.33333 3.33333 3.33333 3.33333 3.33333
⋮ ⋱
-3.33333 -3.33333 -3.33333 -3.33333 … -3.33333 -3.33333 -3.33333
-3.51852 -3.51852 -3.51852 -3.51852 -3.51852 -3.51852 -3.51852
-3.7037 -3.7037 -3.7037 -3.7037 -3.7037 -3.7037 -3.7037
-3.88889 -3.88889 -3.88889 -3.88889 -3.88889 -3.88889 -3.88889
-4.07407 -4.07407 -4.07407 -4.07407 -4.07407 -4.07407 -4.07407
-4.25926 -4.25926 -4.25926 -4.25926 … -4.25926 -4.25926 -4.25926
-4.44444 -4.44444 -4.44444 -4.44444 -4.44444 -4.44444 -4.44444
-4.62963 -4.62963 -4.62963 -4.62963 -4.62963 -4.62963 -4.62963
-4.81481 -4.81481 -4.81481 -4.81481 -4.81481 -4.81481 -4.81481
[:, :, 53, 1] =
5.0 5.0 5.0 5.0 … 5.0 5.0 5.0
4.81481 4.81481 4.81481 4.81481 4.81481 4.81481 4.81481
4.62963 4.62963 4.62963 4.62963 4.62963 4.62963 4.62963
4.44444 4.44444 4.44444 4.44444 4.44444 4.44444 4.44444
4.25926 4.25926 4.25926 4.25926 4.25926 4.25926 4.25926
4.07407 4.07407 4.07407 4.07407 … 4.07407 4.07407 4.07407
3.88889 3.88889 3.88889 3.88889 3.88889 3.88889 3.88889
3.7037 3.7037 3.7037 3.7037 3.7037 3.7037 3.7037
3.51852 3.51852 3.51852 3.51852 3.51852 3.51852 3.51852
3.33333 3.33333 3.33333 3.33333 3.33333 3.33333 3.33333
⋮ ⋱
-3.33333 -3.33333 -3.33333 -3.33333 … -3.33333 -3.33333 -3.33333
-3.51852 -3.51852 -3.51852 -3.51852 -3.51852 -3.51852 -3.51852
-3.7037 -3.7037 -3.7037 -3.7037 -3.7037 -3.7037 -3.7037
-3.88889 -3.88889 -3.88889 -3.88889 -3.88889 -3.88889 -3.88889
-4.07407 -4.07407 -4.07407 -4.07407 -4.07407 -4.07407 -4.07407
-4.25926 -4.25926 -4.25926 -4.25926 … -4.25926 -4.25926 -4.25926
-4.44444 -4.44444 -4.44444 -4.44444 -4.44444 -4.44444 -4.44444
-4.62963 -4.62963 -4.62963 -4.62963 -4.62963 -4.62963 -4.62963
-4.81481 -4.81481 -4.81481 -4.81481 -4.81481 -4.81481 -4.81481
[:, :, 54, 1] =
5.0 5.0 5.0 5.0 … 5.0 5.0 5.0
4.81481 4.81481 4.81481 4.81481 4.81481 4.81481 4.81481
4.62963 4.62963 4.62963 4.62963 4.62963 4.62963 4.62963
4.44444 4.44444 4.44444 4.44444 4.44444 4.44444 4.44444
4.25926 4.25926 4.25926 4.25926 4.25926 4.25926 4.25926
4.07407 4.07407 4.07407 4.07407 … 4.07407 4.07407 4.07407
3.88889 3.88889 3.88889 3.88889 3.88889 3.88889 3.88889
3.7037 3.7037 3.7037 3.7037 3.7037 3.7037 3.7037
3.51852 3.51852 3.51852 3.51852 3.51852 3.51852 3.51852
3.33333 3.33333 3.33333 3.33333 3.33333 3.33333 3.33333
⋮ ⋱
-3.33333 -3.33333 -3.33333 -3.33333 … -3.33333 -3.33333 -3.33333
-3.51852 -3.51852 -3.51852 -3.51852 -3.51852 -3.51852 -3.51852
-3.7037 -3.7037 -3.7037 -3.7037 -3.7037 -3.7037 -3.7037
-3.88889 -3.88889 -3.88889 -3.88889 -3.88889 -3.88889 -3.88889
-4.07407 -4.07407 -4.07407 -4.07407 -4.07407 -4.07407 -4.07407
-4.25926 -4.25926 -4.25926 -4.25926 … -4.25926 -4.25926 -4.25926
-4.44444 -4.44444 -4.44444 -4.44444 -4.44444 -4.44444 -4.44444
-4.62963 -4.62963 -4.62963 -4.62963 -4.62963 -4.62963 -4.62963
-4.81481 -4.81481 -4.81481 -4.81481 -4.81481 -4.81481 -4.81481Then:
\[δρ = χ_0 δV = χ_0 (δV_{\rm ext} + K δρ),\]
which implies
\[δρ = (1-χ_0 K)^{-1} χ_0 δV_{\rm ext}.\]
From this we identify the polarizability operator to be $χ = (1-χ_0 K)^{-1} χ_0$. Numerically, we apply $χ$ to $δV = -x$ by solving a linear equation (the Dyson equation) iteratively.
First we do this using the DFTK.solve_ΩplusK_split function provided by DFTK, which uses an adaptive Krylov subspace algorithm [HS2025]:
# Multiply δVext times the Bloch waves, then solve the Dyson equation:
δVψ = DFTK.multiply_ψ_by_blochwave(scfres.basis, scfres.ψ, δVext)
res = DFTK.solve_ΩplusK_split(scfres, δVψ; verbose=true)(δψ = Matrix{ComplexF64}[[0.003086105290411181 - 0.0011570959982453882im 0.0 + 0.0im 0.0 + 0.0im 0.0 + 0.0im; -0.1061840013058762 - 0.29060480557902063im 0.0 + 0.0im 0.0 + 0.0im 0.0 + 0.0im; … ; -0.014261699542761776 - 0.03850909985461951im 0.0 + 0.0im 0.0 + 0.0im 0.0 + 0.0im; 0.02870521723520972 + 0.07567696075227925im 0.0 + 0.0im 0.0 + 0.0im 0.0 + 0.0im]], δρ = [-3.620694192637767e-7 -2.502014760795852e-7 … -4.880643727687691e-8 -2.502015133490809e-7; -6.293616992167283e-7 -4.793036444734016e-7 … -1.207784494215547e-7 -4.793036721617347e-7; … ; 1.3626262689213901e-7 1.135623144093002e-7 … 4.711990334610914e-8 1.1356230699659275e-7; 1.3305481051766328e-7 1.3632645139727682e-7 … 5.200685619737264e-8 1.3632642429770043e-7;;; -2.502014859287184e-7 -1.7404904568304266e-7 … -3.5564946394605205e-8 -1.740490732859188e-7; -4.793036518666546e-7 -3.679731948919663e-7 … -1.1202895624517107e-7 -3.679732166417522e-7; … ; 1.1356231083901567e-7 1.1211099700212697e-7 … 1.0980683834599898e-7 1.1211099055939183e-7; 1.3632644386541226e-7 1.2550869925628817e-7 … 5.094908932152392e-8 1.255086791745632e-7;;; -4.880642465777934e-8 -3.556493511234006e-8 … -1.1412150521177913e-8 -3.556494149834413e-8; -1.207784394550953e-7 -1.1202894440474266e-7 … -9.78622575903291e-8 -1.1202895109690401e-7; … ; 4.711990555321552e-8 1.098068466486683e-7 … 2.493823176670479e-7 1.0980684104518815e-7; 5.200686473326584e-8 5.094909876314189e-8 … 4.9106280598591004e-8 5.0949093007720866e-8;;; … ;;; 6.475761304515082e-8 4.110333411153984e-8 … -6.681699712281718e-9 4.1103341600903146e-8; 1.6897558284250845e-7 6.65328524235346e-8 … -1.3294452340299122e-7 6.653285695288105e-8; … ; 3.2677650417742836e-8 1.6017370319622466e-7 … 3.847036706605293e-7 1.6017369949577377e-7; -9.141992511263444e-8 -3.583034739793322e-8 … 7.037284955515127e-8 -3.5830343611366475e-8;;; -4.880643261950865e-8 -3.5564941090490256e-8 … -1.1412152204529989e-8 -3.556494717413161e-8; -1.2077844559219403e-7 -1.1202895086688443e-7 … -9.786225983408922e-8 -1.120289550272207e-7; … ; 4.711990573807169e-8 1.0980684307465212e-7 … 2.493823087513038e-7 1.0980683867622185e-7; 5.20068594943182e-8 5.0949093455722647e-8 … 4.9106276390918515e-8 5.094908834658089e-8;;; -2.5020150329360455e-7 -1.7404905876402702e-7 … -3.556494881417711e-8 -1.7404908545334388e-7; -4.793036645260742e-7 -3.6797320527834933e-7 … -1.1202895641700717e-7 -3.6797322443274635e-7; … ; 1.1356230928825263e-7 1.1211099553032628e-7 … 1.0980683565957902e-7 1.1211098878016676e-7; 1.363264315214469e-7 1.25508689892893e-7 … 5.094908645443027e-8 1.2550867009846545e-7;;;;], δHtotψ = Matrix{ComplexF64}[[-0.0017951485785370322 + 0.0006730680548543059im -0.08736970920893486 - 0.03805335525222682im 1.2421138554658384 - 1.2575982306233637im -0.6442802040184978 + 0.8394918336716252im; 0.06792769376257593 + 0.18630463563227254im -0.7356520460955388 + 1.4510310821134929im 0.27415100693777916 - 0.27806402406417624im -0.14214754324660905 + 0.18553666778698302im; … ; 0.016899605218306756 + 0.04601474216746307im 0.008212674135929134 - 0.018986238344880425im -0.02113843487506009 + 0.019802665052085673im 0.012785290253770567 - 0.010188154740791761im; -0.02566353639834422 - 0.06767163080241616im -0.01027138405471991 + 0.023439301537901462im -0.029030977500693198 - 0.0008964683746390126im 0.025354880446641415 + 0.03860012510945076im]], δVind = [0.0019654319362370273 0.00208100048618632 … 0.003153665278031989 0.0020810007955972378; 0.00388913156919915 0.0047792700031459535 … 0.010276672343740205 0.004779270340786578; … ; -0.00615084932745164 -0.008269605136681917 … -0.017302969501390122 -0.00826960456848782; 0.00023868169891399496 -0.0003986981600705272 … -0.003664695127877632 -0.0003986978059559139;;; 0.002081000568446646 0.002223228067788227 … 0.0034201198403682243 0.0022232284153020443; 0.0047792700891658305 0.005734954102263238 … 0.010903598659396435 0.005734954472290274; … ; -0.008269605037670281 -0.009601010076288927 … -0.012623112809056331 -0.009601009474764345; -0.00039869806516888684 -0.0010730139644493332 … -0.004230681128428555 -0.0010730135672022646;;; 0.0031536644922369496 0.003420118816931034 … 0.00483250465521556 0.0034201194080780324; 0.010276671488606302 0.010903597674179297 … 0.012272142061642713 0.010903598237975185; … ; -0.017302970856836956 -0.01262311394226925 … -0.00962026654083625 -0.012623113324580068; -0.0036646960341648683 -0.004230682225741562 … -0.005491108971305244 -0.004230681604263915;;; … ;;; -0.002336560703249231 -0.002815463305943854 … 0.005693318280138861 -0.0028154638169175903; -0.012913094724496121 -0.0178411840663333 … 0.011498341400219315 -0.01784118465029509; … ; -0.025174460733356552 -0.01281088770486523 … -0.008395875659882028 -0.012810887052170294; 0.007597979574452099 0.010225415164555647 … -0.0054497331838359465 0.010225414510024159;;; 0.0031536650050040514 0.003420119400885197 … 0.004832505543863694 0.0034201199768617006; 0.010276672070775752 0.010903598253547276 … 0.012272142813644308 0.010903598896321854; … ; -0.01730296974903605 -0.012623113287899983 … -0.009620266082996886 -0.012623112627273588; -0.0036646954533612434 -0.004230681634506811 … -0.005491108357155015 -0.0042306810111414385;;; 0.002081000714735458 0.002223228235667796 … 0.003420120136691302 0.00222322857804726; 0.0047792702574103695 0.005734954288021839 … 0.010903599057035931 0.005734954668044062; … ; -0.008269604711039474 -0.009601009733832384 … -0.01262311253489984 -0.0096010091578582; -0.00039869790212949453 -0.0010730137780339062 … -0.004230680854787067 -0.0010730133887636448;;;;], δρ0 = [-3.602140063638292e-7 -2.4888426084499396e-7 … -4.8528368315700854e-8 -2.4888426031653604e-7; -5.972250272498438e-7 -4.6069782393102805e-7 … -1.1905768802691744e-7 -4.6069782403445845e-7; … ; 1.2735790687313252e-7 1.0836340241468281e-7 … 4.640718744055608e-8 1.0836340139683404e-7; 1.0352009924139939e-7 1.195796359803443e-7 … 5.069506954355462e-8 1.195796365162754e-7;;; -2.4888425942457604e-7 -1.7310671595225957e-7 … -3.535697804992056e-8 -1.7310671587026763e-7; -4.606978230275844e-7 -3.5721707314756636e-7 … -1.1091440164969129e-7 -3.5721707438485756e-7; … ; 1.0836340081406847e-7 1.0846497839460014e-7 … 1.0862237819621342e-7 1.0846497774844038e-7; 1.1957963667758556e-7 1.1611148821389761e-7 … 5.0208313695791237e-8 1.1611148866313148e-7;;; -4.852836764489596e-8 -3.535697758984945e-8 … -1.1341317148839973e-8 -3.53569776174446e-8; -1.1905768695057613e-7 -1.109143996234097e-7 … -9.770669636997782e-8 -1.1091440054732944e-7; … ; 4.640718666071091e-8 1.0862237831979663e-7 … 2.484524733155184e-7 1.0862237754044131e-7; 5.069506916118041e-8 5.02083133058078e-8 … 4.926859358347597e-8 5.020831333537415e-8;;; … ;;; 6.434570118488023e-8 4.083862049625344e-8 … -6.637933417027731e-9 4.0838619932394704e-8; 1.697503452899619e-7 6.683083704985573e-8 … -1.334048307355267e-7 6.683083844537394e-8; … ; 3.271000194050471e-8 1.6023525995655186e-7 … 3.839014881658955e-7 1.6023526245881496e-7; -9.271367169248466e-8 -3.63362066777745e-8 … 7.126310625420428e-8 -3.6336206036421703e-8;;; -4.852836877337097e-8 -3.535697843944827e-8 … -1.1341317123098915e-8 -3.53569782144e-8; -1.1905768801245707e-7 -1.1091440144423585e-7 … -9.770669478297717e-8 -1.1091439989830326e-7; … ; 4.640718852887717e-8 1.0862237860484239e-7 … 2.4845247253981267e-7 1.0862237899642943e-7; 5.069506901766293e-8 5.0208312659070753e-8 … 4.926859322838845e-8 5.0208313289764306e-8;;; -2.4888426160967323e-7 -1.7310671784832178e-7 … -3.5356977934047947e-8 -1.7310671700869622e-7; -4.606978247176176e-7 -3.5721707531186964e-7 … -1.10914399552573e-7 -3.572170741120817e-7; … ; 1.0836340175639972e-7 1.0846497921532083e-7 … 1.0862237744730401e-7 1.0846497825599781e-7; 1.1957963545402026e-7 1.1611148715858239e-7 … 5.02083131489167e-8 1.1611148780187653e-7;;;;], δeigenvalues = [[0.0004950049229463429, 0.0967985937660379, 0.01749262953153934, 0.04025849587659798]], δoccupation = [[0.0, 0.0, 0.0, 0.0]], δεF = 0.0, ε_adj = DFTK.DielectricAdjoint{Array{Float64, 4}, Vector{Matrix{ComplexF64}}, Vector{Vector{Float64}}, Float64, Vector{Vector{Float64}}, StaticArraysCore.SVector{3, Float64}}(Hamiltonian(PlaneWaveBasis(model = Model(lda_x+lda_c_pw, spin_polarization = :none), Ecut = 30.0 Ha, kgrid = MonkhorstPack([1, 1, 1])), HamiltonianBlock[DFTK.DftHamiltonianBlock(PlaneWaveBasis(model = Model(lda_x+lda_c_pw, spin_polarization = :none), Ecut = 30.0 Ha, kgrid = MonkhorstPack([1, 1, 1])), KPoint([ 0, 0, 0], spin = 1, num. G vectors = 7809), Any[DFTK.FourierMultiplication{Float64, Vector{Float64}}(PlaneWaveBasis(model = Model(lda_x+lda_c_pw, spin_polarization = :none), Ecut = 30.0 Ha, kgrid = MonkhorstPack([1, 1, 1])), KPoint([ 0, 0, 0], spin = 1, num. G vectors = 7809), [0.0, 0.19739208802178715, 0.7895683520871486, 1.7765287921960844, 3.1582734083485944, 4.934802200544679, 7.106115168784338, 9.67221231306757, 12.633093633394378, 15.988759129764759 … 20.133992978222288, 16.38354330580833, 13.027877809437953, 10.066996489111146, 7.5008993448279115, 5.329586376588253, 3.5530575843921683, 2.1713129682396586, 1.184352528130723, 0.5921762640653614]), DFTK.RealSpaceMultiplication{Float64, Array{Float64, 3}}(PlaneWaveBasis(model = Model(lda_x+lda_c_pw, spin_polarization = :none), Ecut = 30.0 Ha, kgrid = MonkhorstPack([1, 1, 1])), KPoint([ 0, 0, 0], spin = 1, num. G vectors = 7809), [0.1603906580688546 0.1602471565738079 … 0.15981840036267833 0.16024715657380792; 0.16024715657380803 0.16010321540768782 … 0.15967315066597979 0.16010321540768782; … ; 0.15981840036267841 0.1596731506659797 … 0.15923916395256382 0.1596731506659797; 0.16024715657380809 0.16010321540768785 … 0.1596731506659798 0.16010321540768785;;; 0.16024715657380786 0.16010321540768765 … 0.15967315066597965 0.16010321540768768; 0.16010321540768782 0.1599588345026577 … 0.15952744714353234 0.1599588345026577; … ; 0.15967315066597973 0.15952744714353229 … 0.15909210945052638 0.15952744714353229; 0.16010321540768785 0.15995883450265772 … 0.15952744714353242 0.15995883450265772;;; 0.1598184003626783 0.15967315066597965 … 0.15923916395256382 0.15967315066597965; 0.15967315066597973 0.15952744714353237 … 0.1590921094505264 0.15952744714353237; … ; 0.15923916395256385 0.15909210945052635 … 0.15865272241895537 0.15909210945052635; 0.1596731506659798 0.15952744714353237 … 0.15909210945052646 0.15952744714353237;;; … ;;; 0.15910963392556968 0.15896221057980836 … 0.1585217272926462 0.15896221057980836; 0.15896221057980842 0.15881431920357322 … 0.15837242838914242 0.15881431920357322; … ; 0.15852172729264621 0.15837242838914237 … 0.157926332522532 0.15837242838914237; 0.15896221057980844 0.15881431920357328 … 0.1583724283891425 0.15881431920357328;;; 0.1598184003626783 0.15967315066597962 … 0.15923916395256382 0.15967315066597962; 0.15967315066597973 0.15952744714353237 … 0.1590921094505264 0.15952744714353237; … ; 0.15923916395256382 0.15909210945052635 … 0.15865272241895534 0.15909210945052635; 0.15967315066597979 0.15952744714353237 … 0.15909210945052646 0.15952744714353237;;; 0.16024715657380786 0.16010321540768765 … 0.15967315066597965 0.16010321540768768; 0.16010321540768782 0.15995883450265772 … 0.15952744714353237 0.15995883450265772; … ; 0.1596731506659797 0.15952744714353229 … 0.15909210945052635 0.15952744714353229; 0.16010321540768785 0.15995883450265772 … 0.15952744714353242 0.15995883450265772]), DFTK.NonlocalOperator{Float64, Matrix{ComplexF64}, Matrix{Float64}}(PlaneWaveBasis(model = Model(lda_x+lda_c_pw, spin_polarization = :none), Ecut = 30.0 Ha, kgrid = MonkhorstPack([1, 1, 1])), KPoint([ 0, 0, 0], spin = 1, num. G vectors = 7809), Matrix{ComplexF64}(undef, 7809, 0), Matrix{Float64}(undef, 0, 0)), DFTK.NoopOperator{Float64}(PlaneWaveBasis(model = Model(lda_x+lda_c_pw, spin_polarization = :none), Ecut = 30.0 Ha, kgrid = MonkhorstPack([1, 1, 1])), KPoint([ 0, 0, 0], spin = 1, num. G vectors = 7809)), DFTK.NoopOperator{Float64}(PlaneWaveBasis(model = Model(lda_x+lda_c_pw, spin_polarization = :none), Ecut = 30.0 Ha, kgrid = MonkhorstPack([1, 1, 1])), KPoint([ 0, 0, 0], spin = 1, num. G vectors = 7809)), DFTK.RealSpaceMultiplication{Float64, Array{Float64, 3}}(PlaneWaveBasis(model = Model(lda_x+lda_c_pw, spin_polarization = :none), Ecut = 30.0 Ha, kgrid = MonkhorstPack([1, 1, 1])), KPoint([ 0, 0, 0], spin = 1, num. G vectors = 7809), [-0.15475965282095328 -0.15461620269656168 … -0.1541875309955631 -0.15461620269606846; -0.1546162026962077 -0.15447229889571576 … -0.15404229467250868 -0.15447229889522493; … ; -0.15418753099613522 -0.15404229467344036 … -0.15360833734366908 -0.1540422946729333; -0.1546162026963786 -0.15447229889589084 … -0.15404229467267314 -0.15447229889539235;;; -0.15461620269617232 -0.15447229889568018 … -0.15404229467247224 -0.1544722988951892; -0.15447229889532887 -0.15432794308685013 … -0.15389660282505388 -0.1543279430863616; … ; -0.15404229467303912 -0.15389660282597772 … -0.1534612870072303 -0.15389660282547277; -0.15447229889549677 -0.15432794308702233 … -0.15389660282521528 -0.15432794308652598;;; -0.1541875309957527 -0.1540422946730491 … -0.15360833734329768 -0.15404229467255737; -0.1540422946726981 -0.153896602825628 … -0.15346128700690168 -0.15389660282513884; … ; -0.15360833734386295 -0.15346128700782377 … -0.15302191313176303 -0.1534612870073177; -0.15404229467286418 -0.1538966028257986 … -0.15346128700706116 -0.15389660282530127;;; … ;;; -0.15347880585224682 -0.15333139196237391 … -0.15289092492464343 -0.15333139196185602; -0.1533313919619985 -0.1531835057401594 … -0.15274162877626715 -0.15318350573964404; … ; -0.15289092492525325 -0.15274162877725933 … -0.15229553958790035 -0.1527416287767257; -0.15333139196218573 -0.15318350574035158 … -0.15274162877644806 -0.15318350573982759;;; -0.15418753099628738 -0.15404229467359304 … -0.15360833734382454 -0.15404229467308642; -0.15404229467322714 -0.15389660282616627 … -0.15346128700742268 -0.15389660282566214; … ; -0.15360833734441814 -0.15346128700838865 … -0.1530219131323108 -0.15346128700786718; -0.15404229467340758 -0.1538966028263514 … -0.15346128700759695 -0.15389660282583897;;; -0.15461620269644136 -0.1544722988959538 … -0.1540422946727373 -0.15447229889545544; -0.15447229889559508 -0.15432794308712092 … -0.15389660282531611 -0.154327943086625; … ; -0.15404229467331837 -0.1538966028262618 … -0.15346128700750594 -0.15389660282574927; -0.1544722988957701 -0.15432794308730036 … -0.15389660282548484 -0.1543279430867965]), DFTK.RealSpaceMultiplication{Float64, SubArray{Float64, 3, Array{Float64, 4}, Tuple{Base.Slice{Base.OneTo{Int64}}, Base.Slice{Base.OneTo{Int64}}, Base.Slice{Base.OneTo{Int64}}, Int64}, true}}(PlaneWaveBasis(model = Model(lda_x+lda_c_pw, spin_polarization = :none), Ecut = 30.0 Ha, kgrid = MonkhorstPack([1, 1, 1])), KPoint([ 0, 0, 0], spin = 1, num. G vectors = 7809), [-0.015163897871257156 -0.012550129343332432 … -0.00496686855889924 -0.012550129360809362; -0.012550129363414285 -0.01036259532951398 … -0.004135159407730051 -0.010362595346263371; … ; -0.004966868515480071 -0.004135159362679929 … -0.0020697807685601593 -0.004135159336175332; -0.012550129335857398 -0.010362595307251952 … -0.0041351593459770936 -0.010362595311789069;;; -0.01255012934835796 -0.010362595312690863 … -0.004135159400965445 -0.01036259533618034; -0.01036259533756688 -0.008622039287500279 … -0.00389751185081905 -0.008622039317424746; … ; -0.004135159319498948 -0.0038975117580871414 … -0.0034841073667349173 -0.003897511741963873; -0.010362595297572554 -0.008622039257698777 … -0.0038975117675891617 -0.008622039265095605;;; -0.004966868531397047 -0.004135159360425766 … -0.0020697808573993132 -0.004135159377017568; -0.004135159389323934 -0.0038975118036918806 … -0.003484107460247089 -0.003897511829369462; … ; -0.002069780748413735 -0.003484107370289544 … -0.005620580291923935 -0.0034841073523840353; -0.004135159316644481 -0.0038975117500088074 … -0.003484107376074801 -0.0038975117453514405;;; … ;;; -0.006547897686549051 -0.004876339750342766 … -0.0014866449232484062 -0.004876339746281544; -0.0048763397244044986 -0.0026779881993586056 … -0.004123054710203007 -0.002677988239925989; … ; -0.0014866449048345806 -0.00412305467733273 … -0.007170895568668758 -0.0041230547127148265; -0.0048763397591903394 -0.0026779882891254997 … -0.00412305470671637 -0.0026779882498727658;;; -0.004966868535264666 -0.004135159357523185 … -0.002069780822720476 -0.004135159366263271; -0.004135159391280644 -0.0038975118229741534 … -0.0034841074084658134 -0.003897511796225835; … ; -0.002069780806209806 -0.003484107386512212 … -0.005620580294957866 -0.0034841073909213584; -0.004135159344643837 -0.003897511752403645 … -0.00348410738195672 -0.003897511771283527;;; -0.012550129347418747 -0.010362595314499859 … -0.004135159366441559 -0.010362595324189712; -0.010362595333859196 -0.008622039296144624 … -0.0038975117952779765 -0.008622039293153354; … ; -0.004135159350228692 -0.003897511784546709 … -0.0034841073559809775 -0.003897511759725165; -0.010362595307630448 -0.00862203926515315 … -0.003897511754135778 -0.00862203926873414])], DFTK.FourierMultiplication{Float64, Vector{Float64}}(PlaneWaveBasis(model = Model(lda_x+lda_c_pw, spin_polarization = :none), Ecut = 30.0 Ha, kgrid = MonkhorstPack([1, 1, 1])), KPoint([ 0, 0, 0], spin = 1, num. G vectors = 7809), [0.0, 0.19739208802178715, 0.7895683520871486, 1.7765287921960844, 3.1582734083485944, 4.934802200544679, 7.106115168784338, 9.67221231306757, 12.633093633394378, 15.988759129764759 … 20.133992978222288, 16.38354330580833, 13.027877809437953, 10.066996489111146, 7.5008993448279115, 5.329586376588253, 3.5530575843921683, 2.1713129682396586, 1.184352528130723, 0.5921762640653614]), DFTK.RealSpaceMultiplication{Float64, Array{Float64, 3}}(PlaneWaveBasis(model = Model(lda_x+lda_c_pw, spin_polarization = :none), Ecut = 30.0 Ha, kgrid = MonkhorstPack([1, 1, 1])), KPoint([ 0, 0, 0], spin = 1, num. G vectors = 7809), [-0.009532892623355823 -0.006919175466086218 … 0.0006640008082159973 -0.006919175483069903; -0.006919175485813965 -0.004731678817541919 … 0.0014956965857410498 -0.00473167883380048; … ; 0.0006640008510631241 0.0014956966298594176 … 0.0035610458403345863 0.0014956966568710813; -0.006919175458427913 -0.0047316787954549445 … 0.0014956966473295835 -0.004731678799493572;;; -0.006919175470722419 -0.0047316788006833855 … 0.0014956965925419603 -0.004731678823681865; -0.004731678825207933 -0.0029911478716927245 … 0.0017333324676594118 -0.002991147901128666; … ; 0.0014956966734416604 0.0017333325594674205 … 0.0021467150765611567 0.0017333325760956462; -0.004731678785381473 -0.002991147842063391 … 0.0017333325507279847 -0.0029911478489638658;;; 0.0006640008355285644 0.0014956966325047677 … 0.0035610457518668297 0.0014956966164047115; 0.0014956966039577071 0.0017333325142124848 … 0.002146714983377639 0.0017333324890240675; … ; 0.003561045860287166 0.0021467150724130326 … 1.0228995268405164e-5 0.0021467150908246087; 0.0014956966764711546 0.0017333325677249723 … 0.002146715067390499 0.0017333325728796636;;; … ;;; -0.0009170696132261985 0.0007544788670916808 … 0.004144157444754351 0.000754478871670794; 0.0007544788934054254 0.002952825264055214 … 0.0015077449026722634 0.0029528252240031963; … ; 0.0041441574625583864 0.0015077449345503062 … -0.0015401026340371125 0.0015077448997018383; 0.0007544788584323732 0.002952825174096196 … 0.001507744905978073 0.0029528252138729272;;; 0.0006640008311262621 0.0014956966348633948 … 0.0035610457860188105 0.0014956966266299315; 0.0014956966014719481 0.0017333324943919481 … 0.0021467150346379145 0.001733332521644391; … ; 0.003561045801935872 0.002146715055625483 … 1.0228991686661755e-5 0.0021467150517378085; 0.001495696647928372 0.001733332564777327 … 0.0021467150609727863 0.0017333325464098685;;; -0.00691917547005224 -0.004731678802765996 … 0.001495696626800781 -0.00473167881195747; -0.004731678821766453 -0.002991147880607826 … 0.0017333325229382784 -0.0029911478771206465; … ; 0.0014956966424326404 0.0017333325327237746 … 0.002146715087039428 0.0017333325580578532; -0.004731678795712704 -0.0029911478497957906 … 0.0017333325639118063 -0.002991147852872935]), DFTK.NonlocalOperator{Float64, Matrix{ComplexF64}, Matrix{Float64}}(PlaneWaveBasis(model = Model(lda_x+lda_c_pw, spin_polarization = :none), Ecut = 30.0 Ha, kgrid = MonkhorstPack([1, 1, 1])), KPoint([ 0, 0, 0], spin = 1, num. G vectors = 7809), Matrix{ComplexF64}(undef, 7809, 0), Matrix{Float64}(undef, 0, 0)), nothing, @NamedTuple{ψ_reals::Array{ComplexF64, 3}}[(ψ_reals = [-2.3894850630358637e-6 + 1.2301039032956905e-6im -0.0034691840125510437 + 0.003218268298298im … -0.0024408342683031586 + 0.002225347835158621im 0.0035427070285445348 - 0.0032508093614383657im; 0.005446036002492062 + 0.014621225688492256im -0.0032797963181838037 - 0.01609763126705992im … 0.005095922915645382 + 0.007581066317841097im -0.007388368047900649 - 0.01230640697894386im; … ; 0.0071043300433178145 + 0.019090842422228364im -0.007057918340761343 - 0.01643277370737064im … 0.003381959547015621 + 0.011978496011586917im -0.005647910611889299 - 0.017734140564302735im; -0.005513429244330443 - 0.014608444470437902im 0.007373895605055171 + 0.012318617370514763im … -0.0017355415277248872 - 0.010639637724915938im 0.0032657318174173656 + 0.016109208071864174im;;; -0.0024183427707357265 + 0.0017239999961991873im 0.0034452793507601344 - 0.002894647062549234im … 0.0011954205506182053 - 0.0011837456011370723im -0.0006628086323776488 + 0.0008963664542009486im; -0.0038948803562302947 - 0.015220859296052804im 0.0024760174833716258 + 0.014398670282577034im … -0.0037730578518784938 - 0.006851189123928511im 0.005120873666843435 + 0.011957653476711643im; … ; -0.006847444867776014 - 0.016734825212318706im 0.006468395821761871 + 0.014090921464650394im … -0.0032036214394775765 - 0.010288073511720728im 0.005361644243810327 + 0.015112535964660159im; 0.0067590535981711455 + 0.013195758825339243im -0.006922328812321617 - 0.010669655798419017im … 0.0020177976668783477 + 0.008594439961496761im -0.004277677177735901 - 0.013110409365145577im;;; 0.0016233548150387377 - 0.0011677630406808784im -0.0018062836985181913 + 0.0014472485412943535im … -0.0002789902054406848 + 0.000352623189434957im -0.0003963728397903436 + 0.00014567994869362335im; 0.002304559396409417 + 0.009910447387903969im -0.0015242369097623837 - 0.008826400316422662im … 0.0021018272870632734 + 0.004674074512639892im -0.002630976063917714 - 0.007804737408555024im; … ; 0.004723478569204795 + 0.01086031218555348im -0.004422839489040025 - 0.009214614510536545im … 0.002400038578083634 + 0.007042650797405621im -0.003786744295222202 - 0.009801879250888788im; -0.004526766829674342 - 0.008310933030951427im 0.0042666149611757884 + 0.006619597547530127im … -0.0016257199729525626 - 0.005268291121080872im 0.0031598836037309425 + 0.007641264580009149im;;; … ;;; 0.0007850861457931422 - 0.0005341652217499027im -0.0003908082463009247 + 0.0001971191134623233im … 0.0008375239087787751 - 0.0006592053388792963im -0.0009449293209043178 + 0.0007087364384832401im; -0.0022987921443275005 - 0.0038989709339311194im 0.0017844602637924695 + 0.0036873010265148594im … -0.001796691108225046 - 0.0021958797057001402im 0.0022653189379703105 + 0.003243266300455668im; … ; -0.0018493910067256835 - 0.006444489812256023im 0.0018411351739562529 + 0.0056611433090626075im … -0.001011609483434994 - 0.0046516265926289406im 0.0015027885976934808 + 0.005973555223794978im; 0.0009643674567199424 + 0.004804376102377546im -0.0011075275289896313 - 0.004026128358957064im … 0.00029282410821076995 + 0.0033772243304117976im -0.0006266816580163977 - 0.0044700928326065256im;;; -0.0017371884517122897 + 0.001193437284194355im 0.0005002160037910507 - 0.00017410077634066423im … -0.0013868212264502955 + 0.0011317544536823465im 0.001910408780584605 - 0.001475757250064867im; 0.004611313875448865 + 0.008289016858702047im -0.0032455069688972476 - 0.007616155006595925im … 0.0030302738173847805 + 0.004021029136921367im -0.00435241093584092 - 0.006594320419984147im; … ; 0.0036156817777015238 + 0.01163908832500956im -0.0034945490599062605 - 0.009867340197710926im … 0.0017942870000669869 + 0.007468632166220042im -0.0028583663430799787 - 0.010454576134179286im; -0.002220036729226765 - 0.009931992497190956im 0.0025453682411769253 + 0.007829570562210458im … -0.0006973576179346286 - 0.005921157522955288im 0.0014385334362177963 + 0.008851252743398385im;;; 0.002511068753052544 - 0.0017384300458820334im 0.0005586767526668014 - 0.0008657499640161654im … 0.0021857075727001727 - 0.0018802394009952601im -0.0035499266853782916 + 0.0029253021295430606im; -0.006781766509868342 - 0.013191818535135421im 0.004333851119076613 + 0.013092259135501627im … -0.004550414162093908 - 0.006304459179758479im 0.006978967250346267 + 0.010651165703846402im; … ; -0.005857093078885048 - 0.01743087916224579im 0.00569113212761467 + 0.014637397292825618im … -0.002756911493863602 - 0.010602265618125578im 0.004584281541452611 + 0.015658977944596513im; 0.0038719273100939603 + 0.015224030714702155im -0.005064388873141863 - 0.011975570713156948im … 0.0012405697957762682 + 0.009141085077133872im -0.0024195003850727434 - 0.014416346758449296im],)])]), [8.519873717954785e-7 4.627265091959304e-7 … 2.4177784152388973e-8 4.627265112702247e-7; 4.627265115793934e-7 2.501447670603227e-7 … 1.3582998744488567e-8 2.5014476835578976e-7; … ; 2.4177783486167918e-8 1.3582998279725298e-8 … 1.565922751217056e-9 1.3582998006286933e-8; 4.627265083087323e-7 2.5014476533850884e-7 … 1.3582998107410813e-8 2.50144765689421e-7;;; 4.62726509792396e-7 2.501447657591625e-7 … 1.3582998674702355e-8 2.501447675759258e-7; 2.501447676831651e-7 1.389560863463933e-7 … 1.1279991472622962e-8 1.3895608788494388e-7; … ; 1.3582997834247522e-8 1.1279990630654264e-8 … 7.937755945307066e-9 1.1279990484260395e-8; 2.501447645898595e-7 1.3895608481418388e-7 … 1.1279990716927468e-8 1.3895608519448126e-7;;; 2.4177783730399885e-8 1.3582998256473004e-8 … 1.5659229598316851e-9 1.3582998427639328e-8; 1.3582998554602958e-8 1.1279991044723043e-8 … 7.937756612297423e-9 1.1279991277869365e-8; … ; 1.5659227039087679e-9 7.937755970659126e-9 … 3.571856348695979e-8 7.937755842944924e-9; 1.3582997804799157e-8 1.1279990557303446e-8 … 7.937756011925534e-9 1.1279990515015223e-8;;; … ;;; 5.79184176847483e-8 2.2816022891189136e-8 … 5.620378839611626e-10 2.281602283131793e-8; 2.2816022508785538e-8 3.4892768771919997e-9 … 1.3458513770111748e-8 3.489277041918545e-9; … ; 5.620378624536113e-10 1.3458513433144405e-8 … 7.728315929597009e-8 1.345851379586336e-8; 2.281602302163141e-8 3.4892772416997605e-9 … 1.3458513734370498e-8 3.489277082310219e-9;;; 2.4177783789741158e-8 1.3582998226524247e-8 … 1.5659228783972343e-9 1.3582998316696202e-8; 1.3582998574787645e-8 1.1279991219802482e-8 … 7.937756242960553e-9 1.1279990976938685e-8; … ; 1.565922839627317e-9 7.937756086373033e-9 … 3.571856354788626e-8 7.937756117821259e-9; 1.3582998093655665e-8 1.1279990579048064e-8 … 7.937756053879729e-9 1.1279990750471632e-8;;; 4.6272650968092497e-7 2.501447658990846e-7 … 1.3582998318533371e-8 2.501447666485367e-7; 2.501447673964001e-7 1.3895608679084648e-7 … 1.1279990968332262e-8 1.3895608663704722e-7; … ; 1.3582998151269956e-8 1.1279990870895478e-8 … 7.93775586859973e-9 1.1279990645523236e-8; 2.5014476536777955e-7 1.3895608519743568e-7 … 1.1279990594777553e-8 1.389560853815573e-7;;;;], Matrix{ComplexF64}[[-0.24393085966210568 + 0.09145876669983086im -0.8847995164545162 - 0.38536920872561103im -1.8606782228309391e-6 - 8.884767594114089e-7im -1.9823771882122014e-7 + 3.803329816397014e-7im; 0.1655019426458059 - 0.062052843936172615im -0.03898773596217184 - 0.016980879901088194im -0.39803692088975334 - 0.3930551488725451im 0.2657508339374689 + 0.20390863982191745im; … ; -0.052581203794761534 + 0.019714652163732176im 0.014101362618809706 + 0.006141765665499026im -0.006918084303946594 - 0.00529399428713321im 0.006508078581161764 + 0.002283686691070428im; 0.09427845853944222 - 0.03534850635263006im -0.02793587841584223 - 0.012167308259109775im 0.017031713286869088 + 0.008282855003658588im -0.021903428950301024 - 0.0017361828032670956im]], [[2.0, 0.0, 0.0, 0.0]], -0.28239440687962375, [[-0.5541104213357841, -0.010678392423463358, 0.1449565831100493, 0.14495671101348204]], 1.0e-6, DFTK.BandtolBalanced{Float64, Vector{Vector{Float64}}}([[0.6300127939257145]], 1.1102230246251565e-16, Inf, 1.0e-6), 100, [0.0, 0.0, 0.0]), info_gmres = (x = [-3.620694192637767e-7, -6.293616992167283e-7, -2.2951212275426885e-7, 4.994026508824494e-7, 7.923313110392089e-7, 5.851072998267395e-7, 6.516971927960387e-7, 1.0189310190271308e-6, 3.638592066626193e-7, -2.201093757990669e-6 … 2.7810187288648138e-6, 1.2016957644471427e-6, 2.445642505661889e-8, 1.0078549369247888e-7, 5.53530892251943e-7, 4.5341735510985634e-7, 4.63193384141252e-9, -1.0437052946901531e-7, 1.1211098878016676e-7, 1.2550867009846545e-7], resid_history = [0.2493920933892154, 0.0037665019476162726, 0.00028511771575586467, 4.682614238453956e-6, 1.6478156526067865e-8, 5.2178509035511174e-11, 3.218652794928696e-8, 1.0663782204138613e-9], converged = true, n_iter = 8, residual_norm = 1.0663782204138613e-9, maxiter = 100, tol = 1.0e-8, s = 0.9346703461464724, residual = [7.409463827276296e-7, 4.943967329583944e-7, 9.96504602163138e-8, -1.3125103879691567e-7, -1.3268980137851054e-7, -6.501758758779599e-8, -4.829897474617991e-8, -4.714216210231375e-8, -8.82473490998314e-9, 1.7489834484352504e-8 … 9.546665532508982e-7, 4.184181277898035e-7, 8.810920764183554e-9, 3.800121959837265e-8, 2.1944544106511103e-7, 1.9252995459475444e-7, 2.220442984960636e-9, -6.258322833246825e-8, 1.025714044160999e-7, 3.270104858166171e-7], restart_history = [6], stage = :finalize, krylovdim = 20, y = [2.760858083312157e-7, -3.17124504180164e-8, 1.0663782204138613e-9, -0.00020434709440023867, 4.6883329503964664e-6, -1.6649727459896717e-8, 5.2178509035511174e-11, 0.0, 0.0, 0.0 … 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0], V = KrylovKit.OrthonormalBasis{Vector{Float64}}([[6.378570950215358e-9, 3.7916907128740044e-7, 5.435196371880748e-8, 1.535397500955027e-7, 6.005745474456382e-7, 6.310162753028188e-7, 8.658226931028641e-7, 1.5849231711862642e-6, 6.391421999527092e-7, -4.194059372527151e-6 … 4.859977006691685e-6, 1.8989454547827627e-6, 3.3708019367489594e-8, 1.2517992955807747e-7, 6.495271468334723e-7, 4.890268603835215e-7, 4.049391722272796e-9, -5.842990090343335e-8, 2.0791475898344933e-8, -4.611382976227923e-8], [1.2758654963499304e-6, 7.031831733677273e-7, 2.3397568638680056e-7, -6.627632975520895e-7, -1.2971845737361725e-6, -1.1053816794111627e-6, -1.3809855172837563e-6, -2.3807360415918326e-6, -9.144182940883154e-7, 5.794805679660865e-6 … -5.086566709877729e-6, -1.9403597106517417e-6, -3.324174973179486e-8, -1.1837777291683872e-7, -5.825305855670762e-7, -3.9632569230781886e-7, -2.4574183260171652e-9, -2.5601271365636006e-9, 1.0464191828935867e-7, 5.320706425432275e-7]]), H = [1.124432851905294 0.005577947493276441 … 0.0 0.0; 0.11727946427356396 1.0199166416050651 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0], R = [1.1305324901055107 0.11135222193130294 … 0.0 0.0; 0.0 1.01439205993577 … 0.0 0.0; … ; 0.0 0.0 … 0.0 0.0; 0.0 0.0 … 0.0 0.0]))From the result of solve_ΩplusK_split we can easily compute the polarisabilities:
println("Non-interacting polarizability: $(dipole(basis, res.δρ0))")
println("Interacting polarizability: $(dipole(basis, res.δρ))")Non-interacting polarizability: 1.9257125290426884
Interacting polarizability: 1.7736548586696415As expected, the interacting polarizability matches the finite difference result. The non-interacting polarizability is higher.
Manual solution of the Dyson equations
To see what goes on under the hood, we also show how to manually solve the Dyson equation using KrylovKit:
using KrylovKit
# Apply ``(1- χ_0 K)``
function dielectric_operator(δρ)
δV = apply_kernel(basis, δρ; scfres.ρ)
χ0δV = apply_χ0(scfres, δV).δρ
δρ - χ0δV
end
# Apply ``χ_0`` once to get non-interacting dipole
δρ_nointeract = apply_χ0(scfres, δVext).δρ
# Solve Dyson equation to get interacting dipole
δρ = linsolve(dielectric_operator, δρ_nointeract, verbosity=3)[1]
println("Non-interacting polarizability: $(dipole(basis, δρ_nointeract))")
println("Interacting polarizability: $(dipole(basis, δρ))")[ Info: GMRES linsolve starts with norm of residual = 4.19e+00
[ Info: GMRES linsolve in iteration 1; step 1: normres = 2.49e-01
[ Info: GMRES linsolve in iteration 1; step 2: normres = 3.77e-03
[ Info: GMRES linsolve in iteration 1; step 3: normres = 2.85e-04
[ Info: GMRES linsolve in iteration 1; step 4: normres = 4.69e-06
[ Info: GMRES linsolve in iteration 1; step 5: normres = 1.09e-08
[ Info: GMRES linsolve in iteration 1; step 6: normres = 6.27e-11
[ Info: GMRES linsolve in iteration 1; step 7: normres = 8.39e-09
[ Info: GMRES linsolve in iteration 2; step 1: normres = 7.81e-10
[ Info: GMRES linsolve in iteration 2; step 2: normres = 3.22e-11
┌ Info: GMRES linsolve converged at iteration 2, step 3:
│ * norm of residual = 3.70e-12
└ * number of operations = 12
Non-interacting polarizability: 1.925712529042633
Interacting polarizability: 1.77365485930914We obtain the identical result to above.
- HS2025M. F. Herbst and B. Sun. Efficient Krylov methods for linear response in plane-wave electronic structure calculations. arXiv 2505.02319