LOOS  v2.3.2
NAMD.py
1 #!/usr/bin/env python
2 
3 import subprocess
4 import sys
5 import os.path
6 import loos
7 
8 # @cond TOOLS_INTERNAL
9 class NAMD:
10 
11  def __init__(self, psf_file, start_pdb, end_pdb, param_file, box,
12  command):
13  self.psf_file = psf_file
14  self.start_pdb = start_pdb
15  self.end_pdb = end_pdb
16  self.param_file = param_file
17 
18  if command is None:
19  self.command = "/opt/bin/namd2"
20  else:
21  self.command = command
22 
23  # in case we're going to do constraints, construct this filename
24  self.cons_k_filename = self.start_pdb[:-4] + ".cons.pdb"
25 
26  self.x_box = box.x()
27  self.y_box = box.y()
28  self.z_box = box.z()
29 
30  def update_box(self, box):
31  self.x_box = box.x()
32  self.y_box = box.y()
33  self.z_box = box.z()
34 
35 
36  def construct_header(self):
37  lines = []
38  lines.append("structure " + self.psf_file)
39  lines.append("coordinates " + self.start_pdb)
40  lines.append("outputname " + self.end_pdb)
41  lines.append("paratypecharmm on")
42  for p in self.param_file:
43  lines.append("parameters " + p)
44  lines.append("exclude scaled1-4")
45  lines.append("1-4scaling 1.0 ")
46  lines.append("binaryoutput no")
47 
48  string = "\n".join(lines)
49  string += self.template() + "\n"
50  return string
51 
52  def construct_box(self):
53  lines = []
54  lines.append("cellBasisVector1 " + str(self.x_box) + " 0 0")
55  lines.append("cellBasisVector2 0 " + str(self.y_box) + " 0")
56  lines.append("cellBasisVector3 0 0 " + str(self.z_box))
57 
58  string = "\n".join(lines)
59  string += "\n"
60  return string
61 
62  def construct_mini(self, num_iter=100):
63  line = "minimize " + str(num_iter) + "\n"
64  return line
65 
66  def construct_dyn(self, num_iter=100):
67  line = "run " + str(num_iter) + "\n"
68  return line
69 
70  def write_inputfile(self, filename, nsteps=100):
71  file = open(filename, "w")
72  file.write(self.construct_header())
73  file.write(self.construct_constraints())
74  file.write(self.construct_box())
75  file.write(self.construct_mini(nsteps))
76  file.write(self.construct_dyn(nsteps))
77  file.close()
78 
79  def write_restraintfile(self, directory, atomicgroup, spring=10.0):
80  pdb = loos.PDB.fromAtomicGroup(atomicgroup.copy())
81  for atom in pdb:
82  atom.coords().zero()
83  heavy = loos.selectAtoms(pdb, '!(name =~ "^H")' )
84  for atom in heavy:
85  atom.coords().x(spring)
86 
87  pdb_file = open(os.path.join(directory,self.cons_k_filename), "w")
88  pdb_file.write(str(pdb))
89  pdb_file.close()
90 
91 
92  def construct_constraints(self):
93  lines = [ "constraints on",
94  "selectConstraints on",
95  "selectConstrZ on",
96  "conskcol X"
97  ]
98  line = "consref " + self.start_pdb
99  lines.append(line)
100 
101  line = "conskfile " + self.cons_k_filename
102  lines.append(line)
103 
104  line = "\n".join(lines)
105  line += "\n"
106  return line
107 
108  def run_namd(self, inputfilename, outfilename):
109  """
110  Run namd on 4 processors, and report any failure to stderr
111  """
112  outfile = open(outfilename, "w")
113  try:
114  subprocess.check_call([self.command, "+p4", inputfilename],
115  stdout=outfile)
116  except subprocess.CalledProcessError:
117  sys.stderr.write("NAMD call failed, inp = %s, out = %s\n" %
118  (inputfilename, outfilename))
119  sys.exit(-1)
120 
121  def template(self):
122  """
123  Template entries for the NAMD input file. Note that these are
124  note settings you'd want to use for any kind of production run
125  (particularly the insanely short cutoff). However, we're mainly
126  using these runs for local minimization and clash-alleviation, so
127  this approach is fine and much faster than doing a realistic cutoff
128  and Ewald summation.
129  """
130  return """
131 temperature 500
132 stepsPerCycle 20
133 switching on
134 switchdist 5
135 cutoff 6
136 pairlistdist 7
137 COMmotion no
138 zeroMomentum yes
139 rigidBonds all
140 rigidTolerance 0.0000000001
141 useSettle on
142 wrapAll on
143 outputenergies 50
144 timestep 2.0
145 langevin on
146 langevinTemp 500
147 langevinDamping 2
148 langevinHydrogen off
149 
150 """
151 # @endcond
152 
153 if __name__ == '__main__':
154  import loos
155 
156  box = loos.GCoord(377, 377, 1000)
157  n = NAMD("generated.psf", "t.pdb", "end", "toppar/par_build.inp", box)
158  print n.construct_header()
159 
160  print n.construct_box()
161 
162  box[0] = 60
163  n.update_box(box)
164  print n.construct_box()
165  box[0] = 377
166  n.update_box(box)
167 
168  print n.construct_mini()
169  #print n.construct_mini(1000)
170 
171  n.write_inputfile("n.inp", 50)
172  print "launching NAMD"
173  n.run_namd("n.inp", "n.out")
174  print "finished NAMD"
Definition: NAMD.py:1
AtomicGroup selectAtoms(const AtomicGroup &source, const std::string selection)
Applies a string-based selection to an atomic group...
Definition: utils.cpp:195
static PDB fromAtomicGroup(const AtomicGroup &)
Class method for creating a PDB from an AtomicGroup.
Definition: pdb.cpp:528