{ "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\n", "--- --------------- --------- --------- ----\n", " 1 -7.234567989945 -0.50 8.0\n", " 2 -7.250283641363 -1.80 -1.41 1.0\n", " 3 -7.251049665235 -3.12 -1.99 1.0\n", " 4 -7.250969370305 + -4.10 -1.88 2.0\n", " 5 -7.251335740565 -3.44 -2.77 1.0\n", " 6 -7.251338430435 -5.57 -3.40 1.0\n", " 7 -7.251338723566 -6.53 -3.57 2.0\n", " 8 -7.251338791097 -7.17 -4.05 1.0\n", " 9 -7.251338797581 -8.19 -4.67 1.0\n", " 10 -7.251338798545 -9.02 -4.85 3.0\n", " 11 -7.251338798659 -9.95 -5.18 1.0\n", " 12 -7.251338798688 -10.53 -5.39 2.0\n", " 13 -7.251338798702 -10.85 -5.99 1.0\n", " 14 -7.251338798704 -11.71 -6.30 3.0\n", " 15 -7.251338798705 -12.51 -6.99 2.0\n", " 16 -7.251338798705 -14.27 -7.55 1.0\n", " 17 -7.251338798705 -15.05 -7.79 2.0\n", " 18 -7.251338798705 -14.45 -8.18 2.0\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.034184863026181254\n", "[ Info: Arnoldi iteration step 2: normres = 0.5319725171416453\n", "[ Info: Arnoldi iteration step 3: normres = 0.2471089116383626\n", "[ Info: Arnoldi iteration step 4: normres = 0.9359558122835486\n", "[ Info: Arnoldi iteration step 5: normres = 0.5808161881623259\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 5: 0 values converged, normres = (4.80e-01, 8.47e-02, 2.55e-01, 1.81e-01, 4.15e-02)\n", "[ Info: Arnoldi iteration step 6: normres = 0.4648591697894539\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 6: 0 values converged, normres = (2.06e-01, 8.56e-02, 3.86e-01, 1.03e-01, 4.19e-02)\n", "[ Info: Arnoldi iteration step 7: normres = 0.09346727050181008\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 7: 0 values converged, normres = (1.05e-02, 1.46e-02, 3.14e-02, 4.53e-02, 6.35e-02)\n", "[ Info: Arnoldi iteration step 8: normres = 0.10816648518912077\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 8: 0 values converged, normres = (4.82e-04, 1.08e-03, 2.58e-03, 1.37e-02, 3.79e-02)\n", "[ Info: Arnoldi iteration step 9: normres = 0.06660828421859441\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 9: 0 values converged, normres = (1.41e-05, 5.29e-05, 1.41e-04, 3.59e-03, 3.52e-02)\n", "[ Info: Arnoldi iteration step 10: normres = 0.09780626001875309\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 10: 0 values converged, normres = (6.10e-07, 3.81e-06, 1.13e-05, 1.39e-03, 4.00e-02)\n", "[ Info: Arnoldi iteration step 11: normres = 0.08491258147977278\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 11: 0 values converged, normres = (2.25e-08, 2.32e-07, 7.65e-07, 4.19e-04, 2.97e-02)\n", "[ Info: Arnoldi iteration step 12: normres = 0.060500297194820876\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 12: 0 values converged, normres = (5.80e-10, 9.77e-09, 3.57e-08, 7.70e-05, 1.04e-02)\n", "[ Info: Arnoldi iteration step 13: normres = 0.04926179791327646\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 13: 0 values converged, normres = (1.22e-11, 3.35e-10, 1.36e-09, 1.12e-05, 2.85e-03)\n", "[ Info: Arnoldi iteration step 14: normres = 0.6045664878144748\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 14: 0 values converged, normres = (3.68e-12, 1.85e-10, 8.59e-10, 1.47e-03, 6.00e-01)\n", "[ Info: Arnoldi iteration step 15: normres = 0.043419063951638845\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 15: 1 values converged, normres = (1.26e-13, 4.18e-11, 3.09e-02, 1.69e-02, 6.65e-06)\n", "[ Info: Arnoldi iteration step 16: normres = 0.6589960751281293\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 16: 1 values converged, normres = (4.20e-14, 2.60e-11, 2.44e-02, 4.68e-03, 6.23e-01)\n", "[ Info: Arnoldi iteration step 17: normres = 0.03601224861167034\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 17: 1 values converged, normres = (1.35e-15, 1.47e-09, 2.90e-02, 2.93e-03, 3.04e-03)\n", "[ Info: Arnoldi iteration step 18: normres = 0.040189387203134024\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 18: 1 values converged, normres = (2.24e-17, 6.30e-04, 4.50e-04, 1.42e-07, 1.27e-04)\n", "[ Info: Arnoldi iteration step 19: normres = 0.20102482168660496\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 19: 1 values converged, normres = (1.99e-18, 8.12e-10, 1.14e-04, 7.23e-06, 1.97e-05)\n", "[ Info: Arnoldi iteration step 20: normres = 0.030932288064949163\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 20: 1 values converged, normres = (2.88e-20, 2.71e-06, 1.06e-06, 2.66e-08, 6.05e-07)\n", "[ Info: Arnoldi iteration step 21: normres = 0.027080167763569436\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 21: 1 values converged, normres = (3.22e-22, 6.20e-10, 5.18e-08, 3.25e-09, 1.14e-08)\n", "[ Info: Arnoldi iteration step 22: normres = 0.0300385123945135\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 22: 1 values converged, normres = (4.02e-24, 2.83e-11, 1.04e-09, 7.13e-11, 2.53e-10)\n", "[ Info: Arnoldi iteration step 23: normres = 0.6262110934723216\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 23: 1 values converged, normres = (1.26e-24, 1.65e-11, 5.93e-10, 4.65e-11, 1.65e-10)\n", "[ Info: Arnoldi iteration step 24: normres = 0.028318595517572452\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 24: 1 values converged, normres = (2.81e-26, 2.39e-12, 8.71e-11, 1.14e-03, 5.96e-04)\n", "[ Info: Arnoldi iteration step 25: normres = 0.08872559254406792\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 25: 2 values converged, normres = (1.04e-27, 1.43e-13, 5.16e-12, 5.89e-09, 1.48e-08)\n", "[ Info: Arnoldi iteration step 26: normres = 0.03678896700914451\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 26: 3 values converged, normres = (1.71e-29, 3.96e-15, 1.43e-13, 7.07e-06, 3.50e-06)\n", "[ Info: Arnoldi iteration step 27: normres = 0.026607909859730012\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 27: 3 values converged, normres = (1.87e-31, 6.91e-17, 2.50e-15, 5.34e-10, 1.93e-09)\n", "[ Info: Arnoldi iteration step 28: normres = 0.09400716675525227\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 28: 3 values converged, normres = (7.63e-33, 4.65e-18, 1.68e-16, 5.87e-10, 6.68e-09)\n", "[ Info: Arnoldi iteration step 29: normres = 0.01555190737422291\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 29: 3 values converged, normres = (4.99e-35, 4.93e-20, 1.78e-18, 1.87e-11, 8.01e-10)\n", "[ Info: Arnoldi iteration step 30: normres = 0.037211996714061585\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 30: 4 values converged, normres = (7.65e-37, 1.20e-21, 4.36e-20, 5.02e-13, 2.15e-11)\n", "[ Info: Arnoldi schursolve in iter 2, krylovdim = 19: 4 values converged, normres = (7.65e-37, 1.20e-21, 4.36e-20, 5.02e-13, 2.15e-11)\n", "[ Info: Arnoldi iteration step 20: normres = 0.06764629636847294\n", "[ Info: Arnoldi schursolve in iter 2, krylovdim = 20: 4 values converged, normres = (2.24e-38, 5.81e-23, 2.10e-21, 2.69e-14, 1.15e-12)\n", "[ Info: Arnoldi iteration step 21: normres = 0.1672286709018542\n", "┌ Info: Arnoldi eigsolve finished after 2 iterations:\n", "│ * 7 eigenvalues converged\n", "│ * norm of residuals = (1.5978801693399404e-39, 6.757836523065391e-24, 2.380094231807992e-22, 3.4681782943249763e-15, 1.4866012021222274e-13, 2.1485774037847807e-14, 7.149966498805555e-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.4" }, "kernelspec": { "name": "julia-1.8", "display_name": "Julia 1.8.4", "language": "julia" } }, "nbformat": 4 }