Coverage for /builds/kinetik161/ase/ase/calculators/turbomole/parameters.py: 38.93%

1# type: ignore 

2"""turbomole parameters management classes and functions""" 

3 

4import os 

5import re 

6from math import floor, log10 

7 

8import numpy as np 

9 

10from ase.calculators.turbomole.reader import parse_data_group, read_data_group 

11from ase.calculators.turbomole.writer import add_data_group, delete_data_group 

12from ase.units import Bohr, Ha 

13 

14 

15class TurbomoleParameters(dict): 

16 """class to manage turbomole parameters""" 

17 

18 available_functionals = [ 

19 'slater-dirac-exchange', 's-vwn', 'vwn', 's-vwn_Gaussian', 'pwlda', 

20 'becke-exchange', 'b-lyp', 'b-vwn', 'lyp', 'b-p', 'pbe', 'tpss', 

21 'bh-lyp', 'b3-lyp', 'b3-lyp_Gaussian', 'pbe0', 'tpssh', 'lhf', 'oep', 

22 'b97-d', 'b2-plyp' 

23 ] 

24 

25 # nested dictionary with parameters attributes 

26 parameter_spec = { 

27 'automatic orbital shift': { 

28 'comment': None, 

29 'default': 0.1, 

30 'group': 'scforbitalshift', 

31 'key': 'automatic', 

32 'mapping': { 

33 'to_control': lambda a: a / Ha, 

34 'from_control': lambda a: a * Ha 

35 }, 

36 'type': float, 

37 'units': 'eV', 

38 'updateable': True 

39 }, 

40 'basis set definition': { 

41 'comment': 'used only in restart', 

42 'default': None, 

43 'group': 'basis', 

44 'key': None, 

45 'type': list, 

46 'units': None, 

47 'updateable': False 

48 }, 

49 'basis set name': { 

50 'comment': 'current default from module "define"', 

51 'default': 'def-SV(P)', 

52 'group': 'basis', 

53 'key': None, 

54 'type': str, 

55 'units': None, 

56 'updateable': False 

57 }, 

58 'closed-shell orbital shift': { 

59 'comment': 'does not work with automatic', 

60 'default': None, 

61 'group': 'scforbitalshift', 

62 'key': 'closedshell', 

63 'mapping': { 

64 'to_control': lambda a: a / Ha, 

65 'from_control': lambda a: a * Ha 

66 }, 

67 'type': float, 

68 'units': 'eV', 

69 'updateable': True 

70 }, 

71 'damping adjustment step': { 

72 'comment': None, 

73 'default': None, 

74 'group': 'scfdamp', 

75 'key': 'step', 

76 'type': float, 

77 'units': None, 

78 'updateable': True 

79 }, 

80 'default eht atomic orbitals': { 

81 'comment': None, 

82 'default': None, 

83 'group': None, 

84 'key': None, 

85 'type': bool, 

86 'units': None, 

87 'updateable': False 

88 }, 

89 'density convergence': { 

90 'comment': None, 

91 'default': None, 

92 'group': 'denconv', 

93 'key': 'denconv', 

94 'mapping': { 

95 'to_control': lambda a: int(-log10(a)), 

96 'from_control': lambda a: 10**(-a) 

97 }, 

98 'non-define': True, 

99 'type': float, 

100 'units': None, 

101 'updateable': True 

102 }, 

103 'density functional': { 

104 'comment': None, 

105 'default': 'b-p', 

106 'group': 'dft', 

107 'key': 'functional', 

108 'type': str, 

109 'units': None, 

110 'updateable': True 

111 }, 

112 'energy convergence': { 

113 'comment': 'jobex -energy <int>', 

114 'default': None, 

115 'group': None, 

116 'key': None, 

117 'mapping': { 

118 'to_control': lambda a: a / Ha, 

119 'from_control': lambda a: a * Ha 

120 }, 

121 'type': float, 

122 'units': 'eV', 

123 'updateable': True 

124 }, 

125 'esp fit': { 

126 'comment': 'ESP fit', 

127 'default': None, 

128 'group': 'esp_fit', 

129 'key': 'esp_fit', 

130 'type': str, 

131 'units': None, 

132 'updateable': True, 

133 'non-define': True 

134 }, 

135 'fermi annealing factor': { 

136 'comment': None, 

137 'default': 0.95, 

138 'group': 'fermi', 

139 'key': 'tmfac', 

140 'type': float, 

141 'units': None, 

142 'updateable': True 

143 }, 

144 'fermi final temperature': { 

145 'comment': None, 

146 'default': 300., 

147 'group': 'fermi', 

148 'key': 'tmend', 

149 'type': float, 

150 'units': 'Kelvin', 

151 'updateable': True 

152 }, 

153 'fermi homo-lumo gap criterion': { 

154 'comment': None, 

155 'default': 0.1, 

156 'group': 'fermi', 

157 'key': 'hlcrt', 

158 'mapping': { 

159 'to_control': lambda a: a / Ha, 

160 'from_control': lambda a: a * Ha 

161 }, 

162 'type': float, 

163 'units': 'eV', 

164 'updateable': True 

165 }, 

166 'fermi initial temperature': { 

167 'comment': None, 

168 'default': 300., 

169 'group': 'fermi', 

170 'key': 'tmstrt', 

171 'type': float, 

172 'units': 'Kelvin', 

173 'updateable': True 

174 }, 

175 'fermi stopping criterion': { 

176 'comment': None, 

177 'default': 0.001, 

178 'group': 'fermi', 

179 'key': 'stop', 

180 'mapping': { 

181 'to_control': lambda a: a / Ha, 

182 'from_control': lambda a: a * Ha 

183 }, 

184 'type': float, 

185 'units': 'eV', 

186 'updateable': True 

187 }, 

188 'force convergence': { 

189 'comment': 'jobex -gcart <int>', 

190 'default': None, 

191 'group': None, 

192 'key': None, 

193 'mapping': { 

194 'to_control': lambda a: a / Ha * Bohr, 

195 'from_control': lambda a: a * Ha / Bohr 

196 }, 

197 'type': float, 

198 'units': 'eV/Angstrom', 

199 'updateable': True 

200 }, 

201 'geometry optimization iterations': { 

202 'comment': 'jobex -c <int>', 

203 'default': None, 

204 'group': None, 

205 'key': None, 

206 'type': int, 

207 'units': None, 

208 'updateable': True 

209 }, 

210 'grid size': { 

211 'comment': None, 

212 'default': 'm3', 

213 'group': 'dft', 

214 'key': 'gridsize', 

215 'type': str, 

216 'units': None, 

217 'updateable': True 

218 }, 

219 'ground state': { 

220 'comment': 'only this is currently supported', 

221 'default': True, 

222 'group': None, 

223 'key': None, 

224 'type': bool, 

225 'units': None, 

226 'updateable': False 

227 }, 

228 'initial damping': { 

229 'comment': None, 

230 'default': None, 

231 'group': 'scfdamp', 

232 'key': 'start', 

233 'type': float, 

234 'units': None, 

235 'updateable': True 

236 }, 

237 'initial guess': { 

238 'comment': '"eht", "hcore" or {"use": "<path/to/control>"}', 

239 'default': 'eht', 

240 'group': None, 

241 'key': None, 

242 'type': (str, dict), 

243 'units': None, 

244 'updateable': False 

245 }, 

246 'minimal damping': { 

247 'comment': None, 

248 'default': None, 

249 'group': 'scfdamp', 

250 'key': 'min', 

251 'type': float, 

252 'units': None, 

253 'updateable': True 

254 }, 

255 'multiplicity': { 

256 'comment': None, 

257 'default': None, 

258 'group': None, 

259 'key': None, 

260 'type': int, 

261 'units': None, 

262 'updateable': False 

263 }, 

264 'non-automatic orbital shift': { 

265 'comment': None, 

266 'default': False, 

267 'group': 'scforbitalshift', 

268 'key': 'noautomatic', 

269 'type': bool, 

270 'units': None, 

271 'updateable': True 

272 }, 

273 'numerical hessian': { 

274 'comment': 'NumForce will be used if dictionary exists', 

275 'default': None, 

276 'group': None, 

277 'key': None, 

278 'type': dict, 

279 'units': None, 

280 'updateable': True 

281 }, 

282 'point group': { 

283 'comment': 'only c1 supported', 

284 'default': 'c1', 

285 'group': 'symmetry', 

286 'key': 'symmetry', 

287 'type': str, 

288 'units': None, 

289 'updateable': False 

290 }, 

291 'ri memory': { 

292 'comment': None, 

293 'default': 1000, 

294 'group': 'ricore', 

295 'key': 'ricore', 

296 'type': int, 

297 'units': 'Megabyte', 

298 'updateable': True 

299 }, 

300 'rohf': { 

301 'comment': 'used only in restart', 

302 'default': None, 

303 'group': None, 

304 'key': None, 

305 'type': bool, 

306 'units': None, 

307 'updateable': False 

308 }, 

309 'scf energy convergence': { 

310 'comment': None, 

311 'default': None, 

312 'group': 'scfconv', 

313 'key': 'scfconv', 

314 'mapping': { 

315 'to_control': lambda a: int(floor(-log10(a / Ha))), 

316 'from_control': lambda a: 10**(-a) * Ha 

317 }, 

318 'type': float, 

319 'units': 'eV', 

320 'updateable': True 

321 }, 

322 'scf iterations': { 

323 'comment': None, 

324 'default': 60, 

325 'group': 'scfiterlimit', 

326 'key': 'scfiterlimit', 

327 'type': int, 

328 'units': None, 

329 'updateable': True 

330 }, 

331 'task': { 

332 'comment': '"energy calculation" = "energy", ' 

333 '"gradient calculation" = "gradient", ' 

334 '"geometry optimization" = "optimize", ' 

335 '"normal mode analysis" = "frequencies"', 

336 'default': 'energy', 

337 'group': None, 

338 'key': None, 

339 'type': str, 

340 'units': None, 

341 'updateable': True 

342 }, 

343 'title': { 

344 'comment': None, 

345 'default': '', 

346 'group': 'title', 

347 'key': 'title', 

348 'type': str, 

349 'units': None, 

350 'updateable': False 

351 }, 

352 'total charge': { 

353 'comment': None, 

354 'default': 0, 

355 'group': None, 

356 'key': None, 

357 'type': int, 

358 'units': None, 

359 'updateable': False 

360 }, 

361 'transition vector': { 

362 'comment': 'vector for transition state optimization', 

363 'default': None, 

364 'group': 'statpt', 

365 'key': 'itrvec', 

366 'type': int, 

367 'units': None, 

368 'updateable': True, 

369 'non-define': True 

370 }, 

371 'uhf': { 

372 'comment': None, 

373 'default': None, 

374 'group': 'uhf', 

375 'key': 'uhf', 

376 'type': bool, 

377 'units': None, 

378 'updateable': False 

379 }, 

380 'use basis set library': { 

381 'comment': 'only true implemented', 

382 'default': True, 

383 'group': 'basis', 

384 'key': None, 

385 'type': bool, 

386 'units': None, 

387 'updateable': False 

388 }, 

389 'use dft': { 

390 'comment': None, 

391 'default': True, 

392 'group': 'dft', 

393 'key': 'dft', 

394 'type': bool, 

395 'units': None, 

396 'updateable': False 

397 }, 

398 'use fermi smearing': { 

399 'comment': None, 

400 'default': False, 

401 'group': 'fermi', 

402 'key': 'fermi', 

403 'type': bool, 

404 'units': None, 

405 'updateable': True 

406 }, 

407 'use redundant internals': { 

408 'comment': None, 

409 'default': False, 

410 'group': 'redundant', 

411 'key': None, 

412 'type': bool, 

413 'units': None, 

414 'updateable': False 

415 }, 

416 'use resolution of identity': { 

417 'comment': None, 

418 'default': False, 

419 'group': 'rij', 

420 'key': 'rij', 

421 'type': bool, 

422 'units': None, 

423 'updateable': False 

424 } 

425 } 

426 

427 spec_names = { 

428 'default': 'default_parameters', 

429 'comment': 'parameter_comment', 

430 'updateable': 'parameter_updateable', 

431 'type': 'parameter_type', 

432 'key': 'parameter_key', 

433 'group': 'parameter_group', 

434 'units': 'parameter_units', 

435 'mapping': 'parameter_mapping', 

436 'non-define': 'parameter_no_define' 

437 } 

438 # flat dictionaries with parameters attributes 

439 default_parameters = {} 

440 parameter_group = {} 

441 parameter_type = {} 

442 parameter_key = {} 

443 parameter_units = {} 

444 parameter_comment = {} 

445 parameter_updateable = {} 

446 parameter_mapping = {} 

447 parameter_no_define = {} 

448 

449 def __init__(self, **kwargs): 

450 # construct flat dictionaries with parameter attributes 

451 for p in self.parameter_spec: 

452 for k in self.spec_names: 

453 if k in list(self.parameter_spec[p].keys()): 

454 subdict = getattr(self, self.spec_names[k]) 

455 subdict.update({p: self.parameter_spec[p][k]}) 

456 super().__init__(**self.default_parameters) 

457 self.update(kwargs) 

458 

459 def update(self, dct): 

460 """check the type of parameters in dct and then update""" 

461 for par in dct.keys(): 

462 if par not in self.parameter_spec: 

463 raise ValueError('invalid parameter: ' + par) 

464 

465 for key, val in dct.items(): 

466 correct_type = self.parameter_spec[key]['type'] 

467 if not isinstance(val, (correct_type, type(None))): 

468 msg = str(key) + ' has wrong type: ' + str(type(val)) 

469 raise TypeError(msg) 

470 self[key] = val 

471 

472 def update_data_groups(self, params_update): 

473 """updates data groups in the control file""" 

474 # construct a list of data groups to update 

475 grps = [] 

476 for p in list(params_update.keys()): 

477 if self.parameter_group[p] is not None: 

478 grps.append(self.parameter_group[p]) 

479 

480 # construct a dictionary of data groups and update params 

481 dgs = {} 

482 for g in grps: 

483 dgs[g] = {} 

484 for p in self.parameter_key: 

485 if g == self.parameter_group[p]: 

486 if self.parameter_group[p] == self.parameter_key[p]: 

487 if p in list(params_update.keys()): 

488 val = params_update[p] 

489 pmap = list(self.parameter_mapping.keys()) 

490 if val is not None and p in pmap: 

491 fun = self.parameter_mapping[p]['to_control'] 

492 val = fun(params_update[p]) 

493 dgs[g] = val 

494 else: 

495 if p in list(self.params_old.keys()): 

496 val = self.params_old[p] 

497 pmap = list(self.parameter_mapping.keys()) 

498 if val is not None and p in pmap: 

499 fun = self.parameter_mapping[p]['to_control'] 

500 val = fun(self.params_old[p]) 

501 dgs[g][self.parameter_key[p]] = val 

502 if p in list(params_update.keys()): 

503 val = params_update[p] 

504 pmap = list(self.parameter_mapping.keys()) 

505 if val is not None and p in pmap: 

506 fun = self.parameter_mapping[p]['to_control'] 

507 val = fun(params_update[p]) 

508 dgs[g][self.parameter_key[p]] = val 

509 

510 # write dgs dictionary to a data group 

511 for g in dgs: 

512 delete_data_group(g) 

513 if isinstance(dgs[g], dict): 

514 string = '' 

515 for key in list(dgs[g].keys()): 

516 if dgs[g][key] is None: 

517 continue 

518 elif isinstance(dgs[g][key], bool): 

519 if dgs[g][key]: 

520 string += ' ' + key 

521 else: 

522 string += ' ' + key + '=' + str(dgs[g][key]) 

523 add_data_group(g, string=string) 

524 else: 

525 if isinstance(dgs[g], bool): 

526 if dgs[g]: 

527 add_data_group(g, string='') 

528 else: 

529 add_data_group(g, string=str(dgs[g])) 

530 

531 def update_no_define_parameters(self): 

532 """process key parameters that are not written with define""" 

533 for p, v in self.items(): 

534 if p in self.parameter_no_define and self.parameter_no_define[p]: 

535 if v: 

536 if p in self.parameter_mapping: 

537 fun = self.parameter_mapping[p]['to_control'] 

538 val = fun(v) 

539 else: 

540 val = v 

541 delete_data_group(self.parameter_group[p]) 

542 if self.parameter_group[p] != self.parameter_key[p]: 

543 val = '\n ' + self.parameter_key[p] + ' ' + str(val) 

544 add_data_group(self.parameter_group[p], str(val)) 

545 else: 

546 delete_data_group(self.parameter_group[p]) 

547 

548 def verify(self): 

549 """detect wrong or not implemented parameters""" 

550 

551 if getattr(self, 'define_str', None) is not None: 

552 assert isinstance(self.define_str, str), 'define_str must be str' 

553 assert len(self.define_str) != 0, 'define_str may not be empty' 

554 else: 

555 for par in self: 

556 assert par in self.parameter_spec, 'invalid parameter: ' + par 

557 

558 if self.get('use dft'): 

559 func_list = [x.lower() for x in self.available_functionals] 

560 func = self['density functional'] 

561 assert func.lower() in func_list, ( 

562 'density functional not available / not supported' 

563 ) 

564 

565 assert self['multiplicity'] is not None, 'multiplicity not defined' 

566 assert self['multiplicity'] > 0, 'multiplicity has wrong value' 

567 

568 if self.get('rohf'): 

569 raise NotImplementedError('ROHF not implemented') 

570 if self['initial guess'] not in ['eht', 'hcore']: 

571 if not (isinstance(self['initial guess'], dict) and 

572 'use' in self['initial guess'].keys()): 

573 raise ValueError('Wrong input for initial guess') 

574 if not self['use basis set library']: 

575 raise NotImplementedError('Explicit basis set definition') 

576 if self['point group'] != 'c1': 

577 raise NotImplementedError('Point group not impemeneted') 

578 

579 def get_define_str(self, natoms): 

580 """construct a define string from the parameters dictionary""" 

581 

582 if getattr(self, 'define_str', None): 

583 return self.define_str 

584 

585 define_str_tpl = ( 

586 '\n__title__\na coord\n__inter__\n' 

587 'bb all __basis_set__\n*\neht\n__eht_aos_str__y\n__charge_str__' 

588 '__occ_str____single_atom_str____norb_str____dft_str____ri_str__' 

589 '__scfiterlimit____fermi_str____damp_str__q\n' 

590 ) 

591 

592 params = self 

593 

594 if params['use redundant internals']: 

595 internals_str = 'ired\n*' 

596 else: 

597 internals_str = '*\nno' 

598 charge_str = str(params['total charge']) + '\n' 

599 

600 if params['multiplicity'] == 1: 

601 if params['uhf']: 

602 occ_str = 'n\ns\n*\n' 

603 else: 

604 occ_str = 'y\n' 

605 elif params['multiplicity'] == 2: 

606 occ_str = 'y\n' 

607 elif params['multiplicity'] == 3: 

608 occ_str = 'n\nt\n*\n' 

609 else: 

610 unpaired = params['multiplicity'] - 1 

611 if params['use fermi smearing']: 

612 occ_str = 'n\nuf ' + str(unpaired) + '\n*\n' 

613 else: 

614 occ_str = 'n\nu ' + str(unpaired) + '\n*\n' 

615 

616 if natoms != 1: 

617 single_atom_str = '' 

618 else: 

619 single_atom_str = '\n' 

620 

621 if params['multiplicity'] == 1 and not params['uhf']: 

622 norb_str = '' 

623 else: 

624 norb_str = 'n\n' 

625 

626 if params['use dft']: 

627 dft_str = 'dft\non\n*\n' 

628 else: 

629 dft_str = '' 

630 

631 if params['density functional']: 

632 dft_str += 'dft\nfunc ' + params['density functional'] + '\n*\n' 

633 

634 if params['grid size']: 

635 dft_str += 'dft\ngrid ' + params['grid size'] + '\n*\n' 

636 

637 if params['use resolution of identity']: 

638 ri_str = 'ri\non\nm ' + str(params['ri memory']) + '\n*\n' 

639 else: 

640 ri_str = '' 

641 

642 if params['scf iterations']: 

643 scfmaxiter = params['scf iterations'] 

644 scfiter_str = 'scf\niter\n' + str(scfmaxiter) + '\n\n' 

645 else: 

646 scfiter_str = '' 

647 if params['scf energy convergence']: 

648 conv = floor(-log10(params['scf energy convergence'] / Ha)) 

649 scfiter_str += 'scf\nconv\n' + str(int(conv)) + '\n\n' 

650 

651 fermi_str = '' 

652 if params['use fermi smearing']: 

653 fermi_str = 'scf\nfermi\n' 

654 if params['fermi initial temperature']: 

655 par = str(params['fermi initial temperature']) 

656 fermi_str += '1\n' + par + '\n' 

657 if params['fermi final temperature']: 

658 par = str(params['fermi final temperature']) 

659 fermi_str += '2\n' + par + '\n' 

660 if params['fermi annealing factor']: 

661 par = str(params['fermi annealing factor']) 

662 fermi_str += '3\n' + par + '\n' 

663 if params['fermi homo-lumo gap criterion']: 

664 par = str(params['fermi homo-lumo gap criterion']) 

665 fermi_str += '4\n' + par + '\n' 

666 if params['fermi stopping criterion']: 

667 par = str(params['fermi stopping criterion']) 

668 fermi_str += '5\n' + par + '\n' 

669 fermi_str += '\n\n' 

670 

671 damp_str = '' 

672 damp_keys = ('initial damping', 'damping adjustment step', 

673 'minimal damping') 

674 damp_pars = [params[k] for k in damp_keys] 

675 if any(damp_pars): 

676 damp_str = 'scf\ndamp\n' 

677 for par in damp_pars: 

678 par_str = str(par) if par else '' 

679 damp_str += par_str + '\n' 

680 damp_str += '\n' 

681 

682 eht_aos_str = 'y\n' if params['default eht atomic orbitals'] else '' 

683 

684 define_str = define_str_tpl 

685 define_str = re.sub('__title__', params['title'], define_str) 

686 define_str = re.sub('__basis_set__', params['basis set name'], 

687 define_str) 

688 define_str = re.sub('__charge_str__', charge_str, define_str) 

689 define_str = re.sub('__occ_str__', occ_str, define_str) 

690 define_str = re.sub('__norb_str__', norb_str, define_str) 

691 define_str = re.sub('__dft_str__', dft_str, define_str) 

692 define_str = re.sub('__ri_str__', ri_str, define_str) 

693 define_str = re.sub('__single_atom_str__', single_atom_str, 

694 define_str) 

695 define_str = re.sub('__inter__', internals_str, define_str) 

696 define_str = re.sub('__scfiterlimit__', scfiter_str, define_str) 

697 define_str = re.sub('__fermi_str__', fermi_str, define_str) 

698 define_str = re.sub('__damp_str__', damp_str, define_str) 

699 define_str = re.sub('__eht_aos_str__', eht_aos_str, define_str) 

700 

701 return define_str 

702 

703 def read_restart(self, atoms, results): 

704 """read parameters from control file""" 

705 

706 params = {} 

707 pdgs = {} 

708 for p in self.parameter_group: 

709 if self.parameter_group[p] and self.parameter_key[p]: 

710 pdgs[p] = parse_data_group( 

711 read_data_group(self.parameter_group[p]), 

712 self.parameter_group[p] 

713 ) 

714 

715 for p in self.parameter_key: 

716 if self.parameter_key[p]: 

717 if self.parameter_key[p] == self.parameter_group[p]: 

718 if pdgs[p] is None: 

719 if self.parameter_type[p] is bool: 

720 params[p] = False 

721 else: 

722 params[p] = None 

723 else: 

724 if self.parameter_type[p] is bool: 

725 params[p] = True 

726 else: 

727 typ = self.parameter_type[p] 

728 val = typ(pdgs[p]) 

729 mapping = self.parameter_mapping 

730 if p in list(mapping.keys()): 

731 fun = mapping[p]['from_control'] 

732 val = fun(val) 

733 params[p] = val 

734 else: 

735 if pdgs[p] is None: 

736 params[p] = None 

737 elif isinstance(pdgs[p], str): 

738 if self.parameter_type[p] is bool: 

739 params[p] = (pdgs[p] == self.parameter_key[p]) 

740 else: 

741 if self.parameter_key[p] not in list(pdgs[p].keys()): 

742 if self.parameter_type[p] is bool: 

743 params[p] = False 

744 else: 

745 params[p] = None 

746 else: 

747 typ = self.parameter_type[p] 

748 val = typ(pdgs[p][self.parameter_key[p]]) 

749 mapping = self.parameter_mapping 

750 if p in list(mapping.keys()): 

751 fun = mapping[p]['from_control'] 

752 val = fun(val) 

753 params[p] = val 

754 

755 # non-group or non-key parameters 

756 

757 # per-element and per-atom basis sets not implemented in calculator 

758 basis_sets = {bs['nickname'] for bs in results['basis set']} 

759 assert len(basis_sets) == 1 

760 params['basis set name'] = list(basis_sets)[0] 

761 params['basis set definition'] = results['basis set'] 

762 

763 # rohf, multiplicity and total charge 

764 orbs = results['molecular orbitals'] 

765 params['rohf'] = (bool(len(read_data_group('rohf'))) or 

766 bool(len(read_data_group('roothaan')))) 

767 core_charge = 0 

768 if results['ecps']: 

769 for ecp in results['ecps']: 

770 for symbol in atoms.get_chemical_symbols(): 

771 if symbol.lower() == ecp['element'].lower(): 

772 core_charge -= ecp['number of core electrons'] 

773 if params['uhf']: 

774 alpha_occ = [o['occupancy'] for o in orbs if o['spin'] == 'alpha'] 

775 beta_occ = [o['occupancy'] for o in orbs if o['spin'] == 'beta'] 

776 spin = (np.sum(alpha_occ) - np.sum(beta_occ)) * 0.5 

777 params['multiplicity'] = int(2 * spin + 1) 

778 nuclear_charge = int(sum(atoms.numbers)) 

779 electron_charge = -int(sum(alpha_occ) + sum(beta_occ)) 

780 electron_charge += core_charge 

781 params['total charge'] = nuclear_charge + electron_charge 

782 elif not params['rohf']: # restricted HF (closed shell) 

783 params['multiplicity'] = 1 

784 nuclear_charge = int(sum(atoms.numbers)) 

785 electron_charge = -int(sum(o['occupancy'] for o in orbs)) 

786 electron_charge += core_charge 

787 params['total charge'] = nuclear_charge + electron_charge 

788 else: 

789 raise NotImplementedError('ROHF not implemented') 

790 

791 # task-related parameters 

792 if os.path.exists('job.start'): 

793 with open('job.start') as log: 

794 lines = log.readlines()