Coverage for /builds/debichem-team/python-ase/ase/build/surfaces_with_termination.py: 98.73%

1import numpy as np 

2 

3from ase.build.general_surface import surface 

4from ase.geometry import get_layers 

5from ase.symbols import string2symbols 

6 

7 

8def surfaces_with_termination(lattice, indices, layers, vacuum=None, tol=1e-10, 

9 termination=None, return_all=False, 

10 verbose=False): 

11 """Create surface from a given lattice and Miller indices with a given 

12 termination 

13 

14 Parameters 

15 ========== 

16 lattice: Atoms object or str 

17 Bulk lattice structure of alloy or pure metal. Note that the 

18 unit-cell must be the conventional cell - not the primitive cell. 

19 One can also give the chemical symbol as a string, in which case the 

20 correct bulk lattice will be generated automatically. 

21 indices: sequence of three int 

22 Surface normal in Miller indices (h,k,l). 

23 layers: int 

24 Number of equivalent layers of the slab. (not the same as the layers 

25 you choose from for terminations) 

26 vacuum: float 

27 Amount of vacuum added on both sides of the slab. 

28 termination: str 

29 the atoms you wish to be in the top layer. There may be many such 

30 terminations, this function returns all terminations with the same 

31 atomic composition. 

32 e.g. 'O' will return oxygen terminated surfaces. 

33 e.g.'TiO' returns surfaces terminated with layers containing both 

34 O and Ti 

35 Returns: 

36 return_surfs: List 

37 a list of surfaces that match the specifications given 

38 

39 """ 

40 lats = translate_lattice(lattice, indices) 

41 return_surfs = [] 

42 check = [] 

43 check2 = [] 

44 for item in lats: 

45 too_similar = False 

46 surf = surface(item, indices, layers, vacuum=vacuum, tol=tol) 

47 surf.wrap(pbc=[True] * 3) # standardize slabs 

48 

49 positions = surf.get_scaled_positions().flatten() 

50 for i, value in enumerate(positions): 

51 if value >= 1 - tol: # move things closer to zero within tol 

52 positions[i] -= 1 

53 surf.set_scaled_positions(np.reshape(positions, (len(surf), 3))) 

54 # rep = find_z_layers(surf) 

55 z_layers, _hs = get_layers(surf, (0, 0, 1)) # just z layers matter 

56 # get the indicies of the atoms in the highest layer 

57 top_layer = [ 

58 i for i, val in enumerate( 

59 z_layers == max(z_layers)) if val] 

60 

61 if termination is not None: 

62 comp = [surf.get_chemical_symbols()[a] for a in top_layer] 

63 term = string2symbols(termination) 

64 # list atoms in top layer and not in requested termination 

65 check = [a for a in comp if a not in term] 

66 # list of atoms in requested termination and not in top layer 

67 check2 = [a for a in term if a not in comp] 

68 if len(return_surfs) > 0: 

69 pos_diff = [a.get_positions() - surf.get_positions() 

70 for a in return_surfs] 

71 for i, su in enumerate(pos_diff): 

72 similarity_test = su.flatten() < tol * 1000 

73 if similarity_test.all(): 

74 # checks if surface is too similar to another surface 

75 too_similar = True 

76 if too_similar: 

77 continue 

78 if return_all is True: 

79 pass 

80 elif check != [] or check2 != []: 

81 continue 

82 return_surfs.append(surf) 

83 return return_surfs 

84 

85 

86def translate_lattice(lattice, indices, tol=10**-3): 

87 """translates a bulk unit cell along a normal vector given by the a set of 

88 miller indices to the next symetric position. This is used to control the 

89 termination of the surface in the smart_surface command 

90 Parameters: 

91 ========== 

92 lattice: Atoms object 

93 atoms object of the bulk unit cell 

94 indices: 1x3 list,tuple, or numpy array 

95 the miller indices you wish to cut along. 

96 returns: 

97 lattice_list: list of Atoms objects 

98 a list of all the different translations of the unit cell that will 

99 yield different terminations of a surface cut along the miller 

100 indices provided. 

101 """ 

102 lattice_list = [] 

103 cell = lattice.get_cell() 

104 pt = [0, 0, 0] 

105 h, k, l = indices # noqa (E741 ambiguous name 'l') 

106 millers = list(indices) 

107 for index, item in enumerate(millers): 

108 if item == 0: 

109 millers[index] = 10**9 # make zeros large numbers 

110 elif pt == [0, 0, 0]: # for numerical stability 

111 pt = list(cell[index] / float(item) / np.linalg.norm(cell[index])) 

112 h1, k1, l1 = millers 

113 N = np.array(cell[0] / h1 + cell[1] / k1 + cell[2] / l1) 

114 n = N / np.linalg.norm(N) # making a unit vector normal to cut plane 

115 # finding distance from cut plan vector 

116 d = [np.round(np.dot(n, (a - pt)) * n, 5) for 

117 a in lattice.get_scaled_positions()] 

118 duplicates = [] 

119 for i, item in enumerate(d): 

120 g = [True for a in d[i + 1:] if np.linalg.norm(a - item) < tol] 

121 if g != []: 

122 duplicates.append(i) 

123 duplicates.reverse() 

124 for i in duplicates: 

125 del d[i] 

126 # put distance to the plane at the end of the array 

127 for i, item in enumerate(d): 

128 d[i] = np.append(item, 

129 np.dot(n, (lattice.get_scaled_positions()[i] - pt))) 

130 d = np.array(d) 

131 d = d[d[:, 3].argsort()] # sort by distance to the plane 

132 d = [a[:3] for a in d] # remove distance 

133 d = list(d) # make it a list again 

134 for i in d: 

135 """ 

136 The above method gives you the boundries of between terminations that 

137 will allow you to build a complete set of terminations. However, it 

138 does not return all the boundries. Thus you must check both above and 

139 below the boundary, and not stray too far from the boundary. If you move 

140 too far away, you risk hitting another boundary you did not find. 

141 """ 

142 lattice1 = lattice.copy() 

143 displacement = (h * cell[0] + k * cell[1] + l * cell[2]) \ 

144 * (i + 10 ** -8) 

145 lattice1.positions -= displacement 

146 lattice_list.append(lattice1) 

147 lattice1 = lattice.copy() 

148 displacement = (h * cell[0] + k * cell[1] + l * cell[2]) \ 

149 * (i - 10 ** -8) 

150 lattice1.positions -= displacement 

151 lattice_list.append(lattice1) 

152 return lattice_list