Coverage for /builds/debichem-team/python-ase/ase/calculators/turbomole/parameters.py: 39.12%

294 statements  

« prev     ^ index     » next       coverage.py v7.5.3, created at 2025-03-06 04:00 +0000

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 'non-define': True, 

95 'type': float, 

96 'units': None, 

97 'updateable': True 

98 }, 

99 'density functional': { 

100 'comment': None, 

101 'default': 'b-p', 

102 'group': 'dft', 

103 'key': 'functional', 

104 'type': str, 

105 'units': None, 

106 'updateable': True 

107 }, 

108 'energy convergence': { 

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

110 'default': None, 

111 'group': None, 

112 'key': None, 

113 'mapping': { 

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

115 'from_control': lambda a: a * Ha 

116 }, 

117 'type': float, 

118 'units': 'eV', 

119 'updateable': True 

120 }, 

121 'esp fit': { 

122 'comment': 'ESP fit', 

123 'default': None, 

124 'group': 'esp_fit', 

125 'key': 'esp_fit', 

126 'type': str, 

127 'units': None, 

128 'updateable': True, 

129 'non-define': True 

130 }, 

131 'fermi annealing factor': { 

132 'comment': None, 

133 'default': 0.95, 

134 'group': 'fermi', 

135 'key': 'tmfac', 

136 'type': float, 

137 'units': None, 

138 'updateable': True 

139 }, 

140 'fermi final temperature': { 

141 'comment': None, 

142 'default': 300., 

143 'group': 'fermi', 

144 'key': 'tmend', 

145 'type': float, 

146 'units': 'Kelvin', 

147 'updateable': True 

148 }, 

149 'fermi homo-lumo gap criterion': { 

150 'comment': None, 

151 'default': 0.1, 

152 'group': 'fermi', 

153 'key': 'hlcrt', 

154 'mapping': { 

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

156 'from_control': lambda a: a * Ha 

157 }, 

158 'type': float, 

159 'units': 'eV', 

160 'updateable': True 

161 }, 

162 'fermi initial temperature': { 

163 'comment': None, 

164 'default': 300., 

165 'group': 'fermi', 

166 'key': 'tmstrt', 

167 'type': float, 

168 'units': 'Kelvin', 

169 'updateable': True 

170 }, 

171 'fermi stopping criterion': { 

172 'comment': None, 

173 'default': 0.001, 

174 'group': 'fermi', 

175 'key': 'stop', 

176 'mapping': { 

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

178 'from_control': lambda a: a * Ha 

179 }, 

180 'type': float, 

181 'units': 'eV', 

182 'updateable': True 

183 }, 

184 'force convergence': { 

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

186 'default': None, 

187 'group': None, 

188 'key': None, 

189 'mapping': { 

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

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

192 }, 

193 'type': float, 

194 'units': 'eV/Angstrom', 

195 'updateable': True 

196 }, 

197 'geometry optimization iterations': { 

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

199 'default': None, 

200 'group': None, 

201 'key': None, 

202 'type': int, 

203 'units': None, 

204 'updateable': True 

205 }, 

206 'grid size': { 

207 'comment': None, 

208 'default': 'm3', 

209 'group': 'dft', 

210 'key': 'gridsize', 

211 'type': str, 

212 'units': None, 

213 'updateable': True 

214 }, 

215 'ground state': { 

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

217 'default': True, 

218 'group': None, 

219 'key': None, 

220 'type': bool, 

221 'units': None, 

222 'updateable': False 

223 }, 

224 'initial damping': { 

225 'comment': None, 

226 'default': None, 

227 'group': 'scfdamp', 

228 'key': 'start', 

229 'type': float, 

230 'units': None, 

231 'updateable': True 

232 }, 

233 'initial guess': { 

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

235 'default': 'eht', 

236 'group': None, 

237 'key': None, 

238 'type': (str, dict), 

239 'units': None, 

240 'updateable': False 

241 }, 

242 'minimal damping': { 

243 'comment': None, 

244 'default': None, 

245 'group': 'scfdamp', 

246 'key': 'min', 

247 'type': float, 

248 'units': None, 

249 'updateable': True 

250 }, 

251 'multiplicity': { 

252 'comment': None, 

253 'default': None, 

254 'group': None, 

255 'key': None, 

256 'type': int, 

257 'units': None, 

258 'updateable': False 

259 }, 

260 'non-automatic orbital shift': { 

261 'comment': None, 

262 'default': False, 

263 'group': 'scforbitalshift', 

264 'key': 'noautomatic', 

265 'type': bool, 

266 'units': None, 

267 'updateable': True 

268 }, 

269 'numerical hessian': { 

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

271 'default': None, 

272 'group': None, 

273 'key': None, 

274 'type': dict, 

275 'units': None, 

276 'updateable': True 

277 }, 

278 'point group': { 

279 'comment': 'only c1 supported', 

280 'default': 'c1', 

281 'group': 'symmetry', 

282 'key': 'symmetry', 

283 'type': str, 

284 'units': None, 

285 'updateable': False 

286 }, 

287 'ri memory': { 

288 'comment': None, 

289 'default': 1000, 

290 'group': 'ricore', 

291 'key': 'ricore', 

292 'type': int, 

293 'units': 'Megabyte', 

294 'updateable': True 

295 }, 

296 'rohf': { 

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

298 'default': None, 

299 'group': None, 

300 'key': None, 

301 'type': bool, 

302 'units': None, 

303 'updateable': False 

304 }, 

305 'scf energy convergence': { 

306 'comment': None, 

307 'default': None, 

308 'group': 'scfconv', 

309 'key': 'scfconv', 

310 'mapping': { 

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

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

313 }, 

314 'type': float, 

315 'units': 'eV', 

316 'updateable': True 

317 }, 

318 'scf iterations': { 

319 'comment': None, 

320 'default': 60, 

321 'group': 'scfiterlimit', 

322 'key': 'scfiterlimit', 

323 'type': int, 

324 'units': None, 

325 'updateable': True 

326 }, 

327 'task': { 

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

329 '"gradient calculation" = "gradient", ' 

330 '"geometry optimization" = "optimize", ' 

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

332 'default': 'energy', 

333 'group': None, 

334 'key': None, 

335 'type': str, 

336 'units': None, 

337 'updateable': True 

338 }, 

339 'title': { 

340 'comment': None, 

341 'default': '', 

342 'group': 'title', 

343 'key': 'title', 

344 'type': str, 

345 'units': None, 

346 'updateable': False 

347 }, 

348 'total charge': { 

349 'comment': None, 

350 'default': 0, 

351 'group': None, 

352 'key': None, 

353 'type': int, 

354 'units': None, 

355 'updateable': False 

356 }, 

357 'transition vector': { 

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

359 'default': None, 

360 'group': 'statpt', 

361 'key': 'itrvec', 

362 'type': int, 

363 'units': None, 

364 'updateable': True, 

365 'non-define': True 

366 }, 

367 'uhf': { 

368 'comment': None, 

369 'default': None, 

370 'group': 'uhf', 

371 'key': 'uhf', 

372 'type': bool, 

373 'units': None, 

374 'updateable': False 

375 }, 

376 'use basis set library': { 

377 'comment': 'only true implemented', 

378 'default': True, 

379 'group': 'basis', 

380 'key': None, 

381 'type': bool, 

382 'units': None, 

383 'updateable': False 

384 }, 

385 'use dft': { 

386 'comment': None, 

387 'default': True, 

388 'group': 'dft', 

389 'key': 'dft', 

390 'type': bool, 

391 'units': None, 

392 'updateable': False 

393 }, 

394 'use fermi smearing': { 

395 'comment': None, 

396 'default': False, 

397 'group': 'fermi', 

398 'key': 'fermi', 

399 'type': bool, 

400 'units': None, 

401 'updateable': True 

402 }, 

403 'use redundant internals': { 

404 'comment': None, 

405 'default': False, 

406 'group': 'redundant', 

407 'key': None, 

408 'type': bool, 

409 'units': None, 

410 'updateable': False 

411 }, 

412 'use resolution of identity': { 

413 'comment': None, 

414 'default': False, 

415 'group': 'rij', 

416 'key': 'rij', 

417 'type': bool, 

418 'units': None, 

419 'updateable': False 

420 } 

421 } 

422 

423 spec_names = { 

424 'default': 'default_parameters', 

425 'comment': 'parameter_comment', 

426 'updateable': 'parameter_updateable', 

427 'type': 'parameter_type', 

428 'key': 'parameter_key', 

429 'group': 'parameter_group', 

430 'units': 'parameter_units', 

431 'mapping': 'parameter_mapping', 

432 'non-define': 'parameter_no_define' 

433 } 

434 # flat dictionaries with parameters attributes 

435 default_parameters = {} 

436 parameter_group = {} 

437 parameter_type = {} 

438 parameter_key = {} 

439 parameter_units = {} 

440 parameter_comment = {} 

441 parameter_updateable = {} 

442 parameter_mapping = {} 

443 parameter_no_define = {} 

444 

445 def __init__(self, **kwargs): 

446 # construct flat dictionaries with parameter attributes 

447 for p in self.parameter_spec: 

448 for k in self.spec_names: 

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

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

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

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

453 self.update(kwargs) 

454 

455 def update(self, dct): 

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

457 for par in dct.keys(): 

458 if par not in self.parameter_spec: 

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

460 

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

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

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

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

465 raise TypeError(msg) 

466 self[key] = val 

467 

468 def update_data_groups(self, params_update): 

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

470 # construct a list of data groups to update 

471 grps = [] 

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

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

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

475 

476 # construct a dictionary of data groups and update params 

477 dgs = {} 

478 for g in grps: 

479 dgs[g] = {} 

480 for p in self.parameter_key: 

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

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

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

484 val = params_update[p] 

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

486 if val is not None and p in pmap: 

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

488 val = fun(params_update[p]) 

489 dgs[g] = val 

490 else: 

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

492 val = self.params_old[p] 

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

494 if val is not None and p in pmap: 

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

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

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

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

499 val = params_update[p] 

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

501 if val is not None and p in pmap: 

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

503 val = fun(params_update[p]) 

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

505 

506 # write dgs dictionary to a data group 

507 for g in dgs: 

508 delete_data_group(g) 

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

510 string = '' 

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

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

513 continue 

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

515 if dgs[g][key]: 

516 string += ' ' + key 

517 else: 

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

519 add_data_group(g, string=string) 

520 else: 

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

522 if dgs[g]: 

523 add_data_group(g, string='') 

524 else: 

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

526 

527 def update_no_define_parameters(self): 

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

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

530 if self.parameter_no_define.get(p): 

531 if v: 

532 if p in self.parameter_mapping: 

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

534 val = fun(v) 

535 else: 

536 val = v 

537 delete_data_group(self.parameter_group[p]) 

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

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

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

541 else: 

542 delete_data_group(self.parameter_group[p]) 

543 

544 def verify(self): 

545 """detect wrong or not implemented parameters""" 

546 

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

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

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

550 else: 

551 for par in self: 

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

553 

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

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

556 func = self['density functional'] 

557 assert func.lower() in func_list, ( 

558 'density functional not available / not supported' 

559 ) 

560 

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

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

563 

564 if self.get('rohf'): 

565 raise NotImplementedError('ROHF not implemented') 

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

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

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

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

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

571 raise NotImplementedError('Explicit basis set definition') 

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

573 raise NotImplementedError('Point group not impemeneted') 

574 

575 def get_define_str(self, natoms): 

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

577 

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

579 return self.define_str 

580 

581 define_str_tpl = ( 

582 '\n__title__\na coord\n__inter__\n' 

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

584 '__occ_str____single_atom_str____norb_str____dft_str____ri_str__' 

585 '__scfiterlimit____fermi_str____damp_str__q\n' 

586 ) 

587 

588 params = self 

589 

590 if params['use redundant internals']: 

591 internals_str = 'ired\n*' 

592 else: 

593 internals_str = '*\nno' 

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

595 

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

597 if params['uhf']: 

598 occ_str = 'n\ns\n*\n' 

599 else: 

600 occ_str = 'y\n' 

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

602 occ_str = 'y\n' 

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

604 occ_str = 'n\nt\n*\n' 

605 else: 

606 unpaired = params['multiplicity'] - 1 

607 if params['use fermi smearing']: 

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

609 else: 

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

611 

612 if natoms != 1: 

613 single_atom_str = '' 

614 else: 

615 single_atom_str = '\n' 

616 

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

618 norb_str = '' 

619 else: 

620 norb_str = 'n\n' 

621 

622 if params['use dft']: 

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

624 else: 

625 dft_str = '' 

626 

627 if params['density functional']: 

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

629 

630 if params['grid size']: 

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

632 

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

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

635 else: 

636 ri_str = '' 

637 

638 if params['scf iterations']: 

639 scfmaxiter = params['scf iterations'] 

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

641 else: 

642 scfiter_str = '' 

643 if params['scf energy convergence']: 

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

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

646 

647 fermi_str = '' 

648 if params['use fermi smearing']: 

649 fermi_str = 'scf\nfermi\n' 

650 if params['fermi initial temperature']: 

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

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

653 if params['fermi final temperature']: 

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

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

656 if params['fermi annealing factor']: 

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

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

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

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

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

662 if params['fermi stopping criterion']: 

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

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

665 fermi_str += '\n\n' 

666 

667 damp_str = '' 

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

669 'minimal damping') 

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

671 if any(damp_pars): 

672 damp_str = 'scf\ndamp\n' 

673 for par in damp_pars: 

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

675 damp_str += par_str + '\n' 

676 damp_str += '\n' 

677 

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

679 

680 define_str = define_str_tpl 

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

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

683 define_str) 

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

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

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

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

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

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

690 define_str) 

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

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

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

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

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

696 

697 return define_str 

698 

699 def read_restart(self, atoms, results): 

700 """read parameters from control file""" 

701 

702 params = {} 

703 pdgs = { 

704 p: parse_data_group( 

705 read_data_group( 

706 self.parameter_group[p]), self.parameter_group[p] 

707 ) 

708 for p in self.parameter_group 

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

710 } 

711 for p in self.parameter_key: 

712 if self.parameter_key[p]: 

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

714 if pdgs[p] is None: 

715 if self.parameter_type[p] is bool: 

716 params[p] = False 

717 else: 

718 params[p] = None 

719 else: 

720 if self.parameter_type[p] is bool: 

721 params[p] = True 

722 else: 

723 typ = self.parameter_type[p] 

724 val = typ(pdgs[p]) 

725 mapping = self.parameter_mapping 

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

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

728 val = fun(val) 

729 params[p] = val 

730 else: 

731 if pdgs[p] is None: 

732 params[p] = None 

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

734 if self.parameter_type[p] is bool: 

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

736 else: 

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

738 if self.parameter_type[p] is bool: 

739 params[p] = False 

740 else: 

741 params[p] = None 

742 else: 

743 typ = self.parameter_type[p] 

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

745 mapping = self.parameter_mapping 

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

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

748 val = fun(val) 

749 params[p] = val 

750 

751 # non-group or non-key parameters 

752 

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

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

755 assert len(basis_sets) == 1 

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

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

758 

759 # rohf, multiplicity and total charge 

760 orbs = results['molecular orbitals'] 

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

762 bool(len(read_data_group('roothaan')))) 

763 core_charge = 0 

764 if results['ecps']: 

765 for ecp in results['ecps']: 

766 for symbol in atoms.get_chemical_symbols(): 

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

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

769 if params['uhf']: 

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

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

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

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

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

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

776 electron_charge += core_charge 

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

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

779 params['multiplicity'] = 1 

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

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

782 electron_charge += core_charge 

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

784 else: 

785 raise NotImplementedError('ROHF not implemented') 

786 

787 # task-related parameters 

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

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

790 lines = log.readlines() 

791 for line in lines: 

792 if 'CRITERION FOR TOTAL SCF-ENERGY' in line: 

793 en = int(re.search(r'10\*{2}\((-\d+)\)', line).group(1)) 

794 mapp = self.parameter_mapping['energy convergence'] 

795 params['energy convergence'] = mapp['from_control'](10**en) 

796 if 'CRITERION FOR MAXIMUM NORM OF SCF-ENERGY GRADIENT' in line: 

797 gr = int(re.search(r'10\*{2}\((-\d+)\)', line).group(1)) 

798 mapp = self.parameter_mapping['force convergence'] 

799 params['force convergence'] = mapp['from_control'](10**gr) 

800 if 'AN OPTIMIZATION WITH MAX' in line: 

801 cy = int(re.search(r'MAX. (\d+) CYCLES', line).group(1)) 

802 params['geometry optimization iterations'] = cy 

803 self.update(params) 

804 self.params_old = params 

805 

806 def update_restart(self, dct): 

807 """update parameters after a restart""" 

808 nulst = [k for k in dct.keys() if not self.parameter_updateable[k]] 

809 if len(nulst) != 0: 

810 raise ValueError(f'parameters {nulst} cannot be changed') 

811 self.update(dct) 

812 self.update_data_groups(dct)