Coverage for /builds/kinetik161/ase/ase/calculators/aims.py: 39.67%

121 statements  

« prev     ^ index     » next       coverage.py v7.2.7, created at 2023-12-10 11:04 +0000

1"""This module defines an ASE interface to FHI-aims. 

2 

3Felix Hanke hanke@liverpool.ac.uk 

4Jonas Bjork j.bjork@liverpool.ac.uk 

5Simon P. Rittmeyer simon.rittmeyer@tum.de 

6 

7Edits on (24.11.2021) by Thomas A. R. Purcell purcell@fhi-berlin.mpg.de 

8""" 

9 

10import os 

11import re 

12 

13import numpy as np 

14 

15from ase.calculators.genericfileio import ( 

16 CalculatorTemplate, 

17 GenericFileIOCalculator, 

18 BaseProfile, 

19 read_stdout, 

20) 

21from ase.io.aims import write_aims, write_control 

22 

23 

24def get_aims_version(string): 

25 match = re.search(r'\s*FHI-aims version\s*:\s*(\S+)', string, re.M) 

26 return match.group(1) 

27 

28 

29class AimsProfile(BaseProfile): 

30 def __init__(self, binary, default_species_directory=None, **kwargs): 

31 super().__init__(**kwargs) 

32 self.binary = binary 

33 self.default_species_directory = default_species_directory 

34 

35 def get_calculator_command(self, inputfile): 

36 return [self.binary] 

37 

38 def version(self): 

39 return get_aims_version(read_stdout(self.binary)) 

40 

41 

42class AimsTemplate(CalculatorTemplate): 

43 def __init__(self): 

44 super().__init__( 

45 'aims', 

46 [ 

47 'energy', 

48 'free_energy', 

49 'forces', 

50 'stress', 

51 'stresses', 

52 'dipole', 

53 'magmom', 

54 ], 

55 ) 

56 

57 self.outputname = 'aims.out' 

58 

59 def update_parameters(self, properties, parameters): 

60 """Check and update the parameters to match the desired calculation 

61 

62 Parameters 

63 ---------- 

64 properties: list of str 

65 The list of properties to calculate 

66 parameters: dict 

67 The parameters used to perform the calculation. 

68 

69 Returns 

70 ------- 

71 dict 

72 The updated parameters object 

73 """ 

74 parameters = dict(parameters) 

75 property_flags = { 

76 'forces': 'compute_forces', 

77 'stress': 'compute_analytical_stress', 

78 'stresses': 'compute_heat_flux', 

79 } 

80 # Ensure FHI-aims will calculate all desired properties 

81 for property in properties: 

82 aims_name = property_flags.get(property, None) 

83 if aims_name is not None: 

84 parameters[aims_name] = True 

85 

86 if 'dipole' in properties: 

87 if 'output' in parameters and 'dipole' not in parameters['output']: 

88 parameters['output'] = list(parameters['output']) 

89 parameters['output'].append('dipole') 

90 elif 'output' not in parameters: 

91 parameters['output'] = ['dipole'] 

92 

93 return parameters 

94 

95 def write_input(self, profile, directory, atoms, parameters, properties): 

96 """Write the geometry.in and control.in files for the calculation 

97 

98 Parameters 

99 ---------- 

100 directory : Path 

101 The working directory to store the input files. 

102 atoms : atoms.Atoms 

103 The atoms object to perform the calculation on. 

104 parameters: dict 

105 The parameters used to perform the calculation. 

106 properties: list of str 

107 The list of properties to calculate 

108 """ 

109 parameters = self.update_parameters(properties, parameters) 

110 

111 ghosts = parameters.pop('ghosts', None) 

112 geo_constrain = parameters.pop('geo_constrain', None) 

113 scaled = parameters.pop('scaled', None) 

114 write_velocities = parameters.pop('write_velocities', None) 

115 

116 if scaled is None: 

117 scaled = np.all(atoms.pbc) 

118 if write_velocities is None: 

119 write_velocities = atoms.has('momenta') 

120 

121 if geo_constrain is None: 

122 geo_constrain = scaled and 'relax_geometry' in parameters 

123 

124 have_lattice_vectors = atoms.pbc.any() 

125 have_k_grid = ( 

126 'k_grid' in parameters 

127 or 'kpts' in parameters 

128 or 'k_grid_density' in parameters 

129 ) 

130 if have_lattice_vectors and not have_k_grid: 

131 raise RuntimeError('Found lattice vectors but no k-grid!') 

132 if not have_lattice_vectors and have_k_grid: 

133 raise RuntimeError('Found k-grid but no lattice vectors!') 

134 

135 geometry_in = directory / 'geometry.in' 

136 

137 write_aims( 

138 geometry_in, 

139 atoms, 

140 scaled, 

141 geo_constrain, 

142 write_velocities=write_velocities, 

143 ghosts=ghosts, 

144 ) 

145 

146 control = directory / 'control.in' 

147 

148 if ( 

149 'species_dir' not in parameters 

150 and profile.default_species_directory is not None 

151 ): 

152 parameters['species_dir'] = profile.default_species_directory 

153 

154 write_control(control, atoms, parameters) 

155 

156 def execute(self, directory, profile): 

157 profile.run(directory, inputfile=None, outputfile=self.outputname) 

158 

159 def read_results(self, directory): 

160 from ase.io.aims import read_aims_results 

161 

162 dst = directory / self.outputname 

163 return read_aims_results(dst, index=-1) 

164 

165 def load_profile(self, cfg, **kwargs): 

166 return AimsProfile.from_config(cfg, self.name, **kwargs) 

167 

168 def socketio_argv(self, profile, unixsocket, port): 

169 return [profile.binary] 

170 

171 def socketio_parameters(self, unixsocket, port): 

172 if port: 

173 use_pimd_wrapper = ('localhost', port) 

174 else: 

175 # (INET port number should be unused.) 

176 use_pimd_wrapper = (f'UNIX:{unixsocket}', 31415) 

177 

178 return dict(use_pimd_wrapper=use_pimd_wrapper, compute_forces=True) 

179 

180 

181class Aims(GenericFileIOCalculator): 

182 def __init__( 

183 self, 

184 profile=None, 

185 directory='.', 

186 parallel_info=None, 

187 parallel=True, 

188 **kwargs, 

189 ): 

190 """Construct the FHI-aims calculator. 

191 

192 The keyword arguments (kwargs) can be one of the ASE standard 

193 keywords: 'xc', 'kpts' and 'smearing' or any of FHI-aims' 

194 native keywords. 

195 

196 

197 Arguments: 

198 

199 cubes: AimsCube object 

200 Cube file specification. 

201 

202 tier: int or array of ints 

203 Set basis set tier for all atomic species. 

204 

205 plus_u : dict 

206 For DFT+U. Adds a +U term to one specific shell of the species. 

207 

208 kwargs : dict 

209 Any of the base class arguments. 

210 

211 """ 

212 

213 super().__init__( 

214 template=AimsTemplate(), 

215 profile=profile, 

216 parameters=kwargs, 

217 parallel_info=parallel_info, 

218 parallel=parallel, 

219 directory=directory, 

220 ) 

221 

222 

223class AimsCube: 

224 'Object to ensure the output of cube files, can be attached to Aims object' 

225 

226 def __init__( 

227 self, 

228 origin=(0, 0, 0), 

229 edges=[(0.1, 0.0, 0.0), (0.0, 0.1, 0.0), (0.0, 0.0, 0.1)], 

230 points=(50, 50, 50), 

231 plots=(), 

232 ): 

233 """parameters: 

234 

235 origin, edges, points: 

236 Same as in the FHI-aims output 

237 plots: 

238 what to print, same names as in FHI-aims""" 

239 

240 self.name = 'AimsCube' 

241 self.origin = origin 

242 self.edges = edges 

243 self.points = points 

244 self.plots = plots 

245 

246 def ncubes(self): 

247 """returns the number of cube files to output""" 

248 return len(self.plots) 

249 

250 def move_to_base_name(self, basename): 

251 """when output tracking is on or the base namem is not standard, 

252 this routine will rename add the base to the cube file output for 

253 easier tracking""" 

254 for plot in self.plots: 

255 found = False 

256 cube = plot.split() 

257 if ( 

258 cube[0] == 'total_density' 

259 or cube[0] == 'spin_density' 

260 or cube[0] == 'delta_density' 

261 ): 

262 found = True 

263 old_name = cube[0] + '.cube' 

264 new_name = basename + '.' + old_name 

265 if cube[0] == 'eigenstate' or cube[0] == 'eigenstate_density': 

266 found = True 

267 state = int(cube[1]) 

268 s_state = cube[1] 

269 for i in [10, 100, 1000, 10000]: 

270 if state < i: 

271 s_state = '0' + s_state 

272 old_name = cube[0] + '_' + s_state + '_spin_1.cube' 

273 new_name = basename + '.' + old_name 

274 if found: 

275 # XXX Should not use platform dependent commands! 

276 os.system('mv ' + old_name + ' ' + new_name) 

277 

278 def add_plot(self, name): 

279 """in case you forgot one ...""" 

280 self.plots += [name] 

281 

282 def write(self, file): 

283 """write the necessary output to the already opened control.in""" 

284 file.write('output cube ' + self.plots[0] + '\n') 

285 file.write(' cube origin ') 

286 for ival in self.origin: 

287 file.write(str(ival) + ' ') 

288 file.write('\n') 

289 for i in range(3): 

290 file.write(' cube edge ' + str(self.points[i]) + ' ') 

291 for ival in self.edges[i]: 

292 file.write(str(ival) + ' ') 

293 file.write('\n') 

294 if self.ncubes() > 1: 

295 for i in range(self.ncubes() - 1): 

296 file.write('output cube ' + self.plots[i + 1] + '\n')