Coverage for tbkit/graphene.py: 100%

102 statements  

« prev     ^ index     » next       coverage.py v7.16.0, created at 2026-09-22 13:16 +0100

1from __future__ import annotations 

2 

3from tbkit.lattice import * 

4from tbkit.plot import * 

5from tbkit.system import * 

6from math import sqrt 

7 

8PI = np.pi 

9ATOL = 1e-3 

10DX = 0.5 * sqrt(3) 

11DY = 0.5 

12 

13 

14################################# 

15# CLASS GRAPHENE 

16################################# 

17 

18 

19class GrapheneLattice(Lattice): 

20 def __init__(self) -> None: 

21 unit_cell = [{'tag': 'a', 'r0': (0, 0)}, 

22 {'tag': 'b', 'r0': (DX, DY)}] 

23 prim_vec = [(2*DX, 0.), (DX, 1.5)] 

24 Lattice.__init__(self, unit_cell=unit_cell, prim_vec=prim_vec) 

25 self.butterfly = np.array([]) 

26 self.betas = np.array([]) 

27 

28 def triangle_zigzag(self, n: int) -> None: 

29 ''' 

30 Triangular flake with zigzag terminations. 

31 

32 :param: n. Int. Number of plackets along the edges.  

33 ''' 

34 error_handling.positive_int(n, 'n') 

35 self.get_lattice(n1=n+2, n2=n+2) 

36 self.boundary_line(cx=-sqrt(3), cy=-1, co=-3*n+2) 

37 

38 def hexagon_zigzag(self, n: int) -> None: 

39 ''' 

40 Hexagonal flake with zigzag terminations. 

41 

42 :param: n. Int. Number of plackets along the edges.  

43 ''' 

44 error_handling.positive_int(n, 'n') 

45 self.get_lattice(n1=2*n, n2=2*n) 

46 self.boundary_line(cx=sqrt(3), cy=1, co=3*(n-1)) 

47 self.boundary_line(cx=-sqrt(3), cy=-1, co=-9*n+2.5) 

48 

49 def triangle_armchair(self, n: int) -> None: 

50 ''' 

51 Triangular flake with armchair terminations. 

52 

53 :param: n. Int. Number of plackets along the edges.  

54 ''' 

55 error_handling.positive_int(n, 'n') 

56 self.get_lattice(n1=2*n, n2=2*n) 

57 self.boundary_line(cx=1, cy=0, co=sqrt(3)/2*(2*n-1)-0.1) 

58 self.boundary_line(cx=-1/sqrt(3), cy=1, co=-n-0.1) 

59 self.boundary_line(cx=-1/sqrt(3), cy=-1, co=-4*n+0.1) 

60 

61 def hexagon_armchair(self, n: int) -> None: 

62 ''' 

63 Hexagonal flake with armchair terminations. 

64 

65 :param: n. Int. Number of plackets along each edge.  

66 ''' 

67 error_handling.positive_int(n, 'n') 

68 nn = 3 * n - 2 

69 self.get_lattice(n1=2*nn, n2=2*nn) 

70 self.boundary_line(cx=1, cy=0, co=sqrt(3)/2* (2*nn-1)-.1) 

71 self.boundary_line(cx=-1/sqrt(3), cy=1, co=-nn -.1) 

72 self.boundary_line(cx=-1/sqrt(3), cy=-1, co=-4*nn +1-.1) 

73 self.coor['x'] -= self.coor['x'].min() 

74 self.boundary_line(cx=-1, cy=0, co=-DX * 2*nn-0.1) 

75 self.boundary_line(cx=1/sqrt(3), cy=1, co=3*n -3.) 

76 self.boundary_line(cx=1/sqrt(3), cy=-1, co=-6*n+4) 

77 

78 def square(self, n: int) -> None: 

79 ''' 

80 Squared flake. 

81 

82 :param: n. Int. Number of plackets along x.  

83 ''' 

84 error_handling.positive_int(n, 'n') 

85 n2 = int(1.5*DX*n) 

86 self.get_lattice(n1=2*n, n2=n2) 

87 self.boundary_line(cx=1, cy=0, co=DX*(2*n-2)+0.1) 

88 self.boundary_line(cx=-1, cy=0, co=-DX*(4*n)+0.5) 

89 

90 def circle(self, n: int) -> None: 

91 ''' 

92 Circular flake. 

93 

94 :param: n. Int. Number of plackets along the diameter.  

95 ''' 

96 error_handling.positive_int(n, 'n') 

97 self.get_lattice(n1=2*n, n2=2*n) 

98 self.sites = len(self.coor) 

99 self.center() 

100 if n % 2 == 0: 

101 self.shift_x(shift=-DX) 

102 self.ellipse_in(rx=DX*(n+1), ry=DX*(n+1), x0=0., y0=0.) 

103 self.remove_dangling() 

104 

105class GrapheneSystem(System): 

106 def __init__(self, lat: Lattice) -> None: 

107 System.__init__(self, lat) 

108 

109 def set_hop_linear_strain(self, t: complex, beta: float) -> None: 

110 ''' 

111 Set nearest neighbors hoppings according to the linear trixial strain.  

112  

113 :param t: Hopping value without strain. 

114 :param beta: Strength of the strain. 

115 ''' 

116 error_handling.number(t, 't') 

117 error_handling.real_number(beta, 'beta') 

118 self.get_distances() 

119 ind = np.argwhere(np.isclose(self.dist_uni[1], self.vec_hop['dis'], atol=ATOL)) 

120 ind_up = ind[ind[:, 1] > ind[:, 0]] 

121 self.hop = np.zeros(len(ind_up), dtype=HOP_DTYPE) 

122 self.hop['n'] = 1 

123 self.hop['i'] = ind_up[:, 0] 

124 self.hop['j'] = ind_up[:, 1] 

125 self.hop['ang'] = self.vec_hop['ang'][ind_up[:, 0], ind_up[:, 1]] 

126 # change angle (to get the correct strain) 

127 self.hop['ang'][np.isclose(30., self.hop['ang'], ATOL)] = -150. 

128 self.hop['ang'][np.isclose(150., self.hop['ang'], ATOL)] = - 30. 

129 x_center = .5 * (self.lat.coor['x'][ind_up[:, 0]] + self.lat.coor['x'][ind_up[:, 1]]) 

130 y_center = .5 * (self.lat.coor['y'][ind_up[:, 0]] + self.lat.coor['y'][ind_up[:, 1]]) 

131 self.hop['t'] = t * (1. + 0.25 * beta * (np.cos(PI / 180 * self.hop['ang']) * x_center + 

132 np.sin(PI / 180 * self.hop['ang']) * y_center)) 

133 # back to the former angle 

134 self.hop['ang'][np.isclose(-150., self.hop['ang'])] = 30. 

135 self.hop['ang'][np.isclose(-30., self.hop['ang'])] = 150. 

136 

137 def get_butterfly(self, t: complex, N: int) -> None: 

138 '''' 

139 Get energies depending on strain. 

140 

141 :param t: Unstrained hopping value. 

142 :param N: number of strain values between min and max strains. 

143 ''' 

144 error_handling.number(t, 't') 

145 error_handling.positive_int(N, 'N') 

146 beta_lims = self.get_beta_lims() 

147 self.betas = np.linspace(beta_lims[0], beta_lims[1], N) 

148 self.butterfly = np.zeros((N, self.lat.sites)) 

149 for i, beta in enumerate(self.betas): 

150 self.set_hop_linear_strain(t=1, beta=beta) 

151 self.get_ham() 

152 self.butterfly[i] = LA.eigvalsh(self.ham.toarray()) 

153 

154 def get_beta_lims(self) -> NDArray[np.float64]: 

155 ''' 

156 Get the extremal values of strain keeping positive hoppings. 

157 ''' 

158 beta_lims = np.zeros(2) 

159 yb_min_val = self.lat.coor['y'][self.lat.coor['tag'] == 'b'].min() 

160 yb_min = self.lat.coor['y'][self.lat.coor['y'] == yb_min_val][0] 

161 ym = 0.5 * (2 * yb_min + 1) 

162 beta_lims[1] = -4. / ym + 1e-6 

163 yb_max_val = self.lat.coor['y'][self.lat.coor['tag'] == 'a'].max() 

164 yb_max = self.lat.coor['y'][self.lat.coor['y'] == yb_max_val][0] 

165 ym = 0.5 * (2 * yb_max - 1) 

166 beta_lims[0] = -4. / ym + 1e-6 

167 print('Strain limits: {}'.format(beta_lims)) 

168 return beta_lims 

169 

170 

171# Backward-compatible camelCase aliases (pre-0.2 API). 

172grapheneLat = GrapheneLattice 

173grapheneSys = GrapheneSystem