forked from CCQC/optavc
-
Notifications
You must be signed in to change notification settings - Fork 0
/
main.py
227 lines (182 loc) · 8.93 KB
/
main.py
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
"""
Wrapper's for running optavc with a single method call.
run_optavc() currently supports:
running finite differences of singlepoints to calculate Hessians
Running optimizations by performing finite differences of gradients.
run_optavc_mpi()
has never been used or tested all code is legacy for use on NERSC. Use WILL require work
"""
import psi4
from optavc.options import Options
from . import optimize
from . import findifcalcs
from . import xtpl
from .template import TemplateFileProcessor
from .calculations import AnalyticHessian
def run_optavc(jobtype,
options_dict,
restart_iteration=0,
xtpl_restart=None,
sow=True,
path=".",
test_input=False,
molecule=None):
"""Optavc driver to unify and simplify input. This is the only internal method of optavc a user should need to
interact with. Creates any needed internal class objects based on the provided dictionary.
Performs an Optimization or Hessian calculation
Parameters
----------
jobtype: str
upper and lowercase variants of HESS, FREQUENCY, FREQUENCIES, HESSIAN, OPT, and OPTIMIZATION are allowed
options_dict: dict
should contain BOTH optavc and psi4 options for optking and the finite difference proecedure.
should not contain ANY options relevant for the calculations that optking will submit.
restart_iteration: int
For optimizations only. orresponds to the iteration we should begin trying to sow gradients for
sow: bool
For hessians only. if False optavc will reap the hessian.
path: str
for Hessians only. Specifies where to write displacement directories.
test_input: bool
Check all calculation in STEP00, HESS, or otherwise. This is a good way to test that the energy regex
is working as expected. Raises an AssertionError if optavc considers any calculations to have failed.
Returns
-------
result : list[int] or list[list[int]]
the gradient after optimization the hessian after a hessian calculation
energy : int
final energy or energy of non displaced point
molecule : molecule.Molecule
final molecule.
"""
options_obj, input_obj, molecule = initialize_optavc(options_dict, molecule)
xtpl_inputs = None
if test_input:
try:
test_singlepoints(jobtype, molecule, options_obj, input_obj, xtpl_inputs, path)
except AssertionError as e:
print("test_singlepoints has failed. Please check regexes after running pytest on "
"test_gradobjects.py")
print(str(e))
raise
calc_obj, calc_type = create_calc_objects(jobtype, molecule, options_obj, input_obj, path)
if calc_type == 'OPT':
result, energy, molecule = calc_obj.run(restart_iteration, xtpl_restart)
return result, energy, molecule
elif calc_type == 'HESS':
if sow:
calc_obj.compute_result()
else:
calc_obj.get_result(force_resub=True)
final_hess_printout(calc_obj)
return calc_obj.result, calc_obj.energy, calc_obj.molecule
def initialize_optavc(options_dict, molecule=None):
"""Create the three basic objects required for running optavc. Options, InputFile, Molecule """
options_obj = Options(**options_dict)
if options_obj.xtpl:
options_obj.template_file_path = options_obj._xtpl_templates[0][0]
template_file_string = open(options_obj.template_file_path).read()
tfp = TemplateFileProcessor(template_file_string, options_obj)
input_obj = tfp.input_file_object
if not molecule:
molecule = tfp.molecule
return options_obj, input_obj, molecule
def create_calc_objects(jobtype, molecule, options_obj, input_obj, path='.'):
"""Creates the internal classes optavc uses to keep track of displacements and calculations. Unfifies possible jobtype"""
if jobtype.upper() in ['OPT', "OPTIMIZATION"]:
calc_obj = optimize.Optimization(molecule, input_obj, options_obj, path=path)
calc_type = 'OPT'
elif jobtype.upper() in ["HESS", "FREQUENCY", "FREQUENCIES", "HESSIAN"]:
calc_type = 'HESS'
if path == '.':
path = './HESS'
use_procedure, calc_obj = xtpl.xtpl_delta_wrapper("HESSIAN", molecule, options_obj, path)
if not use_procedure:
if options_obj.dertype == 'HESSIAN':
calc_obj = AnalyticHessian(molecule, input_obj, options_obj, path)
else:
calc_obj = findifcalcs.Hessian(molecule, input_obj, options_obj, path)
else:
raise ValueError(
"""Can only run hessians or optimizations. For gradients see findifcalcs.py to run or
add wrapper here""")
return calc_obj, calc_type
def final_hess_printout(hess_obj):
""" Here because we don't want printouts for each hessian in the XtpleDelta procedure """
psi_version = float(psi4.__version__[:3])
psi4.core.print_out("\n\n\n============================================================\n")
psi4.core.print_out("========================= OPTAVC ===========================\n")
psi4.core.print_out("============================================================\n")
psi4.core.print_out("\nThe final computed hessian is:\n\n")
psi4_mol_obj = hess_obj.molecule.cast_to_psi4_molecule_object()
wfn = psi4.core.Wavefunction.build(psi4_mol_obj, 'def2-tzvp')
wfn.set_hessian(psi4.core.Matrix.from_array(hess_obj.result))
wfn.set_energy(hess_obj.energy)
psi4.core.set_variable("CURRENT ENERGY", hess_obj.energy)
psi4.driver.vibanal_wfn(wfn)
if psi_version <= 1.4:
psi4.driver._hessian_write(wfn)
else:
psi4.driver.driver_findif.hessian_write(wfn)
def test_singlepoints(jobtype, molecule, options_obj, input_obj, xtpl_inputs=None, path=""):
""" Method to test that the singlepoints which will be created in the course
`jobtype` and their regexes will work with the output generated by from the given templates.
Can be called from run_optavc by specifying test_input=True """
if jobtype.upper() in ['OPT', "OPTIMIZATION"]:
opt_obj = optimize.Optimization(molecule, input_obj, options_obj, xtpl_inputs)
if opt_obj.options.xtpl:
for gradient, _ in opt_obj.create_xtpl_gradients(iteration=0, restart_iteration=0,
user_xtpl_restart=False):
ref_singlepoint = gradient.calculations[0]
assert ref_singlepoint.check_status(ref_singlepoint.options.energy_regex)
assert ref_singlepoint.get_result()
else:
# single gradient reap only
gradient = opt_obj.create_opt_gradient(iteration=0)
ref_singlepoint = gradient.calculations[0]
assert ref_singlepoint.check_status(ref_singlepoint.options.energy_regex)
assert ref_singlepoint.get_result()
elif jobtype.upper() in ["HESS", "FREQUENCY", "FREQUENCIES", "HESSIAN"]:
if path == '.':
path = './HESS'
use_procedure, xtpl_hess_obj = xtpl.xtpl_delta_wrapper("HESSIAN", molecule, options_obj,
path, iteration=0)
if use_procedure:
for hess_obj in xtpl_hess_obj.calc_objects:
ref_singlepoint = hess_obj.calculations[0]
assert ref_singlepoint.check_status(ref_singlepoint.options.energy_regex)
assert ref_singlepoint.get_result()
else:
hess_obj = findifcalcs.Hessian(molecule, input_obj, options_obj, path)
ref_singlepoint = hess_obj.calculations[0]
assert ref_singlepoint.check_status(ref_singlepoint.options.energy_regex)
assert ref_singlepoint.get_result()
def run_optavc_mpi():
from optavc.mpi4py_iface import mpirun
options_kwargs = {
'template_file_path': "template.dat",
'queue': "shared",
'command': "qchem input.dat output.dat",
'time_limit': "48:00:00", # has no effect
# 'memory' : "10 GB", #calculator uses memory, number of nodes, and numcores to
# distribute resources
'energy_regex': r"CCSD\(T\) total energy[=\s]+(-\d+\.\d+)",
'success_regex': r"beer",
'fail_regex': r"coffee",
'input_name': "input.dat",
'output_name': "output.dat",
# 'readhess' : False,
'mpi': True, # set to true to use mpi, false to not
# 'submitter' : None,
'maxiter': 20,
'findif': {
'points': 3
}, # set to 5 if you want slightly better frequencies
# 'job_array' : True,
'optking': {
'max_force_g_convergence': 1e-7, # tighter than this is not recommended
'rms_force_g_convergence': 1e-7,
}
}
options_obj = Options(**options_kwargs)
mpirun(options_obj)