Coverage for gpaw/response/heisenberg.py: 98%

1"""Methods to calculate material properties in the Heisenberg model. 

2Primarily based on the magnon dispersion relations.""" 

3 

4import numpy as np 

5 

6 

7def calculate_single_site_magnon_energies(J_qx, q_qc, mm): 

8 """Compute the magnon energies from the isotropic exchange constants of a 

9 system with a single magnetic site in the unit cell, as a function of the 

10 wave vector q: 

11 

12 ħω(q) = g μ_B / M [J(0) - J(q)] 

13 

14 Parameters 

15 ---------- 

16 J_qx : np.ndarray 

17 Isotropic exchange constants as a function of q. 

18 J_qx can have any number of additional dimensions x, which will be 

19 treated independently. 

20 q_qc : np.ndarray 

21 q-vectors in relative coordinates. Has to include q=0. 

22 mm : float 

23 Magnetic moment of the site in μ_B. 

24 

25 Returns 

26 ------- 

27 E_qx : np.ndarray 

28 Magnon energies as a function of q and x. Same shape as input J_qx. 

29 """ 

30 assert J_qx.shape[0] == q_qc.shape[0] 

31 assert np.allclose(J_qx.imag, 0.), \ 

32 "The exchange constants of an isotropic system with a single "\ 

33 "magnetic sublattice are always real." 

34 

35 q0 = get_q0_index(q_qc) 

36 J0_x = J_qx[q0] 

37 

38 # Compute energies 

39 E_qx = 2. / mm * (J0_x[np.newaxis, ...] - J_qx) 

40 

41 return E_qx.real 

42 

43 

44def calculate_fm_magnon_energies(J_qabx, q_qc, mm_ax): 

45 """Compute the magnon eigenmode energies from the isotropic exchange 

46 constants of a ferromagnetic system with an arbitrary number of magnetic 

47 sites in the unit cell, as a function of the wave vector q. 

48 

49 The eigenmodes are calculated as the eigenvalues to the dynamic spin wave 

50 matrix: 

51 

52 H^ab(q) = g μ_B / sqrt(M_a M_b) [Σ_c J^ac(0) δ_ab - J^ab(q)] 

53 

54 Parameters 

55 ---------- 

56 J_qabx : np.ndarray 

57 Isotropic exchange constants as a function of q and sublattice indices 

58 a and b. J_qabx can have any number of additional dimensions x, which 

59 will be treated independently. 

60 q_qc : np.ndarray 

61 q-vectors in relative coordinates. Has to include q=0. 

62 mm_ax : np.ndarray 

63 Magnetic moments of the sublattice sites in μ_B. 

64 

65 Returns 

66 ------- 

67 E_qnx : np.ndarray 

68 Magnon eigenmode energies as a function of q, mode index n and x. 

69 """ 

70 H_qabx = generate_fm_dynamic_spin_wave_matrix(J_qabx, q_qc, mm_ax) 

71 

72 # Move magnetic site axes in order to prepare for np.linalg.eig 

73 H_qbxa = np.moveaxis(H_qabx, 1, -1) 

74 H_qxab = np.moveaxis(H_qbxa, 1, -1) 

75 

76 # Diagonalize the matrix 

77 E_qxn, _ = np.linalg.eig(H_qxab) 

78 

79 # Transpose to output format 

80 E_qnx = np.moveaxis(E_qxn, -1, 1) 

81 

82 # Eigenvalues should be real, otherwise input J_qabx has been invalid 

83 assert np.allclose(E_qnx.imag, 0.) 

84 E_qnx = E_qnx.real 

85 

86 return E_qnx 

87 

88 

89def generate_fm_dynamic_spin_wave_matrix(J_qabx, q_qc, mm_ax): 

90 """Generate the dynamic spin wave matrix from the isotropic exchange 

91 constants of a ferromagnet: 

92 

93 H^ab(q) = g μ_B / sqrt(M_a M_b) [Σ_c J^ac(0) δ_ab - J^ab(q)] 

94 

95 Parameters 

96 ---------- 

97 J_qabx : np.ndarray 

98 Isotropic exchange constants as a function of q and sublattice indices 

99 a and b. J_qabx can have any number of additional dimensions x, which 

100 will be treated independently. 

101 q_qc : np.ndarray 

102 q-vectors in relative coordinates. Has to include q=0. 

103 mm_ax : np.ndarray 

104 Magnetic moments of the sublattice sites in μ_B. 

105 

106 Returns 

107 ------- 

108 H_qabx : np.ndarray 

109 Dynamic spin wave matrix. Has the same shape as the input J_qabx 

110 """ 

111 assert len(J_qabx.shape) >= 3 

112 assert J_qabx.shape[1] == J_qabx.shape[2] 

113 assert J_qabx.shape[1] == mm_ax.shape[0] 

114 assert J_qabx.shape[0] == q_qc.shape[0] 

115 assert J_qabx.shape[3:] == mm_ax.shape[1:] 

116 assert np.allclose(J_qabx, np.conj(np.moveaxis(J_qabx, 2, 1))), \ 

117 "The isotropic exchange constants J^ab(q) are Hermitian by definition." 

118 na = mm_ax.shape[0] 

119 

120 # Get J^ab(0) 

121 q0 = get_q0_index(q_qc) 

122 J0_acx = J_qabx[q0] 

123 assert np.allclose(J0_acx.imag, 0.), \ 

124 "For a collinear system without spin-orbit coupling, the exchange "\ 

125 "constants are not only isotropic, but also reciprocal. In "\ 

126 "particular, this implies that [J^ab(q)]^*=J^ab(-q)." 

127 

128 # Set up magnetic moment prefactor as outer product 

129 mm_inv_abx = 2. / np.sqrt(mm_ax[:, np.newaxis, ...] 

130 * mm_ax[np.newaxis, ...]) 

131 

132 # Calculate diagonal component Σ_c J^ac(0) δ_ab 

133 J0_ax = np.sum(J0_acx, axis=1) 

134 diagonal_mapping = np.zeros((na, na, na)) 

135 np.fill_diagonal(diagonal_mapping, 1) 

136 J0_abx = np.tensordot(diagonal_mapping, J0_ax, axes=(2, 0)) 

137 

138 # Calculate the dynamic spin wave matrix 

139 H_qabx = mm_inv_abx[np.newaxis, ...] * (J0_abx[np.newaxis, ...] - 

140 J_qabx) 

141 

142 return H_qabx 

143 

144 

145def get_q0_index(q_qc): 

146 """Find index corresponding q=0 in q-vector array.""" 

147 q0_indices = np.argwhere(np.all(q_qc == 0, axis=1)) 

148 

149 if len(q0_indices) >= 1: 

150 return int(q0_indices[0, 0]) 

151 else: 

152 raise ValueError('q_qc has to include q=0, i.e. q_c = [0., 0., 0.]')