{ "cells": [ { "cell_type": "markdown", "source": [ "# Eigenvalues of the dielectric matrix\n", "\n", "We compute a few eigenvalues of the dielectric matrix ($q=0$, $ω=0$) iteratively." ], "metadata": {} }, { "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "n Energy log10(ΔE) log10(Δρ) Diag Δtime\n", "--- --------------- --------- --------- ---- ------\n", " 1 -7.232797979692 -0.50 6.0 \n", " 2 -7.249367036029 -1.78 -1.39 1.0 10.0ms\n", " 3 -7.250972694135 -2.79 -1.82 2.0 11.6ms\n", " 4 -7.250953560013 + -4.72 -1.87 2.0 12.1ms\n", " 5 -7.251315577644 -3.44 -2.48 1.0 10.1ms\n", " 6 -7.251337763936 -4.65 -3.06 1.0 10.0ms\n", " 7 -7.251338657416 -6.05 -3.59 1.0 10.0ms\n", " 8 -7.251338783761 -6.90 -3.89 2.0 11.8ms\n", " 9 -7.251338796349 -7.90 -4.32 1.0 10.4ms\n", " 10 -7.251338798482 -8.67 -4.97 2.0 11.9ms\n", " 11 -7.251338798634 -9.82 -5.05 3.0 13.8ms\n", " 12 -7.251338798692 -10.23 -5.55 1.0 10.5ms\n", " 13 -7.251338798698 -11.23 -5.60 3.0 14.0ms\n", " 14 -7.251338798704 -11.23 -6.22 1.0 10.5ms\n", " 15 -7.251338798705 -12.16 -6.98 3.0 14.0ms\n", " 16 -7.251338798705 -13.91 -7.16 2.0 12.0ms\n", " 17 -7.251338798705 -14.45 -7.44 1.0 10.3ms\n", " 18 -7.251338798705 + -Inf -8.05 2.0 12.3ms\n" ] } ], "cell_type": "code", "source": [ "using DFTK\n", "using Plots\n", "using KrylovKit\n", "using Printf\n", "\n", "# Calculation parameters\n", "kgrid = [1, 1, 1]\n", "Ecut = 5\n", "\n", "# Silicon lattice\n", "a = 10.26\n", "lattice = a / 2 .* [[0 1 1.]; [1 0 1.]; [1 1 0.]]\n", "Si = ElementPsp(:Si, psp=load_psp(\"hgh/lda/Si-q4\"))\n", "atoms = [Si, Si]\n", "positions = [ones(3)/8, -ones(3)/8]\n", "\n", "# Compute the dielectric operator without symmetries\n", "model = model_LDA(lattice, atoms, positions, symmetries=false)\n", "basis = PlaneWaveBasis(model; Ecut, kgrid)\n", "scfres = self_consistent_field(basis, tol=1e-8);" ], "metadata": {}, "execution_count": 1 }, { "cell_type": "markdown", "source": [ "Applying $ε^† ≔ (1- χ_0 K)$ …" ], "metadata": {} }, { "outputs": [], "cell_type": "code", "source": [ "function eps_fun(δρ)\n", " δV = apply_kernel(basis, δρ; ρ=scfres.ρ)\n", " χ0δV = apply_χ0(scfres, δV)\n", " δρ - χ0δV\n", "end;" ], "metadata": {}, "execution_count": 2 }, { "cell_type": "markdown", "source": [ "… eagerly diagonalizes the subspace matrix at each iteration" ], "metadata": {} }, { "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "[ Info: Arnoldi iteration step 1: normres = 0.06198435834205509\n", "[ Info: Arnoldi iteration step 2: normres = 0.5158645562562821\n", "[ Info: Arnoldi iteration step 3: normres = 0.8459982776164054\n", "[ Info: Arnoldi iteration step 4: normres = 0.19587794323692703\n", "[ Info: Arnoldi iteration step 5: normres = 0.40686045264409537\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 5: 0 values converged, normres = (1.47e-02, 3.02e-02, 3.49e-01, 2.06e-01, 1.02e-02)\n", "[ Info: Arnoldi iteration step 6: normres = 0.37411436293354805\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 6: 0 values converged, normres = (4.99e-03, 2.61e-01, 2.25e-01, 1.03e-01, 3.41e-02)\n", "[ Info: Arnoldi iteration step 7: normres = 0.09537329598050488\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 7: 0 values converged, normres = (2.53e-04, 2.84e-02, 6.29e-03, 6.85e-02, 5.49e-02)\n", "[ Info: Arnoldi iteration step 8: normres = 0.12263535557645752\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 8: 0 values converged, normres = (1.33e-05, 2.44e-03, 6.13e-04, 2.49e-02, 3.67e-02)\n", "[ Info: Arnoldi iteration step 9: normres = 0.07015692665996354\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 9: 0 values converged, normres = (4.11e-07, 1.26e-04, 3.53e-05, 7.16e-03, 3.61e-02)\n", "[ Info: Arnoldi iteration step 10: normres = 0.09549714432971323\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 10: 0 values converged, normres = (1.72e-08, 8.79e-06, 2.75e-06, 2.61e-03, 3.31e-02)\n", "[ Info: Arnoldi iteration step 11: normres = 0.07107279716594946\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 11: 0 values converged, normres = (5.24e-10, 4.37e-07, 1.52e-07, 5.94e-04, 1.60e-02)\n", "[ Info: Arnoldi iteration step 12: normres = 0.066254625964104\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 12: 0 values converged, normres = (1.47e-11, 2.00e-08, 7.67e-09, 1.11e-04, 5.40e-03)\n", "[ Info: Arnoldi iteration step 13: normres = 0.060141829679643835\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 13: 1 values converged, normres = (3.78e-13, 8.35e-10, 3.55e-10, 1.96e-05, 1.75e-03)\n", "[ Info: Arnoldi iteration step 14: normres = 0.6752175667957129\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 14: 1 values converged, normres = (1.61e-13, 8.30e-10, 4.44e-10, 6.64e-01, 4.36e-02)\n", "[ Info: Arnoldi iteration step 15: normres = 0.06571931169706294\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 15: 1 values converged, normres = (6.68e-15, 1.89e-10, 4.10e-02, 1.29e-03, 3.63e-05)\n", "[ Info: Arnoldi iteration step 16: normres = 0.6601214999636914\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 16: 1 values converged, normres = (3.44e-15, 3.37e-10, 1.07e-01, 1.12e-02, 6.49e-01)\n", "[ Info: Arnoldi iteration step 17: normres = 0.03540273724803503\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 17: 1 values converged, normres = (7.09e-17, 6.45e-09, 1.72e-02, 7.25e-06, 5.38e-03)\n", "[ Info: Arnoldi iteration step 18: normres = 0.01599247542469627\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 18: 1 values converged, normres = (4.70e-19, 1.81e-08, 1.82e-04, 3.23e-05, 5.21e-05)\n", "[ Info: Arnoldi iteration step 19: normres = 0.19244883663341256\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 19: 1 values converged, normres = (3.89e-20, 9.86e-06, 2.27e-05, 9.77e-09, 9.25e-06)\n", "[ Info: Arnoldi iteration step 20: normres = 0.05594731736472959\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 20: 1 values converged, normres = (1.04e-21, 9.71e-09, 1.17e-06, 3.59e-07, 3.43e-07)\n", "[ Info: Arnoldi iteration step 21: normres = 0.028143217504210106\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 21: 1 values converged, normres = (1.22e-23, 7.30e-09, 2.10e-08, 9.18e-10, 1.04e-08)\n", "[ Info: Arnoldi iteration step 22: normres = 0.02605433776547058\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 22: 1 values converged, normres = (1.31e-25, 8.02e-11, 3.73e-10, 2.36e-11, 1.96e-10)\n", "[ Info: Arnoldi iteration step 23: normres = 0.6633721844226979\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 23: 1 values converged, normres = (5.21e-26, 7.00e-11, 3.25e-10, 2.50e-11, 2.08e-10)\n", "[ Info: Arnoldi iteration step 24: normres = 0.01921959552375744\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 24: 1 values converged, normres = (6.53e-28, 4.74e-12, 2.20e-11, 2.46e-10, 2.06e-09)\n", "[ Info: Arnoldi iteration step 25: normres = 0.03686838831214523\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 25: 3 values converged, normres = (9.90e-30, 1.15e-13, 5.34e-13, 2.98e-05, 2.37e-04)\n", "[ Info: Arnoldi iteration step 26: normres = 0.06820801126164587\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 26: 3 values converged, normres = (3.05e-31, 5.99e-15, 2.78e-14, 1.35e-05, 1.02e-05)\n", "[ Info: Arnoldi iteration step 27: normres = 0.03017217897859313\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 27: 3 values converged, normres = (3.81e-33, 1.20e-16, 5.58e-16, 2.68e-09, 2.18e-08)\n", "[ Info: Arnoldi iteration step 28: normres = 0.09352966133290364\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 28: 3 values converged, normres = (1.55e-34, 8.09e-18, 3.76e-17, 4.31e-10, 1.69e-08)\n", "[ Info: Arnoldi iteration step 29: normres = 0.03653691232545321\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 29: 3 values converged, normres = (2.38e-36, 2.00e-19, 9.29e-19, 7.52e-12, 8.57e-10)\n", "[ Info: Arnoldi iteration step 30: normres = 0.19178918475219967\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 30: 3 values converged, normres = (2.14e-37, 3.15e-20, 1.46e-19, 1.32e-12, 1.52e-10)\n", "[ Info: Arnoldi schursolve in iter 2, krylovdim = 19: 3 values converged, normres = (2.14e-37, 3.15e-20, 1.46e-19, 1.32e-12, 1.52e-10)\n", "[ Info: Arnoldi iteration step 20: normres = 0.04558371531091134\n", "[ Info: Arnoldi schursolve in iter 2, krylovdim = 20: 4 values converged, normres = (4.26e-39, 1.04e-21, 4.84e-21, 4.89e-14, 5.61e-12)\n", "[ Info: Arnoldi iteration step 21: normres = 0.07873596414034613\n", "┌ Info: Arnoldi eigsolve finished after 2 iterations:\n", "│ * 6 eigenvalues converged\n", "│ * norm of residuals = (1.428561320805377e-40, 5.711665972977044e-23, 2.63581390873658e-22, 2.9658733101081714e-15, 3.404235781471831e-13, 1.6997031325785676e-13)\n", "└ * number of operations = 32\n" ] } ], "cell_type": "code", "source": [ "eigsolve(eps_fun, randn(size(scfres.ρ)), 5, :LM; eager=true, verbosity=3);" ], "metadata": {}, "execution_count": 3 } ], "nbformat_minor": 3, "metadata": { "language_info": { "file_extension": ".jl", "mimetype": "application/julia", "name": "julia", "version": "1.8.5" }, "kernelspec": { "name": "julia-1.8", "display_name": "Julia 1.8.5", "language": "julia" } }, "nbformat": 4 }