1import pytest
2from gpaw.mpi import world
3from gpaw.utilities import compiled_with_sl
4import numpy as np
5from ase.build import bulk
6from gpaw import GPAW, FermiDirac
7from gpaw.response.bse import BSE
8from gpaw.test import findpeak, equal
9
10pytestmark = pytest.mark.skipif(
11    world.size != 4 or not compiled_with_sl(),
12    reason='world.size != 4 or not compiled_with_sl()')
13
14
15def test_response_bse_silicon(in_tmp_dir):
16    GS = 1
17    nosym = 1
18    bse = 1
19    check = 1
20
21    if GS:
22        a = 5.431  # From PRB 73,045112 (2006)
23        atoms = bulk('Si', 'diamond', a=a)
24        atoms.positions -= a / 8
25        calc = GPAW(mode='pw',
26                    kpts={'size': (2, 2, 2), 'gamma': True},
27                    occupations=FermiDirac(0.001),
28                    nbands=12,
29                    convergence={'bands': -4})
30        atoms.calc = calc
31        atoms.get_potential_energy()
32        calc.write('Si.gpw', 'all')
33
34    if bse:
35        eshift = 0.8
36        bse = BSE('Si.gpw',
37                  ecut=50.,
38                  valence_bands=range(4),
39                  conduction_bands=range(4, 8),
40                  eshift=eshift,
41                  nbands=8,
42                  write_h=False,
43                  write_v=False)
44        w_w, eps_w = bse.get_dielectric_function(filename=None,
45                                                 eta=0.2,
46                                                 w_w=np.linspace(0, 10, 2001))
47    if check:
48        w_ = 2.552
49        I_ = 421.15
50        w, I = findpeak(w_w, eps_w.imag)
51        equal(w, w_, 0.01)
52        equal(I, I_, 0.1)
53
54    if GS and nosym:
55        atoms = bulk('Si', 'diamond', a=a)
56        calc = GPAW(mode='pw',
57                    kpts={'size': (2, 2, 2), 'gamma': True},
58                    occupations=FermiDirac(0.001),
59                    nbands=12,
60                    symmetry='off',
61                    convergence={'bands': -4})
62        atoms.calc = calc
63        atoms.get_potential_energy()
64        calc.write('Si.gpw', 'all')
65
66    if bse and nosym:
67        bse = BSE('Si.gpw',
68                  ecut=50.,
69                  valence_bands=range(4),
70                  conduction_bands=range(4, 8),
71                  eshift=eshift,
72                  nbands=8,
73                  write_h=False,
74                  write_v=False)
75        w_w, eps_w = bse.get_dielectric_function(filename=None,
76                                                 eta=0.2,
77                                                 w_w=np.linspace(0, 10, 2001))
78
79    if check and nosym:
80        w, I = findpeak(w_w, eps_w.imag)
81        equal(w, w_, 0.01)
82        equal(I, I_, 0.1)
83