Coverage for gpaw/test/hubbard/test_hubbard_U_noncollinear.py: 100%

1import pytest 

2 

3import numpy as np 

4from ase.units import Ha 

5 

6from gpaw.hubbard import hubbard 

7 

8 

9def test_hubbard_U_noncollinear(): 

10 rng = np.random.default_rng(seed=555) 

11 

12 # Generate random data simulating a setup with 

13 # s and p bounded and unbounded partial waves. 

14 

15 D_sii = (rng.random([4, 8, 8], dtype=np.float64) 

16 + 1j * rng.random([4, 8, 8], dtype=np.float64)) 

17 # Make Hermitian 

18 D_sii = (D_sii + np.transpose(D_sii.conj(), (0, 2, 1))) / 2 

19 n_j = [2, 2, -1, -1] 

20 l_j = [0, 1, 0, 1] 

21 

22 # Hubbard correction on the p orbitals 

23 l = 1 

24 U = 4 # eV 

25 

26 # Generate inner products for the p-orbitals 

27 N0_q = np.zeros(len(l_j) * (len(l_j) + 1) // 2) 

28 N0_q[4] = 0.9 # Bounded-bounded 

29 N0_q[6] = 1.1 # Bounded-unbounded 

30 N0_q[9] = 1.2 # Unbounded-unbounded 

31 

32 eU, V_sii = hubbard(D_sii, U / Ha, l, l_j, n_j, N0_q, scale=True) 

33 eU_onlyreal, _ = hubbard(D_sii.real, U / Ha, l, l_j, n_j, N0_q, scale=True) 

34 assert eU.imag == pytest.approx(0, abs=1e-10) 

35 

36 # Assert that disregarding the complex values will alter the energy 

37 assert eU == pytest.approx(-7.825104, abs=1.e-2) 

38 assert eU_onlyreal == pytest.approx(-7.546553, abs=1.e-2) 

39 

40 # Now, we check that the Hubbard Hamiltonian is calculated correctly by 

41 # comparing it with an energy derivative w.r.t. a density matrix element 

42 # calculated through finite difference. 

43 

44 diff = 1e-5 

45 

46 D1_sii = D_sii.copy() 

47 D1_sii[0, 1, 2] += diff / 2 

48 

49 eU1, _ = hubbard(D1_sii, U / Ha, l, l_j, n_j, N0_q, scale=True) 

50 

51 D2_sii = D_sii.copy() 

52 D2_sii[0, 1, 2] -= diff / 2 

53 

54 eU2, _ = hubbard(D2_sii, U / Ha, l, l_j, n_j, N0_q, scale=True) 

55 

56 # print(eU1) 

57 # print(eU2) 

58 # print(V_sii[0, 1:4, 1:4]) 

59 assert (eU1 - eU2) / diff == pytest.approx(V_sii[0, 1, 2], abs=1e-8)