{ "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.234211192012 -0.50 7.0 \n", " 2 -7.249858111261 -1.81 -1.39 1.0 10.7ms\n", " 3 -7.250950812143 -2.96 -1.81 2.0 13.3ms\n", " 4 -7.251035258326 -4.07 -1.92 1.0 10.7ms\n", " 5 -7.251316543840 -3.55 -2.48 1.0 11.5ms\n", " 6 -7.251337748898 -4.67 -3.24 1.0 11.3ms\n", " 7 -7.251338737887 -6.00 -3.55 2.0 14.1ms\n", " 8 -7.251338786528 -7.31 -3.97 1.0 11.6ms\n", " 9 -7.251338794942 -8.08 -4.27 2.0 13.7ms\n", " 10 -7.251338798444 -8.46 -4.82 1.0 12.2ms\n", " 11 -7.251338798652 -9.68 -5.07 2.0 91.4ms\n", " 12 -7.251338798694 -10.37 -5.49 1.0 13.2ms\n", " 13 -7.251338798703 -11.03 -6.05 2.0 14.2ms\n", " 14 -7.251338798704 -11.96 -6.41 2.0 14.9ms\n", " 15 -7.251338798705 -12.98 -6.71 1.0 14.2ms\n", " 16 -7.251338798705 -13.50 -7.30 1.0 12.6ms\n", " 17 -7.251338798705 -14.75 -7.27 3.0 17.1ms\n", " 18 -7.251338798705 -14.45 -7.84 1.0 15.1ms\n", " 19 -7.251338798705 + -14.75 -7.86 3.0 16.6ms\n", " 20 -7.251338798705 -14.75 -8.37 1.0 12.7ms\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.08456233519723523\n", "[ Info: Arnoldi iteration step 2: normres = 0.5179755840173649\n", "[ Info: Arnoldi iteration step 3: normres = 0.9022383685726841\n", "[ Info: Arnoldi iteration step 4: normres = 0.29945840114044103\n", "[ Info: Arnoldi iteration step 5: normres = 0.4327932321046853\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 5: 0 values converged, normres = (3.27e-02, 5.34e-02, 3.05e-01, 2.97e-01, 4.61e-02)\n", "[ Info: Arnoldi iteration step 6: normres = 0.2799077910190465\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 6: 0 values converged, normres = (7.71e-03, 6.63e-02, 2.52e-01, 7.51e-02, 2.96e-02)\n", "[ Info: Arnoldi iteration step 7: normres = 0.0684593419020992\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 7: 0 values converged, normres = (2.54e-04, 5.23e-03, 1.46e-02, 2.18e-02, 5.73e-02)\n", "[ Info: Arnoldi iteration step 8: normres = 0.10264934506905754\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 8: 0 values converged, normres = (1.10e-05, 3.67e-04, 1.13e-03, 6.13e-03, 3.45e-02)\n", "[ Info: Arnoldi iteration step 9: normres = 0.047143873548756675\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 9: 0 values converged, normres = (2.25e-07, 1.24e-05, 4.26e-05, 1.00e-03, 1.89e-02)\n", "[ Info: Arnoldi iteration step 10: normres = 0.09969496360691892\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 10: 0 values converged, normres = (9.76e-09, 8.87e-07, 3.38e-06, 3.33e-04, 2.05e-02)\n", "[ Info: Arnoldi iteration step 11: normres = 0.09886544097188818\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 11: 0 values converged, normres = (4.16e-10, 6.22e-08, 2.63e-07, 1.10e-04, 2.77e-02)\n", "[ Info: Arnoldi iteration step 12: normres = 0.08765763194605695\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 12: 0 values converged, normres = (1.62e-11, 4.05e-09, 1.91e-08, 4.17e-05, 4.83e-02)\n", "[ Info: Arnoldi iteration step 13: normres = 0.13655294900789272\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 13: 1 values converged, normres = (9.49e-13, 3.90e-10, 2.04e-09, 1.86e-05, 4.10e-02)\n", "[ Info: Arnoldi iteration step 14: normres = 0.3869796553350125\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 14: 1 values converged, normres = (3.60e-13, 1.40e-09, 3.84e-01, 1.74e-03, 2.38e-03)\n", "[ Info: Arnoldi iteration step 15: normres = 0.1486236732026821\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 15: 1 values converged, normres = (2.44e-14, 1.20e-09, 3.92e-02, 1.68e-04, 1.18e-06)\n", "[ Info: Arnoldi iteration step 16: normres = 0.1462347982046594\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 16: 1 values converged, normres = (3.40e-15, 2.22e-08, 3.98e-02, 2.79e-04, 1.38e-01)\n", "[ Info: Arnoldi iteration step 17: normres = 0.02222917509155539\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 17: 1 values converged, normres = (3.14e-17, 3.96e-04, 4.75e-04, 2.05e-03, 8.14e-04)\n", "[ Info: Arnoldi iteration step 18: normres = 0.03861433315981629\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 18: 1 values converged, normres = (5.03e-19, 3.41e-08, 1.59e-05, 1.85e-08, 6.26e-05)\n", "[ Info: Arnoldi iteration step 19: normres = 0.11390171382164471\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 19: 1 values converged, normres = (2.39e-20, 8.51e-07, 8.64e-07, 4.83e-06, 2.10e-06)\n", "[ Info: Arnoldi iteration step 20: normres = 0.06766662926495326\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 20: 1 values converged, normres = (7.96e-22, 2.02e-08, 7.04e-08, 3.24e-07, 1.63e-07)\n", "[ Info: Arnoldi iteration step 21: normres = 0.023676060286042317\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 21: 1 values converged, normres = (7.79e-24, 4.42e-10, 1.06e-09, 6.28e-09, 6.53e-12)\n", "[ Info: Arnoldi iteration step 22: normres = 0.049009561202773874\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 22: 1 values converged, normres = (1.58e-25, 1.43e-11, 3.46e-11, 2.22e-10, 4.07e-11)\n", "[ Info: Arnoldi iteration step 23: normres = 0.6363803881591422\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 23: 1 values converged, normres = (6.80e-26, 1.59e-11, 3.82e-11, 3.24e-10, 5.91e-11)\n", "[ Info: Arnoldi iteration step 24: normres = 0.019534648524535653\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 24: 2 values converged, normres = (7.76e-28, 8.44e-13, 2.03e-12, 1.01e-02, 3.76e-03)\n", "[ Info: Arnoldi iteration step 25: normres = 0.06999266782733406\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 25: 3 values converged, normres = (2.24e-29, 3.89e-14, 9.37e-14, 1.07e-08, 6.81e-10)\n", "[ Info: Arnoldi iteration step 26: normres = 0.04679167682916041\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 26: 3 values converged, normres = (4.73e-31, 1.39e-15, 3.34e-15, 5.50e-09, 1.68e-09)\n", "[ Info: Arnoldi iteration step 27: normres = 0.027873412254176694\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 27: 3 values converged, normres = (5.42e-33, 2.54e-17, 6.12e-17, 1.43e-07, 3.20e-07)\n", "[ Info: Arnoldi iteration step 28: normres = 0.09382633291733547\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 28: 3 values converged, normres = (2.23e-34, 1.73e-18, 4.17e-18, 1.42e-08, 2.65e-08)\n", "[ Info: Arnoldi iteration step 29: normres = 0.07578667197829558\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 29: 3 values converged, normres = (7.12e-36, 8.95e-20, 2.16e-19, 1.15e-09, 1.54e-09)\n", "[ Info: Arnoldi iteration step 30: normres = 0.06832169931206256\n", "[ Info: Arnoldi schursolve in iter 1, krylovdim = 30: 3 values converged, normres = (2.40e-37, 5.46e-21, 1.32e-20, 8.68e-11, 1.02e-10)\n", "[ Info: Arnoldi schursolve in iter 2, krylovdim = 19: 3 values converged, normres = (2.41e-37, 5.46e-21, 1.32e-20, 8.68e-11, 1.02e-10)\n", "[ Info: Arnoldi iteration step 20: normres = 0.06582324653140685\n", "[ Info: Arnoldi schursolve in iter 2, krylovdim = 20: 3 values converged, normres = (6.62e-39, 2.43e-22, 5.86e-22, 4.21e-12, 5.08e-12)\n", "[ Info: Arnoldi iteration step 21: normres = 0.03653804150881428\n", "┌ Info: Arnoldi eigsolve finished after 2 iterations:\n", "│ * 6 eigenvalues converged\n", "│ * norm of residuals = (1.040582233206911e-40, 6.2698064102266406e-24, 1.6339024302874763e-23, 1.2050992547309896e-13, 1.4396409442930286e-13, 2.571462452083856e-14)\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 }