Coverage for gpaw/test/test_lvlmm.py: 100%

1import numpy as np 

2import pytest 

3from gpaw.spinorbit import get_L_vlmm 

4 

5 

6class QuantumState: 

7 """Quantum State, representing real Y_lm or complex Y_l^m""" 

8 def __init__(self, l: int, m: int, factor: complex = 1.0): 

9 assert type(l) is int 

10 assert type(m) is int 

11 self.l = l 

12 self.m = m 

13 self.factor = factor 

14 

15 

16def real_to_complex_quantum_state(state: QuantumState): 

17 """represent the real Y_lm in complex Y_l^m""" 

18 if state.m == 0: 

19 return [QuantumState(state.l, state.m, 1.0)] 

20 elif state.m > 0: 

21 return [ 

22 QuantumState(state.l, -1 * state.m, 

23 1 / np.sqrt(2.0) * state.factor), 

24 QuantumState(state.l, state.m, 

25 (-1.0) ** state.m / np.sqrt(2.0) * state.factor)] 

26 else: # state.m < 0: 

27 return [ 

28 QuantumState(state.l, state.m, 1j / np.sqrt(2.0) * state.factor), 

29 QuantumState( 

30 state.l, 

31 -1 * state.m, 

32 -1j * (-1.0) ** state.m / np.sqrt(2.0) * state.factor, 

33 ), 

34 ] 

35 

36 

37def matrix_element(bra: QuantumState, ket: QuantumState, operator): 

38 """<Y_lm' | O | Y_lm >""" 

39 assert bra.l == ket.l 

40 result = 0.0 

41 kets_complex = real_to_complex_quantum_state(ket) 

42 bras_complex = real_to_complex_quantum_state(bra) 

43 for b in bras_complex: 

44 for k in kets_complex: 

45 operator_factor, newk = operator(k) 

46 if newk.m != b.m: 

47 continue 

48 result += np.conj(b.factor) * newk.factor * operator_factor 

49 return result 

50 

51 

52def Lz(ket: QuantumState): 

53 """L_z acting on complex Y_l^m""" 

54 return ket.m, QuantumState(ket.l, ket.m, ket.factor) 

55 

56 

57def Lp(ket: QuantumState): 

58 """L_+ acting on complex Y_l^m""" 

59 if ket.m < ket.l: 

60 result = np.sqrt(ket.l * (ket.l + 1) - ket.m * (ket.m + 1)) 

61 else: 

62 result = 0 

63 return result, QuantumState(ket.l, ket.m + 1, ket.factor) 

64 

65 

66def Lm(ket: QuantumState): 

67 """L_- acting on complex Y_l^m""" 

68 if ket.m > -ket.l: 

69 result = np.sqrt(ket.l * (ket.l + 1) - ket.m * (ket.m - 1)) 

70 else: 

71 result = 0 

72 return result, QuantumState(ket.l, ket.m - 1, ket.factor) 

73 

74 

75def test_lvlmm(): 

76 

77 L_vlmm = get_L_vlmm() 

78 

79 for l in [0, 1, 2, 3]: 

80 # print(f"Doing {'spdf'[l]}:") 

81 result_lz = np.zeros((2 * l + 1, 2 * l + 1), dtype=complex) 

82 result_lp = np.zeros((2 * l + 1, 2 * l + 1), dtype=complex) 

83 result_lm = np.zeros((2 * l + 1, 2 * l + 1), dtype=complex) 

84 for mp in range(-l, l + 1): 

85 bra = QuantumState(l, mp) # real Y_lm 

86 for m in range(-l, l + 1): 

87 ket = QuantumState(l, m) # real Y_lm 

88 result_lz[mp + l, m + l] = matrix_element(bra, ket, Lz) 

89 result_lp[mp + l, m + l] = matrix_element(bra, ket, Lp) 

90 result_lm[mp + l, m + l] = matrix_element(bra, ket, Lm) 

91 

92 result_lx = (result_lp + result_lm) / 2.0 

93 result_ly = (result_lp - result_lm) / 2j 

94 

95 L_vlmm[0][l] == pytest.approx(result_lx) 

96 L_vlmm[1][l] == pytest.approx(result_ly) 

97 L_vlmm[2][l] == pytest.approx(result_lz)