Coverage for gpaw/test/xc/test_xc.py: 98%

1import pytest 

2import numpy as np 

3from gpaw.xc.libxc import LibXC, short_names 

4from gpaw.xc.kernel import XCKernel, codes 

5from gpaw.xc import libraries 

6 

7HAS_LIBXC_7 = int((libraries['libxc'] or '0.0').split('.')[0]) >= 7 

8 

9 

10@pytest.mark.ci 

11@pytest.mark.libxc 

12def test_xc_xc(): 

13 funcs = [] 

14 modes = [] 

15 for name, value in short_names.items(): 

16 try: 

17 if 'MGGA' in value and not HAS_LIBXC_7: 

18 with pytest.warns(UserWarning, 

19 match='should be compiled with'): 

20 LibXC(name) 

21 else: 

22 LibXC(name) 

23 except NameError: 

24 continue 

25 if name == 'SCAN': 

26 continue 

27 funcs.append(name) 

28 modes.append(0) 

29 for name in codes: 

30 funcs.append(name) 

31 modes.append(1) 

32 

33 def create_xc(func, mode): 

34 if mode == 0: 

35 xc = LibXC(func, disable_fhc=False) 

36 else: 

37 xc = XCKernel(func) 

38 return xc 

39 

40 def f1(n_xg, xc): 

41 e_g = np.empty_like(n_xg[0]) 

42 n_sg = n_xg[:1] 

43 sigma_xg = n_xg[1:2] 

44 tau_sg = n_xg[2:] 

45 dedn_sg = np.zeros_like(n_sg) 

46 dedsigma_xg = np.zeros_like(sigma_xg) 

47 dedtau_sg = np.zeros_like(tau_sg) 

48 xc.calculate(e_g, 

49 n_sg, dedn_sg, 

50 sigma_xg, dedsigma_xg, 

51 tau_sg, dedtau_sg) 

52 return e_g, np.concatenate((dedn_sg, dedsigma_xg, dedtau_sg)) 

53 

54 def f2(n_xg, xc): 

55 e_g = np.empty_like(n_xg[0]) 

56 n_sg = n_xg[:2] 

57 sigma_xg = n_xg[2:5] 

58 # Triangle inequality: 

59 assert (sigma_xg[1] <= (sigma_xg[0] + sigma_xg[2]) / 2).all() 

60 assert (sigma_xg[1] >= -(sigma_xg[0] + sigma_xg[2]) / 2).all() 

61 tau_sg = n_xg[5:] 

62 dedn_sg = np.zeros_like(n_sg) 

63 dedsigma_xg = np.zeros_like(sigma_xg) 

64 dedtau_sg = np.zeros_like(tau_sg) 

65 xc.calculate(e_g, 

66 n_sg, dedn_sg, 

67 sigma_xg, dedsigma_xg, 

68 tau_sg, dedtau_sg) 

69 return e_g, np.concatenate((dedn_sg, dedsigma_xg, dedtau_sg)) 

70 

71 eps = 1.0e-6 

72 

73 n_xg = np.array( 

74 [[0.2, 0.01, 0.4], 

75 [0.2, 0.1, 0.5], 

76 [0.01, 0.01, 0.2], 

77 [0.1, 0.3, 0.5]]).T.copy() 

78 

79 for i, func in enumerate(funcs): 

80 print(i, func) 

81 xc = create_xc(funcs[i], modes[i]) 

82 e0_g, d0_xg = f1(n_xg, xc) 

83 d_xg = np.empty_like(d0_xg) 

84 for x, n_g in enumerate(n_xg): 

85 m_xg = n_xg.copy() 

86 m_xg[x] += eps 

87 d_xg[x] = 0.5 * f1(m_xg, xc)[0] / eps 

88 m_xg[x] -= 2 * eps 

89 d_xg[x] -= 0.5 * f1(m_xg, xc)[0] / eps 

90 

91 ns_xg = np.empty((7, len(n_g))) 

92 ns_xg[:2] = n_xg[0] / 2 

93 ns_xg[2:5] = n_xg[1] / 4 

94 ns_xg[5:] = n_xg[2] / 2 

95 es_g, ds_xg = f2(ns_xg, xc) 

96 error = (abs(d0_xg - d_xg).max() + 

97 abs(es_g - e0_g).max() + 

98 abs(ds_xg[:2] - d0_xg[0]).max() + 

99 abs(ds_xg[2:5].sum(0) / 4 - d0_xg[1]).max() + 

100 abs(ds_xg[5:] - d0_xg[2]).max()) 

101 assert error == pytest.approx(0, abs=6e-9) 

102 del xc 

103 

104 # Numbers from old lxc_xc.py test: 

105 na = 2.0 

106 nb = 1.0 

107 sigma0 = 2.0 # (0.0, 1.0, 1.0) 

108 sigma1 = 2.0 

109 sigma2 = 5.0 # (1.0, 2.0, 0.0) 

110 taua = (3 * np.pi**2)**(2. / 3.) * na**(5. / 3.) / 2 * sigma0 

111 taub = (3 * np.pi**2)**(2. / 3.) * nb**(5. / 3.) / 2 * sigma2 

112 

113 n_xg = np.array( 

114 [[na, nb, sigma0, sigma1, sigma2, taua, taub], 

115 [0.1, 0.1, 0.025, 0.02, 0.025, 0.25, 0.25], 

116 [0.1, 0.1, 0.125, 0.12, 0.125, 0.0025, 0.025], 

117 [0.1, 0.1, 0.01, 0.01, 0.015, 0.2, 0.2], 

118 [0.1, 0.2, 0.1, -0.08, 0.10, 0.01, 0.05], 

119 [0.1, 0.1, 0.1, 0.01, 0.01, 0.01, 0.01], 

120 [0.1, 0.1, 0.1, 0.14, 0.20, 0.01, 0.05]]).T.copy() 

121 

122 for i, func in enumerate(funcs): 

123 xc = create_xc(funcs[i], modes[i]) 

124 if xc.type == 'MGGA': 

125 N_xg = n_xg[:, :1].copy() 

126 else: 

127 N_xg = n_xg 

128 e0_g, d0_xg = f2(N_xg, xc) 

129 d_xg = np.empty_like(d0_xg) 

130 for x, n_g in enumerate(N_xg): 

131 m_xg = N_xg.copy() 

132 m_xg[x] += eps 

133 d_xg[x] = 0.5 * f2(m_xg, xc)[0] / eps 

134 m_xg[x] -= 2 * eps 

135 d_xg[x] -= 0.5 * f2(m_xg, xc)[0] / eps 

136 assert abs(d0_xg - d_xg).max() == pytest.approx(0, abs=2e-8) 

137 del xc