Coverage for hiphive/input_output/gpumd.py: 98%

75 statements  

« prev     ^ index     » next       coverage.py v7.10.1, created at 2025-08-01 17:04 +0000

1from itertools import permutations, product 

2import ase 

3import numpy as np 

4 

5 

6def write_fcs_gpumd( 

7 fname_fc: str, 

8 fname_clusters: str, 

9 fcs, 

10 order: int, 

11 tol: float = 1e-10, 

12) -> None: 

13 """ 

14 Writes force constants of given order in GPUMD format. 

15 

16 Parameters 

17 ---------- 

18 fname_fc 

19 Name of file which contains the lookup force constants. 

20 fname_clusters 

21 Name of file which contains the clusters and the force constant lookup index. 

22 fcs 

23 Force constants in :class:`ForceConstants` format. 

24 order 

25 Force constants for this order will be written to file. 

26 tol 

27 If the norm of a force constant term is less than this value it will be excluded 

28 from the output; 

29 if two force-constants differ by this value or less, they are considered equal. 

30 """ 

31 cluster_lookup, fc_lookup = _get_lookup_data_smart(fcs, order, tol) 

32 _write_clusters(fname_clusters, cluster_lookup, order) 

33 _write_fc_lookup(fname_fc, fc_lookup, order) 

34 

35 

36def _write_fc_lookup(fname, fc_lookup, order): 

37 """ Writes the lookup force constants to file """ 

38 fmt = '{}' + ' {}'*order 

39 with open(fname, 'w') as f: 

40 f.write(str(len(fc_lookup)) + '\n\n') 

41 for fc in fc_lookup: 

42 for xyz in product(range(3), repeat=order): 

43 f.write(fmt.format(*xyz, fc[xyz])+'\n') 

44 f.write('\n') 

45 

46 

47def _write_clusters(fname, cluster_lookup, order): 

48 """ Writes the cluster lookup to file """ 

49 fmt = '{}' + ' {}'*order 

50 with open(fname, 'w') as f: 

51 f.write(str(len(cluster_lookup)) + '\n\n') 

52 for c, i in cluster_lookup.items(): 

53 line = fmt.format(*c, i) + '\n' 

54 f.write(line) 

55 

56 

57def _get_clusters(fcs, 

58 order: int, 

59 tol: float): 

60 """ Collect all relevant clusters; for 2nd and 3rd-order force constants 

61 all permutations are included. 

62 """ 

63 if order in [2, 3]: 

64 clusters = [] 

65 for c in fcs._fc_dict.keys(): 

66 if len(c) == order and np.linalg.norm(fcs[c]) > tol: 

67 for ci in permutations(c): 

68 clusters.append(ci) 

69 clusters = list(sorted(set(clusters))) 

70 else: 

71 clusters = [c for c in fcs._fc_dict.keys() if len(c) == order and np.linalg.norm(fcs[c]) > tol] # noqa 

72 return clusters 

73 

74 

75def _get_lookup_data_naive(fcs, 

76 order: int, 

77 tol: float): 

78 """ Groups force constants for a given order into groups for which the 

79 force constant is identical. """ 

80 fc_lookup = [] 

81 cluster_lookup = dict() 

82 

83 clusters = _get_clusters(fcs, order, tol) 

84 

85 for c in clusters: 

86 fc1 = fcs[c] 

87 if np.linalg.norm(fc1) < tol: 87 ↛ 88line 87 didn't jump to line 88 because the condition on line 87 was never true

88 continue 

89 for i, fc2 in enumerate(fc_lookup): 

90 if np.linalg.norm(fc1 - fc2) < tol: 

91 cluster_lookup[c] = i 

92 break 

93 else: 

94 cluster_lookup[c] = len(fc_lookup) 

95 fc_lookup.append(fc1) 

96 return cluster_lookup, fc_lookup 

97 

98 

99def _get_lookup_data_smart(fcs, 

100 order: int, 

101 tol: float): 

102 """ Groups force constants for a given order into groups for which the 

103 force constant is identical. """ 

104 fc_lookup = [] 

105 cluster_lookup = dict() 

106 axis = tuple(range(1, order+1)) 

107 

108 clusters = _get_clusters(fcs, order, tol) 

109 fc_all = np.array([fcs[c] for c in clusters]) 

110 

111 indices = list(range(len(clusters))) 

112 while len(indices) > 0: 

113 i = indices[0] 

114 delta = fc_all[indices] - fc_all[i] 

115 delta_norm = np.sqrt(np.sum(delta**2, axis=axis)) 

116 

117 inds_to_del = [indices[x] for x in np.where(delta_norm < tol)[0]] 

118 assert i in inds_to_del 

119 

120 fc_lookup.append(fc_all[i]) 

121 for j in inds_to_del: 

122 indices.remove(j) 

123 cluster_lookup[clusters[j]] = len(fc_lookup)-1 

124 return cluster_lookup, fc_lookup 

125 

126 

127def write_fcp_txt( 

128 fname: str, 

129 path: str, 

130 n_types: int, 

131 max_order: int, 

132 heat_current_order: int = 2, 

133) -> None: 

134 """ Write driver potential file for GPUMD. 

135 

136 Parameters 

137 ---------- 

138 fname 

139 Output file name. 

140 path 

141 Path to directory with force constant file. 

142 n_types 

143 Number of atom types. 

144 max_order 

145 Maximum order of the force constant potential. 

146 heat_current_order 

147 Heat current order used in thermal conductivity. 

148 

149 

150 Format is a simple file containing the following:: 

151 

152 fcp number_of_atom_types 

153 highest_force_order heat_current_order 

154 path_to_force_constant_files 

155 

156 which in practice for a binary system with a sixth order model, 

157 consider third-order heat-currents, would mean:: 

158 

159 fcp 2 

160 6 3 

161 /path/to/your/folder 

162 """ 

163 

164 with open(fname, 'w') as f: 

165 f.write('fcp {}\n'.format(n_types)) 

166 f.write('{} {}\n'.format(max_order, heat_current_order)) 

167 f.write('{}'.format(path.rstrip('/'))) # without a trailing '/' 

168 

169 

170def write_r0(fname: str, atoms: ase.Atoms): 

171 """ 

172 Write the GPUMD r0 file with the reference atomic positions. 

173 

174 Parameters 

175 ---------- 

176 fname 

177 Name of file to which to write the atomic positions. 

178 atoms 

179 Input structure. 

180 """ 

181 line = '{} {} {}\n' 

182 with open(fname, 'w') as f: 

183 for a in atoms: 

184 f.write(line.format(*a.position))