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

294 statements  

« prev     ^ index     » next       coverage.py v7.5.3, created at 2025-08-02 00:12 +0000

1# fmt: off 

2 

3# type: ignore 

4"""turbomole parameters management classes and functions""" 

5 

6import os 

7import re 

8from math import floor, log10 

9 

10import numpy as np 

11 

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

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

14from ase.units import Bohr, Ha 

15 

16 

17class TurbomoleParameters(dict): 

18 """class to manage turbomole parameters""" 

19 

20 available_functionals = [ 

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

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

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

24 'b97-d', 'b2-plyp' 

25 ] 

26 

27 # nested dictionary with parameters attributes 

28 parameter_spec = { 

29 'automatic orbital shift': { 

30 'comment': None, 

31 'default': 0.1, 

32 'group': 'scforbitalshift', 

33 'key': 'automatic', 

34 'mapping': { 

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

36 'from_control': lambda a: a * Ha 

37 }, 

38 'type': float, 

39 'units': 'eV', 

40 'updateable': True 

41 }, 

42 'basis set definition': { 

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

44 'default': None, 

45 'group': 'basis', 

46 'key': None, 

47 'type': list, 

48 'units': None, 

49 'updateable': False 

50 }, 

51 'basis set name': { 

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

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

54 'group': 'basis', 

55 'key': None, 

56 'type': str, 

57 'units': None, 

58 'updateable': False 

59 }, 

60 'closed-shell orbital shift': { 

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

62 'default': None, 

63 'group': 'scforbitalshift', 

64 'key': 'closedshell', 

65 'mapping': { 

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

67 'from_control': lambda a: a * Ha 

68 }, 

69 'type': float, 

70 'units': 'eV', 

71 'updateable': True 

72 }, 

73 'damping adjustment step': { 

74 'comment': None, 

75 'default': None, 

76 'group': 'scfdamp', 

77 'key': 'step', 

78 'type': float, 

79 'units': None, 

80 'updateable': True 

81 }, 

82 'default eht atomic orbitals': { 

83 'comment': None, 

84 'default': None, 

85 'group': None, 

86 'key': None, 

87 'type': bool, 

88 'units': None, 

89 'updateable': False 

90 }, 

91 'density convergence': { 

92 'comment': None, 

93 'default': None, 

94 'group': 'denconv', 

95 'key': 'denconv', 

96 'non-define': True, 

97 'type': float, 

98 'units': None, 

99 'updateable': True 

100 }, 

101 'density functional': { 

102 'comment': None, 

103 'default': 'b-p', 

104 'group': 'dft', 

105 'key': 'functional', 

106 'type': str, 

107 'units': None, 

108 'updateable': True 

109 }, 

110 'energy convergence': { 

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

112 'default': None, 

113 'group': None, 

114 'key': None, 

115 'mapping': { 

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

117 'from_control': lambda a: a * Ha 

118 }, 

119 'type': float, 

120 'units': 'eV', 

121 'updateable': True 

122 }, 

123 'esp fit': { 

124 'comment': 'ESP fit', 

125 'default': None, 

126 'group': 'esp_fit', 

127 'key': 'esp_fit', 

128 'type': str, 

129 'units': None, 

130 'updateable': True, 

131 'non-define': True 

132 }, 

133 'fermi annealing factor': { 

134 'comment': None, 

135 'default': 0.95, 

136 'group': 'fermi', 

137 'key': 'tmfac', 

138 'type': float, 

139 'units': None, 

140 'updateable': True 

141 }, 

142 'fermi final temperature': { 

143 'comment': None, 

144 'default': 300., 

145 'group': 'fermi', 

146 'key': 'tmend', 

147 'type': float, 

148 'units': 'Kelvin', 

149 'updateable': True 

150 }, 

151 'fermi homo-lumo gap criterion': { 

152 'comment': None, 

153 'default': 0.1, 

154 'group': 'fermi', 

155 'key': 'hlcrt', 

156 'mapping': { 

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

158 'from_control': lambda a: a * Ha 

159 }, 

160 'type': float, 

161 'units': 'eV', 

162 'updateable': True 

163 }, 

164 'fermi initial temperature': { 

165 'comment': None, 

166 'default': 300., 

167 'group': 'fermi', 

168 'key': 'tmstrt', 

169 'type': float, 

170 'units': 'Kelvin', 

171 'updateable': True 

172 }, 

173 'fermi stopping criterion': { 

174 'comment': None, 

175 'default': 0.001, 

176 'group': 'fermi', 

177 'key': 'stop', 

178 'mapping': { 

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

180 'from_control': lambda a: a * Ha 

181 }, 

182 'type': float, 

183 'units': 'eV', 

184 'updateable': True 

185 }, 

186 'force convergence': { 

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

188 'default': None, 

189 'group': None, 

190 'key': None, 

191 'mapping': { 

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

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

194 }, 

195 'type': float, 

196 'units': 'eV/Angstrom', 

197 'updateable': True 

198 }, 

199 'geometry optimization iterations': { 

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

201 'default': None, 

202 'group': None, 

203 'key': None, 

204 'type': int, 

205 'units': None, 

206 'updateable': True 

207 }, 

208 'grid size': { 

209 'comment': None, 

210 'default': 'm3', 

211 'group': 'dft', 

212 'key': 'gridsize', 

213 'type': str, 

214 'units': None, 

215 'updateable': True 

216 }, 

217 'ground state': { 

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

219 'default': True, 

220 'group': None, 

221 'key': None, 

222 'type': bool, 

223 'units': None, 

224 'updateable': False 

225 }, 

226 'initial damping': { 

227 'comment': None, 

228 'default': None, 

229 'group': 'scfdamp', 

230 'key': 'start', 

231 'type': float, 

232 'units': None, 

233 'updateable': True 

234 }, 

235 'initial guess': { 

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

237 'default': 'eht', 

238 'group': None, 

239 'key': None, 

240 'type': (str, dict), 

241 'units': None, 

242 'updateable': False 

243 }, 

244 'minimal damping': { 

245 'comment': None, 

246 'default': None, 

247 'group': 'scfdamp', 

248 'key': 'min', 

249 'type': float, 

250 'units': None, 

251 'updateable': True 

252 }, 

253 'multiplicity': { 

254 'comment': None, 

255 'default': None, 

256 'group': None, 

257 'key': None, 

258 'type': int, 

259 'units': None, 

260 'updateable': False 

261 }, 

262 'non-automatic orbital shift': { 

263 'comment': None, 

264 'default': False, 

265 'group': 'scforbitalshift', 

266 'key': 'noautomatic', 

267 'type': bool, 

268 'units': None, 

269 'updateable': True 

270 }, 

271 'numerical hessian': { 

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

273 'default': None, 

274 'group': None, 

275 'key': None, 

276 'type': dict, 

277 'units': None, 

278 'updateable': True 

279 }, 

280 'point group': { 

281 'comment': 'only c1 supported', 

282 'default': 'c1', 

283 'group': 'symmetry', 

284 'key': 'symmetry', 

285 'type': str, 

286 'units': None, 

287 'updateable': False 

288 }, 

289 'ri memory': { 

290 'comment': None, 

291 'default': 1000, 

292 'group': 'ricore', 

293 'key': 'ricore', 

294 'type': int, 

295 'units': 'Megabyte', 

296 'updateable': True 

297 }, 

298 'rohf': { 

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

300 'default': None, 

301 'group': None, 

302 'key': None, 

303 'type': bool, 

304 'units': None, 

305 'updateable': False 

306 }, 

307 'scf energy convergence': { 

308 'comment': None, 

309 'default': None, 

310 'group': 'scfconv', 

311 'key': 'scfconv', 

312 'mapping': { 

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

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

315 }, 

316 'type': float, 

317 'units': 'eV', 

318 'updateable': True 

319 }, 

320 'scf iterations': { 

321 'comment': None, 

322 'default': 60, 

323 'group': 'scfiterlimit', 

324 'key': 'scfiterlimit', 

325 'type': int, 

326 'units': None, 

327 'updateable': True 

328 }, 

329 'task': { 

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

331 '"gradient calculation" = "gradient", ' 

332 '"geometry optimization" = "optimize", ' 

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

334 'default': 'energy', 

335 'group': None, 

336 'key': None, 

337 'type': str, 

338 'units': None, 

339 'updateable': True 

340 }, 

341 'title': { 

342 'comment': None, 

343 'default': '', 

344 'group': 'title', 

345 'key': 'title', 

346 'type': str, 

347 'units': None, 

348 'updateable': False 

349 }, 

350 'total charge': { 

351 'comment': None, 

352 'default': 0, 

353 'group': None, 

354 'key': None, 

355 'type': int, 

356 'units': None, 

357 'updateable': False 

358 }, 

359 'transition vector': { 

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

361 'default': None, 

362 'group': 'statpt', 

363 'key': 'itrvec', 

364 'type': int, 

365 'units': None, 

366 'updateable': True, 

367 'non-define': True 

368 }, 

369 'uhf': { 

370 'comment': None, 

371 'default': None, 

372 'group': 'uhf', 

373 'key': 'uhf', 

374 'type': bool, 

375 'units': None, 

376 'updateable': False 

377 }, 

378 'use basis set library': { 

379 'comment': 'only true implemented', 

380 'default': True, 

381 'group': 'basis', 

382 'key': None, 

383 'type': bool, 

384 'units': None, 

385 'updateable': False 

386 }, 

387 'use dft': { 

388 'comment': None, 

389 'default': True, 

390 'group': 'dft', 

391 'key': 'dft', 

392 'type': bool, 

393 'units': None, 

394 'updateable': False 

395 }, 

396 'use fermi smearing': { 

397 'comment': None, 

398 'default': False, 

399 'group': 'fermi', 

400 'key': 'fermi', 

401 'type': bool, 

402 'units': None, 

403 'updateable': True 

404 }, 

405 'use redundant internals': { 

406 'comment': None, 

407 'default': False, 

408 'group': 'redundant', 

409 'key': None, 

410 'type': bool, 

411 'units': None, 

412 'updateable': False 

413 }, 

414 'use resolution of identity': { 

415 'comment': None, 

416 'default': False, 

417 'group': 'rij', 

418 'key': 'rij', 

419 'type': bool, 

420 'units': None, 

421 'updateable': False 

422 } 

423 } 

424 

425 spec_names = { 

426 'default': 'default_parameters', 

427 'comment': 'parameter_comment', 

428 'updateable': 'parameter_updateable', 

429 'type': 'parameter_type', 

430 'key': 'parameter_key', 

431 'group': 'parameter_group', 

432 'units': 'parameter_units', 

433 'mapping': 'parameter_mapping', 

434 'non-define': 'parameter_no_define' 

435 } 

436 # flat dictionaries with parameters attributes 

437 default_parameters = {} 

438 parameter_group = {} 

439 parameter_type = {} 

440 parameter_key = {} 

441 parameter_units = {} 

442 parameter_comment = {} 

443 parameter_updateable = {} 

444 parameter_mapping = {} 

445 parameter_no_define = {} 

446 

447 def __init__(self, **kwargs): 

448 # construct flat dictionaries with parameter attributes 

449 for p in self.parameter_spec: 

450 for k in self.spec_names: 

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

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

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

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

455 self.update(kwargs) 

456 

457 def update(self, dct): 

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

459 for par in dct.keys(): 

460 if par not in self.parameter_spec: 

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

462 

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

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

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

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

467 raise TypeError(msg) 

468 self[key] = val 

469 

470 def update_data_groups(self, params_update): 

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

472 # construct a list of data groups to update 

473 grps = [] 

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

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

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

477 

478 # construct a dictionary of data groups and update params 

479 dgs = {} 

480 for g in grps: 

481 dgs[g] = {} 

482 for p in self.parameter_key: 

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

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

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

486 val = params_update[p] 

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

488 if val is not None and p in pmap: 

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

490 val = fun(params_update[p]) 

491 dgs[g] = val 

492 else: 

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

494 val = self.params_old[p] 

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

496 if val is not None and p in pmap: 

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

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

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

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

501 val = params_update[p] 

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

503 if val is not None and p in pmap: 

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

505 val = fun(params_update[p]) 

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

507 

508 # write dgs dictionary to a data group 

509 for g in dgs: 

510 delete_data_group(g) 

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

512 string = '' 

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

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

515 continue 

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

517 if dgs[g][key]: 

518 string += ' ' + key 

519 else: 

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

521 add_data_group(g, string=string) 

522 else: 

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

524 if dgs[g]: 

525 add_data_group(g, string='') 

526 else: 

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

528 

529 def update_no_define_parameters(self): 

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

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

532 if self.parameter_no_define.get(p): 

533 if v: 

534 if p in self.parameter_mapping: 

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

536 val = fun(v) 

537 else: 

538 val = v 

539 delete_data_group(self.parameter_group[p]) 

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

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

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

543 else: 

544 delete_data_group(self.parameter_group[p]) 

545 

546 def verify(self): 

547 """detect wrong or not implemented parameters""" 

548 

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

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

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

552 else: 

553 for par in self: 

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

555 

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

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

558 func = self['density functional'] 

559 assert func.lower() in func_list, ( 

560 'density functional not available / not supported' 

561 ) 

562 

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

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

565 

566 if self.get('rohf'): 

567 raise NotImplementedError('ROHF not implemented') 

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

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

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

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

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

573 raise NotImplementedError('Explicit basis set definition') 

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

575 raise NotImplementedError('Point group not impemeneted') 

576 

577 def get_define_str(self, natoms): 

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

579 

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

581 return self.define_str 

582 

583 define_str_tpl = ( 

584 '\n__title__\na coord\n__inter__\n' 

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

586 '__occ_str____single_atom_str____norb_str____dft_str____ri_str__' 

587 '__scfiterlimit____fermi_str____damp_str__q\n' 

588 ) 

589 

590 params = self 

591 

592 if params['use redundant internals']: 

593 internals_str = 'ired\n*' 

594 else: 

595 internals_str = '*\nno' 

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

597 

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

599 if params['uhf']: 

600 occ_str = 'n\ns\n*\n' 

601 else: 

602 occ_str = 'y\n' 

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

604 occ_str = 'y\n' 

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

606 occ_str = 'n\nt\n*\n' 

607 else: 

608 unpaired = params['multiplicity'] - 1 

609 if params['use fermi smearing']: 

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

611 else: 

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

613 

614 if natoms != 1: 

615 single_atom_str = '' 

616 else: 

617 single_atom_str = '\n' 

618 

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

620 norb_str = '' 

621 else: 

622 norb_str = 'n\n' 

623 

624 if params['use dft']: 

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

626 else: 

627 dft_str = '' 

628 

629 if params['density functional']: 

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

631 

632 if params['grid size']: 

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

634 

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

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

637 else: 

638 ri_str = '' 

639 

640 if params['scf iterations']: 

641 scfmaxiter = params['scf iterations'] 

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

643 else: 

644 scfiter_str = '' 

645 if params['scf energy convergence']: 

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

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

648 

649 fermi_str = '' 

650 if params['use fermi smearing']: 

651 fermi_str = 'scf\nfermi\n' 

652 if params['fermi initial temperature']: 

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

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

655 if params['fermi final temperature']: 

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

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

658 if params['fermi annealing factor']: 

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

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

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

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

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

664 if params['fermi stopping criterion']: 

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

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

667 fermi_str += '\n\n' 

668 

669 damp_str = '' 

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

671 'minimal damping') 

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

673 if any(damp_pars): 

674 damp_str = 'scf\ndamp\n' 

675 for par in damp_pars: 

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

677 damp_str += par_str + '\n' 

678 damp_str += '\n' 

679 

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

681 

682 define_str = define_str_tpl 

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

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

685 define_str) 

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

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

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

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

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

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

692 define_str) 

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

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

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

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

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

698 

699 return define_str 

700 

701 def read_restart(self, atoms, results): 

702 """read parameters from control file""" 

703 

704 params = {} 

705 pdgs = { 

706 p: parse_data_group( 

707 read_data_group( 

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

709 ) 

710 for p in self.parameter_group 

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

712 } 

713 for p in self.parameter_key: 

714 if self.parameter_key[p]: 

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

716 if pdgs[p] is None: 

717 if self.parameter_type[p] is bool: 

718 params[p] = False 

719 else: 

720 params[p] = None 

721 else: 

722 if self.parameter_type[p] is bool: 

723 params[p] = True 

724 else: 

725 typ = self.parameter_type[p] 

726 val = typ(pdgs[p]) 

727 mapping = self.parameter_mapping 

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

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

730 val = fun(val) 

731 params[p] = val 

732 else: 

733 if pdgs[p] is None: 

734 params[p] = None 

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

736 if self.parameter_type[p] is bool: 

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

738 else: 

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

740 if self.parameter_type[p] is bool: 

741 params[p] = False 

742 else: 

743 params[p] = None 

744 else: 

745 typ = self.parameter_type[p] 

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

747 mapping = self.parameter_mapping 

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

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

750 val = fun(val) 

751 params[p] = val 

752 

753 # non-group or non-key parameters 

754 

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

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

757 assert len(basis_sets) == 1 

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

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

760 

761 # rohf, multiplicity and total charge 

762 orbs = results['molecular orbitals'] 

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

764 bool(len(read_data_group('roothaan')))) 

765 core_charge = 0 

766 if results['ecps']: 

767 for ecp in results['ecps']: 

768 for symbol in atoms.get_chemical_symbols(): 

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

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

771 if params['uhf']: 

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

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

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

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

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

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

778 electron_charge += core_charge 

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

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

781 params['multiplicity'] = 1 

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

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

784 electron_charge += core_charge 

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

786 else: 

787 raise NotImplementedError('ROHF not implemented') 

788 

789 # task-related parameters 

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

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

792 lines = log.readlines() 

793 for line in lines: 

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

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

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

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

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

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

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

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

802 if 'AN OPTIMIZATION WITH MAX' in line: 

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

804 params['geometry optimization iterations'] = cy 

805 self.update(params) 

806 self.params_old = params 

807 

808 def update_restart(self, dct): 

809 """update parameters after a restart""" 

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

811 if len(nulst) != 0: 

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

813 self.update(dct) 

814 self.update_data_groups(dct)