Coverage for /builds/kinetik161/ase/ase/io/dmol.py: 94.16%

1""" 

2IO functions for DMol3 file formats. 

3 

4read/write functionality for car, incoor and arc file formats 

5only car format is added to known ase file extensions 

6use format='dmol-arc' or 'dmol-incoor' for others 

7 

8car structure file - Angstrom and cellpar description of cell. 

9incoor structure file - Bohr and cellvector describption of cell. 

10 Note: incoor file not used if car file present. 

11arc multiple-structure file - Angstrom and cellpar description of cell. 

12 

13 

14The formats follow strict formatting 

15 

16car 

17---- 

18col: 1-5 atom name 

19col: 7-20 x Cartesian coordinate of atom in A 

20col: 22-35 y Cartesian coordinate of atom in A 

21col: 37-50 z Cartesian coordinate of atom in A 

22col: 52-55 type of residue containing atom 

23col: 57-63 residue sequence name relative to beginning of current molecule, 

24 left justified 

25col: 64-70 potential type of atom left justified 

26col: 72-73 element symbol 

27col: 75-80 partial charge on atom 

28 

29 

30incoor 

31------- 

32$cell vectors 

33 37.83609647462165 0.00000000000000 0.00000000000000 

34 0.00000000000000 37.60366016124745 0.00000000000000 

35 0.00000000000000 0.00000000000000 25.29020473078921 

36$coordinates 

37Si 15.94182672614820 1.85274838936809 16.01426481346124 

38Si 4.45559370448989 2.68957177851318 -0.05326937257442 

39$end 

40 

41 

42arc 

43---- 

44multiple images of car format separated with $end 

45 

46 

47""" 

48 

49from datetime import datetime 

50 

51import numpy as np 

52 

53from ase import Atom, Atoms 

54from ase.geometry.cell import cell_to_cellpar, cellpar_to_cell 

55from ase.units import Bohr 

56from ase.utils import reader, writer 

57 

58 

59@writer 

60def write_dmol_car(fd, atoms): 

61 """ Write a dmol car-file from an Atoms object 

62 

63 Notes 

64 ----- 

65 The positions written to file are rotated as to align with the cell when 

66 reading (due to cellpar information) 

67 Can not handle multiple images. 

68 Only allows for pbc 111 or 000. 

69 """ 

70 

71 fd.write('!BIOSYM archive 3\n') 

72 dt = datetime.now() 

73 

74 symbols = atoms.get_chemical_symbols() 

75 if np.all(atoms.pbc): 

76 # Rotate positions so they will align with cellpar cell 

77 cellpar = cell_to_cellpar(atoms.cell) 

78 new_cell = cellpar_to_cell(cellpar) 

79 lstsq_fit = np.linalg.lstsq(atoms.cell, new_cell, rcond=-1) 

80 # rcond=-1 silences FutureWarning in numpy 1.14 

81 R = lstsq_fit[0] 

82 positions = np.dot(atoms.positions, R) 

83 

84 fd.write('PBC=ON\n\n') 

85 fd.write('!DATE %s\n' % dt.strftime('%b %d %H:%M:%S %Y')) 

86 fd.write('PBC %9.5f %9.5f %9.5f %9.5f %9.5f %9.5f\n' % tuple(cellpar)) 

87 elif not np.any(atoms.pbc): # [False,False,False] 

88 fd.write('PBC=OFF\n\n') 

89 fd.write('!DATE %s\n' % dt.strftime('%b %d %H:%M:%S %Y')) 

90 positions = atoms.positions 

91 else: 

92 raise ValueError('PBC must be all true or all false for .car format') 

93 

94 for i, (sym, pos) in enumerate(zip(symbols, positions)): 

95 fd.write('%-6s %12.8f %12.8f %12.8f XXXX 1 xx %-2s ' 

96 '0.000\n' % (sym + str(i + 1), pos[0], pos[1], pos[2], sym)) 

97 fd.write('end\nend\n') 

98 

99 

100@reader 

101def read_dmol_car(fd): 

102 """ Read a dmol car-file and return an Atoms object. 

103 

104 Notes 

105 ----- 

106 Cell is constructed from cellpar so orientation of cell might be off. 

107 """ 

108 

109 lines = fd.readlines() 

110 atoms = Atoms() 

111 

112 start_line = 4 

113 

114 if lines[1][4:6] == 'ON': 

115 start_line += 1 

116 cell_dat = np.array([float(fld) for fld in lines[4].split()[1:7]]) 

117 cell = cellpar_to_cell(cell_dat) 

118 pbc = [True, True, True] 

119 else: 

120 cell = np.zeros((3, 3)) 

121 pbc = [False, False, False] 

122 

123 symbols = [] 

124 positions = [] 

125 for line in lines[start_line:]: 

126 if line.startswith('end'): 

127 break 

128 flds = line.split() 

129 symbols.append(flds[7]) 

130 positions.append(flds[1:4]) 

131 atoms.append(Atom(flds[7], flds[1:4])) 

132 atoms = Atoms(symbols=symbols, positions=positions, cell=cell, pbc=pbc) 

133 return atoms 

134 

135 

136@writer 

137def write_dmol_incoor(fd, atoms, bohr=True): 

138 """ Write a dmol incoor-file from an Atoms object 

139 

140 Notes 

141 ----- 

142 Only used for pbc 111. 

143 Can not handle multiple images. 

144 DMol3 expect data in .incoor files to be in bohr, if bohr is false however 

145 the data is written in Angstroms. 

146 """ 

147 

148 if not np.all(atoms.pbc): 

149 raise ValueError('PBC must be all true for .incoor format') 

150 

151 if bohr: 

152 cell = atoms.cell / Bohr 

153 positions = atoms.positions / Bohr 

154 else: 

155 cell = atoms.cell 

156 positions = atoms.positions 

157 

158 fd.write('$cell vectors\n') 

159 fd.write(' {:18.14f} {:18.14f} {:18.14f}\n'.format( 

160 cell[0, 0], cell[0, 1], cell[0, 2])) 

161 fd.write(' {:18.14f} {:18.14f} {:18.14f}\n'.format( 

162 cell[1, 0], cell[1, 1], cell[1, 2])) 

163 fd.write(' {:18.14f} {:18.14f} {:18.14f}\n'.format( 

164 cell[2, 0], cell[2, 1], cell[2, 2])) 

165 

166 fd.write('$coordinates\n') 

167 for a, pos in zip(atoms, positions): 

168 fd.write('%-12s%18.14f %18.14f %18.14f \n' % ( 

169 a.symbol, pos[0], pos[1], pos[2])) 

170 fd.write('$end\n') 

171 

172 

173@reader 

174def read_dmol_incoor(fd, bohr=True): 

175 """ Reads an incoor file and returns an atoms object. 

176 

177 Notes 

178 ----- 

179 If bohr is True then incoor is assumed to be in bohr and the data 

180 is rescaled to Angstrom. 

181 """ 

182 

183 lines = fd.readlines() 

184 symbols = [] 

185 positions = [] 

186 for i, line in enumerate(lines): 

187 if line.startswith('$cell vectors'): 

188 cell = np.zeros((3, 3)) 

189 for j, line in enumerate(lines[i + 1:i + 4]): 

190 cell[j, :] = [float(fld) for fld in line.split()] 

191 if line.startswith('$coordinates'): 

192 j = i + 1 

193 while True: 

194 if lines[j].startswith('$end'): 

195 break 

196 flds = lines[j].split() 

197 symbols.append(flds[0]) 

198 positions.append(flds[1:4]) 

199 j += 1 

200 atoms = Atoms(symbols=symbols, positions=positions, cell=cell, pbc=True) 

201 if bohr: 

202 atoms.cell = atoms.cell * Bohr 

203 atoms.positions = atoms.positions * Bohr 

204 return atoms 

205 

206 

207@writer 

208def write_dmol_arc(fd, images): 

209 """ Writes all images to file filename in arc format. 

210 

211 Similar to the .car format only pbc 111 or 000 is supported. 

212 """ 

213 

214 fd.write('!BIOSYM archive 3\n') 

215 if np.all(images[0].pbc): 

216 fd.write('PBC=ON\n\n') 

217 # Rotate positions so they will align with cellpar cell 

218 elif not np.any(images[0].pbc): 

219 fd.write('PBC=OFF\n\n') 

220 else: 

221 raise ValueError('PBC must be all true or all false for .arc format') 

222 for atoms in images: 

223 dt = datetime.now() 

224 symbols = atoms.get_chemical_symbols() 

225 if np.all(atoms.pbc): 

226 cellpar = cell_to_cellpar(atoms.cell) 

227 new_cell = cellpar_to_cell(cellpar) 

228 lstsq_fit = np.linalg.lstsq(atoms.cell, new_cell, rcond=-1) 

229 R = lstsq_fit[0] 

230 fd.write('!DATE %s\n' % dt.strftime('%b %d %H:%m:%S %Y')) 

231 fd.write('PBC %9.5f %9.5f %9.5f %9.5f %9.5f %9.5f\n' 

232 % tuple(cellpar)) 

233 positions = np.dot(atoms.positions, R) 

234 elif not np.any(atoms.pbc): # [False,False,False] 

235 fd.write('!DATE %s\n' % dt.strftime('%b %d %H:%m:%S %Y')) 

236 positions = atoms.positions 

237 else: 

238 raise ValueError( 

239 'PBC must be all true or all false for .arc format') 

240 for i, (sym, pos) in enumerate(zip(symbols, positions)): 

241 fd.write( 

242 '%-6s %12.8f %12.8f %12.8f XXXX 1 xx %-2s ' 

243 '0.000\n' % (sym + str(i + 1), pos[0], pos[1], pos[2], sym)) 

244 fd.write('end\nend\n') 

245 fd.write('\n') 

246 

247 

248@reader 

249def read_dmol_arc(fd, index=-1): 

250 """ Read a dmol arc-file and return a series of Atoms objects (images). """ 

251 

252 lines = fd.readlines() 

253 images = [] 

254 

255 if lines[1].startswith('PBC=ON'): 

256 pbc = True 

257 elif lines[1].startswith('PBC=OFF'): 

258 pbc = False 

259 else: 

260 raise RuntimeError('Could not read pbc from second line in file') 

261 

262 i = 0 

263 while i < len(lines): 

264 cell = np.zeros((3, 3)) 

265 symbols = [] 

266 positions = [] 

267 # parse single image 

268 if lines[i].startswith('!DATE'): 

269 # read cell 

270 if pbc: 

271 cell_dat = np.array([float(fld) 

272 for fld in lines[i + 1].split()[1:7]]) 

273 cell = cellpar_to_cell(cell_dat) 

274 i += 1 

275 i += 1 

276 # read atoms 

277 while not lines[i].startswith('end'): 

278 flds = lines[i].split() 

279 symbols.append(flds[7]) 

280 positions.append(flds[1:4]) 

281 i += 1 

282 image = Atoms(symbols=symbols, positions=positions, cell=cell, 

283 pbc=pbc) 

284 images.append(image) 

285 if len(images) == index: 

286 return images[-1] 

287 i += 1 

288 

289 # return requested images, code borrowed from ase/io/trajectory.py 

290 if isinstance(index, int): 

291 return images[index] 

292 else: 

293 from ase.io.formats import index2range 

294 indices = index2range(index, len(images)) 

295 return [images[j] for j in indices]