Coverage for gpaw/test/test_jellium.py: 96%

28 statements  

« prev     ^ index     » next       coverage.py v7.7.1, created at 2025-07-20 00:19 +0000

1import pytest 

2import numpy as np 

3from ase import Atoms 

4from ase.units import Bohr 

5from gpaw.jellium import JelliumSlab 

6from gpaw import GPAW, Mixer 

7 

8rs = 5.0 * Bohr # Wigner-Seitz radius 

9h = 0.24 # grid-spacing 

10a = 8 * h # lattice constant 

11v = 3 * a # vacuum 

12L = 8 * a # thickness 

13k = 6 # number of k-points (k*k*1) 

14 

15ne = a**2 * L / (4 * np.pi / 3 * rs**3) 

16 

17 

18@pytest.mark.libxc 

19def test_jellium(in_tmp_dir, gpaw_new): 

20 x = h / 4 # make sure surfaces are between grid-points 

21 bc = JelliumSlab(ne, z1=v - x, z2=v + L + x) 

22 

23 surf = Atoms(pbc=(True, True, False), 

24 cell=(a, a, v + L + v)) 

25 params = dict( 

26 mode='fd', 

27 poissonsolver={'dipolelayer': 'xy'}, 

28 xc='LDA_X+LDA_C_WIGNER', 

29 eigensolver='dav', 

30 kpts=[k, k, 1], 

31 h=h, 

32 maxiter=300, 

33 convergence={'density': 1e-5}, 

34 mixer=Mixer(0.3, 7, 100), 

35 nbands=int(ne / 2) + 15) 

36 if gpaw_new: 

37 surf.calc = GPAW(**params, environment=bc) 

38 else: 

39 surf.calc = GPAW(**params, background_charge=bc) 

40 

41 _ = surf.get_potential_energy() 

42 

43 # Get the work function 

44 v_r = surf.calc.get_electrostatic_potential() 

45 efermi = surf.calc.get_fermi_level() 

46 phi = v_r[:, :, -1].mean() - efermi 

47 assert phi == pytest.approx(2.715, abs=1e-3) 

48 # Reference value: Lang and Kohn, 1971, Theory of Metal Surfaces: 

49 # Work function 

50 # DOI 10.1103/PhysRevB.3.1215 

51 # r_s = 5, work function = 2.73 eV 

52 

53 surf.calc.write('jellium.gpw')