-
Notifications
You must be signed in to change notification settings - Fork 0
/
run_opt.py
63 lines (58 loc) · 1.99 KB
/
run_opt.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
from __future__ import division
import math, os, sys
import time
import iotbx.pdb
from scitbx.array_family import flex
import mmtbx.model
from libtbx.utils import null_out
from libtbx import easy_run
base=[
"qr.refine",
"%s",
"mode=opt",
"stpmax=0.2",
"gradient_only=true", # no problem with running up to upper bound.
"clustering=false",
"use_convergence_test=true",
"rmsd_tolerance=0.001",
"number_of_micro_cycles=25"
]
cmd_cctbx = " ".join(base+["restraints=cctbx","> %s.log"])
cmd_xtb = " ".join(base+["restraints=qm","engine_name=xtb","quantum.nproc=8","> %s.log"])
cmd_terachem = " ".join(base+["restraints=qm","engine_name=terachem",
"basis=6-31g","> %s.log"])
cmd_mopac = " ".join(base+["restraints=qm","engine_name=mopac","> %s.log"])
def run(args):
assert len(args)==1
mode = args[0]
assert mode in ["cctbx", "xtb", "terachem", "mopac"]
#
if( mode=="cctbx"): cmd = cmd_cctbx
elif(mode=="xtb"): cmd = cmd_xtb
elif(mode=="terachem"): cmd = cmd_terachem
elif(mode=="mopac"): cmd = cmd_mopac
else: assert 0
#
results_prefix = "./%s_opt"%mode
perturbed_prefix = "./perturbed/"
rmsd_dirs = ["0.3/","0.6/","0.9/","1.2/","1.5/"]
base_names = ["0","1","2","3","4","5","6","7","8","9"]
easy_run.call("rm -rf %s_opt"%mode)
os.makedirs(results_prefix)
os.chdir(results_prefix)
for rmsd_dir in rmsd_dirs:
print rmsd_dir
os.makedirs(rmsd_dir)
for bn in base_names:
file_name = "../"+perturbed_prefix+rmsd_dir+bn+".pdb"
if(not os.path.exists(file_name)): assert 0
cmd_full = cmd%(file_name,bn)
print "running command:\n%s"%(cmd_full)
easy_run.call(cmd_full)
easy_run.call("cp pdb/%s_refined.pdb %s/%s.pdb"%(bn,rmsd_dir,bn))
easy_run.call("mv %s.log %s"%(bn,rmsd_dir))
easy_run.call("mv pdb %s/%s_pdb"%(rmsd_dir,bn))
if __name__ == "__main__":
t0 = time.time()
run(args=sys.argv[1:])
print "Time: %6.4f s" % (time.time() - t0)