Coverage for /builds/ase/ase/ase/io/dmol.py: 94.16%
154 statements
« prev ^ index » next coverage.py v7.5.3, created at 2025-08-02 00:12 +0000
1# fmt: off
3"""
4IO functions for DMol3 file formats.
6read/write functionality for car, incoor and arc file formats
7only car format is added to known ase file extensions
8use format='dmol-arc' or 'dmol-incoor' for others
10car structure file - Angstrom and cellpar description of cell.
11incoor structure file - Bohr and cellvector describption of cell.
12 Note: incoor file not used if car file present.
13arc multiple-structure file - Angstrom and cellpar description of cell.
16The formats follow strict formatting
18car
19----
20col: 1-5 atom name
21col: 7-20 x Cartesian coordinate of atom in A
22col: 22-35 y Cartesian coordinate of atom in A
23col: 37-50 z Cartesian coordinate of atom in A
24col: 52-55 type of residue containing atom
25col: 57-63 residue sequence name relative to beginning of current molecule,
26 left justified
27col: 64-70 potential type of atom left justified
28col: 72-73 element symbol
29col: 75-80 partial charge on atom
32incoor
33-------
34$cell vectors
35 37.83609647462165 0.00000000000000 0.00000000000000
36 0.00000000000000 37.60366016124745 0.00000000000000
37 0.00000000000000 0.00000000000000 25.29020473078921
38$coordinates
39Si 15.94182672614820 1.85274838936809 16.01426481346124
40Si 4.45559370448989 2.68957177851318 -0.05326937257442
41$end
44arc
45----
46multiple images of car format separated with $end
49"""
51from datetime import datetime
53import numpy as np
55from ase import Atom, Atoms
56from ase.geometry.cell import cell_to_cellpar, cellpar_to_cell
57from ase.units import Bohr
58from ase.utils import reader, writer
61@writer
62def write_dmol_car(fd, atoms):
63 """ Write a dmol car-file from an Atoms object
65 Notes
66 -----
67 The positions written to file are rotated as to align with the cell when
68 reading (due to cellpar information)
69 Can not handle multiple images.
70 Only allows for pbc 111 or 000.
71 """
73 fd.write('!BIOSYM archive 3\n')
74 dt = datetime.now()
76 symbols = atoms.get_chemical_symbols()
77 if np.all(atoms.pbc):
78 # Rotate positions so they will align with cellpar cell
79 cellpar = cell_to_cellpar(atoms.cell)
80 new_cell = cellpar_to_cell(cellpar)
81 lstsq_fit = np.linalg.lstsq(atoms.cell, new_cell, rcond=-1)
82 # rcond=-1 silences FutureWarning in numpy 1.14
83 R = lstsq_fit[0]
84 positions = np.dot(atoms.positions, R)
86 fd.write('PBC=ON\n\n')
87 fd.write('!DATE %s\n' % dt.strftime('%b %d %H:%M:%S %Y'))
88 fd.write('PBC %9.5f %9.5f %9.5f %9.5f %9.5f %9.5f\n' % tuple(cellpar))
89 elif not np.any(atoms.pbc): # [False,False,False]
90 fd.write('PBC=OFF\n\n')
91 fd.write('!DATE %s\n' % dt.strftime('%b %d %H:%M:%S %Y'))
92 positions = atoms.positions
93 else:
94 raise ValueError('PBC must be all true or all false for .car format')
96 for i, (sym, pos) in enumerate(zip(symbols, positions)):
97 fd.write('%-6s %12.8f %12.8f %12.8f XXXX 1 xx %-2s '
98 '0.000\n' % (sym + str(i + 1), pos[0], pos[1], pos[2], sym))
99 fd.write('end\nend\n')
102@reader
103def read_dmol_car(fd):
104 """ Read a dmol car-file and return an Atoms object.
106 Notes
107 -----
108 Cell is constructed from cellpar so orientation of cell might be off.
109 """
111 lines = fd.readlines()
112 atoms = Atoms()
114 start_line = 4
116 if lines[1][4:6] == 'ON':
117 start_line += 1
118 cell_dat = np.array([float(fld) for fld in lines[4].split()[1:7]])
119 cell = cellpar_to_cell(cell_dat)
120 pbc = [True, True, True]
121 else:
122 cell = np.zeros((3, 3))
123 pbc = [False, False, False]
125 symbols = []
126 positions = []
127 for line in lines[start_line:]:
128 if line.startswith('end'):
129 break
130 flds = line.split()
131 symbols.append(flds[7])
132 positions.append(flds[1:4])
133 atoms.append(Atom(flds[7], flds[1:4]))
134 atoms = Atoms(symbols=symbols, positions=positions, cell=cell, pbc=pbc)
135 return atoms
138@writer
139def write_dmol_incoor(fd, atoms, bohr=True):
140 """ Write a dmol incoor-file from an Atoms object
142 Notes
143 -----
144 Only used for pbc 111.
145 Can not handle multiple images.
146 DMol3 expect data in .incoor files to be in bohr, if bohr is false however
147 the data is written in Angstroms.
148 """
150 if not np.all(atoms.pbc):
151 raise ValueError('PBC must be all true for .incoor format')
153 if bohr:
154 cell = atoms.cell / Bohr
155 positions = atoms.positions / Bohr
156 else:
157 cell = atoms.cell
158 positions = atoms.positions
160 fd.write('$cell vectors\n')
161 fd.write(' {:18.14f} {:18.14f} {:18.14f}\n'.format(
162 cell[0, 0], cell[0, 1], cell[0, 2]))
163 fd.write(' {:18.14f} {:18.14f} {:18.14f}\n'.format(
164 cell[1, 0], cell[1, 1], cell[1, 2]))
165 fd.write(' {:18.14f} {:18.14f} {:18.14f}\n'.format(
166 cell[2, 0], cell[2, 1], cell[2, 2]))
168 fd.write('$coordinates\n')
169 for a, pos in zip(atoms, positions):
170 fd.write('%-12s%18.14f %18.14f %18.14f \n' % (
171 a.symbol, pos[0], pos[1], pos[2]))
172 fd.write('$end\n')
175@reader
176def read_dmol_incoor(fd, bohr=True):
177 """ Reads an incoor file and returns an atoms object.
179 Notes
180 -----
181 If bohr is True then incoor is assumed to be in bohr and the data
182 is rescaled to Angstrom.
183 """
185 lines = fd.readlines()
186 symbols = []
187 positions = []
188 for i, line in enumerate(lines):
189 if line.startswith('$cell vectors'):
190 cell = np.zeros((3, 3))
191 for j, line in enumerate(lines[i + 1:i + 4]):
192 cell[j, :] = [float(fld) for fld in line.split()]
193 if line.startswith('$coordinates'):
194 j = i + 1
195 while True:
196 if lines[j].startswith('$end'):
197 break
198 flds = lines[j].split()
199 symbols.append(flds[0])
200 positions.append(flds[1:4])
201 j += 1
202 atoms = Atoms(symbols=symbols, positions=positions, cell=cell, pbc=True)
203 if bohr:
204 atoms.cell = atoms.cell * Bohr
205 atoms.positions = atoms.positions * Bohr
206 return atoms
209@writer
210def write_dmol_arc(fd, images):
211 """ Writes all images to file filename in arc format.
213 Similar to the .car format only pbc 111 or 000 is supported.
214 """
216 fd.write('!BIOSYM archive 3\n')
217 if np.all(images[0].pbc):
218 fd.write('PBC=ON\n\n')
219 # Rotate positions so they will align with cellpar cell
220 elif not np.any(images[0].pbc):
221 fd.write('PBC=OFF\n\n')
222 else:
223 raise ValueError('PBC must be all true or all false for .arc format')
224 for atoms in images:
225 dt = datetime.now()
226 symbols = atoms.get_chemical_symbols()
227 if np.all(atoms.pbc):
228 cellpar = cell_to_cellpar(atoms.cell)
229 new_cell = cellpar_to_cell(cellpar)
230 lstsq_fit = np.linalg.lstsq(atoms.cell, new_cell, rcond=-1)
231 R = lstsq_fit[0]
232 fd.write('!DATE %s\n' % dt.strftime('%b %d %H:%m:%S %Y'))
233 fd.write('PBC %9.5f %9.5f %9.5f %9.5f %9.5f %9.5f\n'
234 % tuple(cellpar))
235 positions = np.dot(atoms.positions, R)
236 elif not np.any(atoms.pbc): # [False,False,False]
237 fd.write('!DATE %s\n' % dt.strftime('%b %d %H:%m:%S %Y'))
238 positions = atoms.positions
239 else:
240 raise ValueError(
241 'PBC must be all true or all false for .arc format')
242 for i, (sym, pos) in enumerate(zip(symbols, positions)):
243 fd.write(
244 '%-6s %12.8f %12.8f %12.8f XXXX 1 xx %-2s '
245 '0.000\n' % (sym + str(i + 1), pos[0], pos[1], pos[2], sym))
246 fd.write('end\nend\n')
247 fd.write('\n')
250@reader
251def read_dmol_arc(fd, index=-1):
252 """ Read a dmol arc-file and return a series of Atoms objects (images). """
254 lines = fd.readlines()
255 images = []
257 if lines[1].startswith('PBC=ON'):
258 pbc = True
259 elif lines[1].startswith('PBC=OFF'):
260 pbc = False
261 else:
262 raise RuntimeError('Could not read pbc from second line in file')
264 i = 0
265 while i < len(lines):
266 cell = np.zeros((3, 3))
267 symbols = []
268 positions = []
269 # parse single image
270 if lines[i].startswith('!DATE'):
271 # read cell
272 if pbc:
273 cell_dat = np.array([float(fld)
274 for fld in lines[i + 1].split()[1:7]])
275 cell = cellpar_to_cell(cell_dat)
276 i += 1
277 i += 1
278 # read atoms
279 while not lines[i].startswith('end'):
280 flds = lines[i].split()
281 symbols.append(flds[7])
282 positions.append(flds[1:4])
283 i += 1
284 image = Atoms(symbols=symbols, positions=positions, cell=cell,
285 pbc=pbc)
286 images.append(image)
287 if len(images) == index:
288 return images[-1]
289 i += 1
291 # return requested images, code borrowed from ase/io/trajectory.py
292 if isinstance(index, int):
293 return images[index]
294 else:
295 from ase.io.formats import index2range
296 indices = index2range(index, len(images))
297 return [images[j] for j in indices]