1import pytest 2import numpy as np 3from ase.units import Bohr 4from ase.build import bulk 5from gpaw import GPAW, FermiDirac 6from gpaw.response.df import DielectricFunction 7from gpaw.test import equal, findpeak 8 9 10@pytest.mark.libxc 11def test_response_diamond_absorption(in_tmp_dir): 12 a = 6.75 * Bohr 13 atoms = bulk('C', 'diamond', a=a) 14 15 calc = GPAW(mode='pw', 16 kpts=(3, 3, 3), 17 eigensolver='rmm-diis', 18 occupations=FermiDirac(0.001)) 19 20 atoms.calc = calc 21 atoms.get_potential_energy() 22 calc.write('C.gpw', 'all') 23 24 eM1_ = 9.727 25 eM2_ = 9.548 26 w0_ = 10.7782 27 I0_ = 5.47 28 w_ = 10.7532 29 I_ = 5.98 30 31 # Macroscopic dielectric constant calculation 32 df = DielectricFunction('C.gpw', frequencies=(0.,), eta=0.001, ecut=50, 33 hilbert=False) 34 eM1, eM2 = df.get_macroscopic_dielectric_constant() 35 equal(eM1, eM1_, 0.01) 36 equal(eM2, eM2_, 0.01) 37 38 # Absorption spectrum calculation RPA 39 df = DielectricFunction('C.gpw', 40 eta=0.25, 41 ecut=50, 42 frequencies=np.linspace(0, 24., 241), 43 hilbert=False) 44 a0, a = df.get_dielectric_function(filename=None) 45 df.check_sum_rule(a.imag) 46 47 equal(a0[0].real, eM1_, 0.01) 48 equal(a[0].real, eM2_, 0.01) 49 w, I = findpeak(np.linspace(0, 24., 241), a0.imag) 50 equal(w, w0_, 0.01) 51 equal(I / (4 * np.pi), I0_, 0.1) 52 w, I = findpeak(np.linspace(0, 24., 241), a.imag) 53 equal(w, w_, 0.01) 54 equal(I / (4 * np.pi), I_, 0.1) 55 56 a0, a = df.get_polarizability(filename=None) 57 58 w, I = findpeak(np.linspace(0, 24., 241), a0.imag) 59 equal(w, w0_, 0.01) 60 equal(I, I0_, 0.01) 61 w, I = findpeak(np.linspace(0, 24., 241), a.imag) 62 equal(w, w_, 0.01) 63 equal(I, I_, 0.01) 64 65 # Absorption spectrum calculation ALDA 66 w0_ = 10.7931 67 I0_ = 5.36 68 w_ = 10.7562 69 I_ = 5.8803 70 71 a0, a = df.get_polarizability(filename=None, xc='ALDA') 72 73 w, I = findpeak(np.linspace(0, 24., 241), a0.imag) 74 equal(w, w0_, 0.01) 75 equal(I, I0_, 0.1) 76 w, I = findpeak(np.linspace(0, 24., 241), a.imag) 77 equal(w, w_, 0.01) 78 equal(I, I_, 0.1) 79 80 # Absorption spectrum calculation long-range kernel 81 w0_ = 10.2189 82 I0_ = 5.14 83 w_ = 10.2906 84 I_ = 5.6955 85 86 a0, a = df.get_polarizability(filename=None, xc='LR0.25') 87 88 w, I = findpeak(np.linspace(0, 24., 241), a0.imag) 89 equal(w, w0_, 0.01) 90 equal(I, I0_, 0.1) 91 w, I = findpeak(np.linspace(0, 24., 241), a.imag) 92 equal(w, w_, 0.01) 93 equal(I, I_, 0.1) 94 95 # Absorption spectrum calculation Bootstrap 96 w0_ = 10.37 97 I0_ = 5.27 98 w_ = 10.4600 99 I_ = 6.0263 100 101 a0, a = df.get_polarizability(filename=None, xc='Bootstrap') 102 103 w, I = findpeak(np.linspace(0, 24., 241), a0.imag) 104 equal(w, w0_, 0.02) 105 equal(I, I0_, 0.2) 106 w, I = findpeak(np.linspace(0, 24., 241), a.imag) 107 equal(w, w_, 0.02) 108 equal(I, I_, 0.2) 109 110 # Absorption spectrum calculation RPA - Wigner-Seitz truncation 111 w0_ = 10.7780 112 I0_ = 5.5472 113 w_ = 10.7532 114 I_ = 6.0750 115 116 df_ws = DielectricFunction('C.gpw', 117 eta=0.25, 118 ecut=50, 119 frequencies=np.linspace(0, 24., 241), 120 hilbert=False, 121 truncation='wigner-seitz') 122 123 a0_ws, a_ws = df_ws.get_polarizability(filename=None) 124 125 w, I = findpeak(np.linspace(0, 24., 241), a0_ws.imag) 126 equal(w, w0_, 0.02) 127 equal(I, I0_, 0.2) 128 w, I = findpeak(np.linspace(0, 24., 241), a_ws.imag) 129 equal(w, w_, 0.1) 130 equal(I, I_, 0.2) 131 # The Wigner-Seitz truncation does not give exactly the 132 # same for kpts=(3,3,3) 133