Coverage for gpaw/test/response/test_pdens_tool.py: 100%

47 statements  

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

1import pytest 

2import numpy as np 

3 

4from gpaw.response.tool import (get_bz_transitions, 

5 get_chi0_integrand, 

6 get_degeneracy_matrix, 

7 get_individual_transition_strengths) 

8 

9 

10@pytest.mark.response 

11def test_response_pdens_tool(in_tmp_dir, gpw_files): 

12 """Calculate optical transition strengths.""" 

13 spins = 'all' 

14 q_c = [0., 0., 0.] 

15 bzk_kc = np.array([[0., 0., 0.]]) 

16 

17 pair, qpd, domain = get_bz_transitions( 

18 gpw_files['silicon_pdens_tool'], q_c, bzk_kc, spins=spins, ecut=10) 

19 

20 nocc1, nocc2 = pair.gs.count_occupied_bands(1e-6) 

21 # XXX should we know 1e-6? 

22 

23 # non-empty bands 

24 n_n = np.arange(0, nocc2) 

25 # not completely filled bands 

26 m_m = np.arange(nocc1, pair.gs.bd.nbands) 

27 

28 nt = len(domain) 

29 nn = len(n_n) 

30 nm = len(m_m) 

31 nG = qpd.ngmax 

32 optical_limit = np.allclose(q_c, 0.) 

33 

34 n_tnmG = np.zeros((nt, nn, nm, nG + 2 * optical_limit), dtype=complex) 

35 df_tnm = np.zeros((nt, nn, nm), dtype=float) 

36 eps_tn = np.zeros((nt, nn), dtype=float) 

37 eps_tm = np.zeros((nt, nm), dtype=float) 

38 

39 for t, point in enumerate(domain): 

40 (n_tnmG[t], df_tnm[t], 

41 eps_tn[t], eps_tm[t]) = get_chi0_integrand(pair, qpd, 

42 n_n, m_m, point) 

43 

44 testNM_ibN = [[[0], [4, 5, 6]], [[0], [7]], 

45 [[1, 2, 3], [4, 5, 6]], [[1, 2, 3], [7]]] 

46 testts_iG = np.array([[0.07, 0.07, 0.07], [0.00, 0.00, 0.00], 

47 [51.34, 51.34, 51.34], [22.69, 22.69, 22.69]]) 

48 

49 for t, (point, n_nmG, 

50 df_nm, eps_n, eps_m) in enumerate(zip(domain, n_tnmG, df_tnm, 

51 eps_tn, eps_tm)): 

52 

53 # Find degeneracies 

54 degmat_Nn, eps_N = get_degeneracy_matrix(eps_n, tol=1.e-2) 

55 degmat_Mm, eps_M = get_degeneracy_matrix(eps_m, tol=1.e-2) 

56 

57 # Find diagonal transition strengths 

58 its_nmG = np.zeros((nn, nm, 1 + 2 * optical_limit)) 

59 for G in range(1 + 2 * optical_limit): 

60 its_nmG[:, :, G] = get_individual_transition_strengths( 

61 n_nmG, df_nm, 

62 G, G) 

63 

64 # Find unique transition strengths 

65 its_NmG = np.tensordot(degmat_Nn, its_nmG, axes=(1, 0)) 

66 ts_MNG = np.tensordot(degmat_Mm, its_NmG, axes=(1, 1)) 

67 ts_NMG = np.transpose(ts_MNG, (1, 0, 2)) 

68 

69 i = 0 

70 for N, ts_MG in enumerate(ts_NMG): 

71 for M, ts_G in enumerate(ts_MG): 

72 degN_n = n_n[np.where(degmat_Nn[N])] 

73 degM_m = m_m[np.where(degmat_Mm[M])] 

74 

75 for testn, n in zip(testNM_ibN[i][0], degN_n): 

76 assert testn == pytest.approx(n, abs=0.5) 

77 for testm, m in zip(testNM_ibN[i][1], degM_m): 

78 assert testm == pytest.approx(m, abs=0.5) 

79 

80 for testts, ts in zip(testts_iG[i], ts_G): 

81 print(ts) 

82 assert testts == pytest.approx(ts, abs=0.01) 

83 

84 i += 1