Coverage for gpaw/test/noncollinear/test_kinetic_energy.py: 100%

1import numpy as np 

2import pytest 

3 

4from gpaw.mpi import world 

5from gpaw.new.ase_interface import GPAW 

6from gpaw.utilities import pack_density 

7 

8 

9@pytest.mark.soc 

10@pytest.mark.skipif(world.size > 1, reason='Gamma-point calculation.') 

11def test_kinetic_energy(gpw_files): 

12 # Test that we calculate the kinetic energy correctly in noncollinear mode. 

13 # We are mostly concerned with the contributions from the complex parts of 

14 # the atomic density matrix and Hamiltonian (spin-orbit). 

15 

16 # Thallium atom in a simulation box, we use a heavy atom with large SOC. 

17 calc = GPAW(gpw_files['Tl_box_pw']) 

18 setup = calc.dft.setups[0] 

19 

20 # Kinetic energy calculated from sum of bands (the standard way): 

21 Ekin1 = (calc.dft.energies._energies['band'] + 

22 calc.dft.energies._energies['kinetic_correction']) 

23 

24 # Kinetic energy calculated directly from second-order derivative: 

25 wfs = calc.dft.ibzwfs.wfs_qs[0][0] 

26 

27 psit = wfs.psit_nX 

28 occ_n = wfs.occ_n 

29 

30 psit_nsG = psit.data 

31 G_plus_k_Gv = psit.desc.G_plus_k_Gv 

32 ucvol = psit.desc.volume 

33 

34 laplacian_on_psit_nsv = np.abs(psit_nsG)**2 @ G_plus_k_Gv**2 

35 laplacian_on_psit_n = - np.sum(np.sum(laplacian_on_psit_nsv, axis=1), 

36 axis=1) * ucvol 

37 Ekin_pseudo = -0.5 * (laplacian_on_psit_n @ occ_n) 

38 

39 D_p = pack_density(calc.dft.density.D_asii[0][0].real) 

40 Ekin_PAW = setup.K_p @ D_p + setup.Kc 

41 

42 Ekin2 = Ekin_pseudo + Ekin_PAW 

43 

44 # We should get the same value for the kinetic 

45 # energy irrespective of the method used. 

46 assert Ekin1 == pytest.approx(Ekin2, abs=1e-5)