Coverage for gpaw/test/gllb/test_ne_disc.py: 100%

43 statements  

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

1import pytest 

2from ase import Atom, Atoms 

3 

4from gpaw import GPAW, restart, Davidson, Mixer 

5from gpaw.atom.generator import Generator 

6from gpaw.atom.configurations import parameters 

7from gpaw.mpi import world 

8 

9# This test calculates the derivative discontinuity of Ne-atom 

10# first on 3D without restart. Then does restart and recalculates. 

11 

12 

13@pytest.mark.slow 

14@pytest.mark.gllb 

15@pytest.mark.libxc 

16def test_gllb_ne_disc(in_tmp_dir, add_cwd_to_setup_paths): 

17 atom = 'Ne' 

18 

19 for xcname in ['GLLB', 'GLLBSC']: 

20 if world.rank == 0: 

21 g = Generator(atom, xcname=xcname, scalarrel=False, nofiles=True) 

22 g.run(**parameters[atom]) 

23 eps = g.e_j[-1] 

24 world.barrier() 

25 

26 a = 10 

27 Ne = Atoms([Atom(atom, (0, 0, 0))], 

28 cell=(a, a, a), pbc=False) 

29 Ne.center() 

30 calc = GPAW(mode='fd', 

31 eigensolver=Davidson(4), 

32 nbands=10, 

33 h=0.18, 

34 xc=xcname, 

35 basis='dzp', 

36 mixer=Mixer(0.6)) 

37 Ne.calc = calc 

38 Ne.get_potential_energy() 

39 homo, lumo = calc.get_homo_lumo() 

40 response = calc.hamiltonian.xc.response 

41 dxc_pot = response.calculate_discontinuity_potential(homo, lumo) 

42 KS, dxc = response.calculate_discontinuity(dxc_pot) 

43 if xcname == 'GLLB': 

44 assert KS + dxc == pytest.approx(24.89, abs=1.5e-1) 

45 else: 

46 assert KS + dxc == pytest.approx(27.70, abs=6.0e-2) 

47 eps3d = calc.wfs.kpt_u[0].eps_n[3] 

48 calc.write('Ne_temp.gpw', mode='all') 

49 

50 atoms, calc = restart('Ne_temp.gpw') 

51 homo, lumo = calc.get_homo_lumo() 

52 response = calc.hamiltonian.xc.response 

53 dxc_pot = response.calculate_discontinuity_potential(homo, lumo) 

54 KS2, dxc2 = response.calculate_discontinuity(dxc_pot) 

55 assert KS == pytest.approx(KS2, abs=1e-5) 

56 assert dxc2 == pytest.approx(dxc, abs=1e-5) 

57 

58 # Hardness of Ne 24.71eV by GLLB+Dxc, experimental I-A = I = 21.56eV 

59 # 

60 # Not sure where 24.71 comes from, but with better grid and better 

61 # stencil, result becomes 24.89. --askhl 

62 

63 if world.rank == 0: 

64 assert eps == pytest.approx(eps3d, abs=1e-3) 

65 if xcname == 'GLLB': 

66 assert 24.89 == pytest.approx(KS2 + dxc2, abs=1.2e-1)