1 #!/bin/env python
2 #
3 # File: OpenMMPerformSimulatedAnnealing.py
4 # Author: Manish Sud <msud@san.rr.com>
5 #
6 # Copyright (C) 2025 Manish Sud. All rights reserved.
7 #
8 # The functionality available in this script is implemented using OpenMM, an
9 # open source molecuar simulation package.
10 #
11 # This file is part of MayaChemTools.
12 #
13 # MayaChemTools is free software; you can redistribute it and/or modify it under
14 # the terms of the GNU Lesser General Public License as published by the Free
15 # Software Foundation; either version 3 of the License, or (at your option) any
16 # later version.
17 #
18 # MayaChemTools is distributed in the hope that it will be useful, but without
19 # any warranty; without even the implied warranty of merchantability of fitness
20 # for a particular purpose. See the GNU Lesser General Public License for more
21 # details.
22 #
23 # You should have received a copy of the GNU Lesser General Public License
24 # along with MayaChemTools; if not, see <http://www.gnu.org/licenses/> or
25 # write to the Free Software Foundation Inc., 59 Temple Place, Suite 330,
26 # Boston, MA, 02111-1307, USA.
27 #
28
29 from __future__ import print_function
30
31 # Add local python path to the global path and import standard library modules...
32 import os
33 import sys; sys.path.insert(0, os.path.join(os.path.dirname(sys.argv[0]), "..", "lib", "Python"))
34 import time
35 import re
36
37 # OpenMM imports...
38 try:
39 import openmm as mm
40 import openmm.app
41 except ImportError as ErrMsg:
42 sys.stderr.write("\nFailed to import OpenMM related module/package: %s\n" % ErrMsg)
43 sys.stderr.write("Check/update your OpenMM environment and try again.\n\n")
44 sys.exit(1)
45
46 # MDTraj import...
47 try:
48 import mdtraj
49 except ImportError as ErrMsg:
50 sys.stderr.write("\nFailed to import MDTraj: %s\n" % ErrMsg)
51 sys.stderr.write("Check/update your MDTraj environment and try again.\n\n")
52 sys.exit(1)
53
54 # MayaChemTools imports...
55 try:
56 from docopt import docopt
57 import MiscUtil
58 import OpenMMUtil
59 except ImportError as ErrMsg:
60 sys.stderr.write("\nFailed to import MayaChemTools module/package: %s\n" % ErrMsg)
61 sys.stderr.write("Check/update your MayaChemTools environment and try again.\n\n")
62 sys.exit(1)
63
64 ScriptName = os.path.basename(sys.argv[0])
65 Options = {}
66 OptionsInfo = {}
67
68 def main():
69 """Start execution of the script."""
70
71 MiscUtil.PrintInfo("\n%s (OpenMM v%s; MayaChemTools v%s; %s): Starting...\n" % (ScriptName, mm.Platform.getOpenMMVersion(), MiscUtil.GetMayaChemToolsVersion(), time.asctime()))
72
73 (WallClockTime, ProcessorTime) = MiscUtil.GetWallClockAndProcessorTime()
74
75 # Retrieve command line arguments and options...
76 RetrieveOptions()
77
78 # Process and validate command line arguments and options...
79 ProcessOptions()
80
81 # Perform actions required by the script...
82 PerformSimulatedAnnealing()
83
84 MiscUtil.PrintInfo("\n%s: Done...\n" % ScriptName)
85 MiscUtil.PrintInfo("Total time: %s" % MiscUtil.GetFormattedElapsedTime(WallClockTime, ProcessorTime))
86
87 def PerformSimulatedAnnealing():
88 """Perform MD simulation."""
89
90 # Prepare system for simulation...
91 System, Barostat, Topology, Positions = PrepareSystem()
92
93 # Freeze and restraint atoms...
94 FreezeRestraintAtoms(System, Topology, Positions)
95
96 # Setup integrator...
97 Integrator = SetupIntegrator()
98
99 # Setup simulation...
100 Simulation = SetupSimulation(System, Integrator, Topology, Positions)
101
102 # Write setup files...
103 WriteSimulationSetupFiles(System, Integrator)
104
105 # Perform minimization...
106 PerformMinimization(Simulation)
107
108 # Set up intial velocities...
109 SetupInitialVelocities(Simulation)
110
111 # Setup reporters for production run...
112 SetupReporters(Simulation)
113
114 # Perform annealing...
115 PerformSimulatedAnnealingWorkflow(Simulation, Integrator, Barostat)
116
117 # Save final state files...
118 WriteFinalStateFiles(Simulation)
119
120 # Reimage and realign trajectory for periodic systems...
121 ProcessTrajectory(System, Topology)
122
123 def PrepareSystem():
124 """Prepare system for simulation."""
125
126 System, Topology, Positions = OpenMMUtil.InitializeSystem(OptionsInfo["Infile"], OptionsInfo["ForcefieldParams"], OptionsInfo["SystemParams"], OptionsInfo["WaterBox"], OptionsInfo["WaterBoxParams"], OptionsInfo["SmallMolFile"], OptionsInfo["SmallMolID"])
127
128 Barostat = None
129 if OptionsInfo["NPTMode"]:
130 if not OpenMMUtil.DoesSystemUsesPeriodicBoundaryConditions(System):
131 MiscUtil.PrintInfo("")
132 MiscUtil.PrintWarning("A barostat needs to be added for NPT simulation. It appears that your system is a non-periodic system and OpenMM might fail during the initialization of the system. You might want to specify a periodic system or add water box to automatically set up a periodic system. ")
133
134 Barostat = OpenMMUtil.InitializeBarostat(OptionsInfo["IntegratorParams"])
135 MiscUtil.PrintInfo("Adding barostat for NPT simulation...")
136 try:
137 System.addForce(Barostat)
138 except Exception as ErrMsg:
139 MiscUtil.PrintInfo("")
140 MiscUtil.PrintError("Failed to add barostat:\n%s\n" % (ErrMsg))
141
142 # Write out a PDB file for the system...
143 PDBFile = OptionsInfo["PDBOutfile"]
144 MiscUtil.PrintInfo("\nWriting PDB file %s..." % PDBFile)
145 OpenMMUtil.WritePDBFile(PDBFile, Topology, Positions, OptionsInfo["OutputParams"]["PDBOutKeepIDs"])
146
147 return (System, Barostat, Topology, Positions)
148
149 def SetupIntegrator():
150 """Setup integrator. """
151
152 Integrator = OpenMMUtil.InitializeIntegrator(OptionsInfo["IntegratorParams"], OptionsInfo["SystemParams"]["ConstraintErrorTolerance"])
153
154 return Integrator
155
156 def SetupSimulation(System, Integrator, Topology, Positions):
157 """Setup simulation. """
158
159 Simulation = OpenMMUtil.InitializeSimulation(System, Integrator, Topology, Positions, OptionsInfo["PlatformParams"])
160
161 return Simulation
162
163 def SetupInitialVelocities(Simulation):
164 """Setup initial velocities."""
165
166 # Set velocities to random values choosen from a Boltzman distribution at a given
167 # temperature...
168 AnnealingParams = OptionsInfo["AnnealingParams"]
169 AnnealingParamsInfo= OpenMMUtil.SetupAnnealingParameters(OptionsInfo["AnnealingParams"])
170
171 MiscUtil.PrintInfo("\nSetting initial velocities to temperature (%s K)..." % AnnealingParams["InitialStart"])
172 Simulation.context.setVelocitiesToTemperature(AnnealingParamsInfo["InitialStart"])
173
174 def PerformMinimization(Simulation):
175 """Perform minimization."""
176
177 SimulationParams = OpenMMUtil.SetupSimulationParameters(OptionsInfo["SimulationParams"])
178
179 if not SimulationParams["Minimization"]:
180 MiscUtil.PrintInfo("\nSkipping energy minimization...")
181 return
182
183 OutputParams = OptionsInfo["OutputParams"]
184
185 # Setup a local minimization reporter...
186 MinimizeReporter = None
187 if OutputParams["MinimizationDataStdout"] or OutputParams["MinimizationDataLog"]:
188 MinimizeReporter = LocalMinimizationReporter()
189
190 if MinimizeReporter is not None:
191 MiscUtil.PrintInfo("\nAdding minimization reporters...")
192 if OutputParams["MinimizationDataLog"]:
193 MiscUtil.PrintInfo("Adding data log minimization reporter (Steps: %s; File: %s)..." % (OutputParams["MinimizationDataSteps"], OutputParams["MinimizationDataLogFile"]))
194 if OutputParams["MinimizationDataStdout"]:
195 MiscUtil.PrintInfo("Adding data stdout minimization reporter (Steps: %s)..." % (OutputParams["MinimizationDataSteps"]))
196 else:
197 MiscUtil.PrintInfo("\nSkipping addition of minimization reporters...")
198
199 MaxSteps = SimulationParams["MinimizationMaxSteps"]
200
201 MaxStepsMsg = "MaxSteps: %s" % ("UntilConverged" if MaxSteps == 0 else MaxSteps)
202 ToleranceMsg = "Tolerance: %.2f kcal/mol/A (%.2f kjoules/mol/nm)" % (SimulationParams["MinimizationToleranceInKcal"], SimulationParams["MinimizationToleranceInJoules"])
203
204 MiscUtil.PrintInfo("\nPerforming energy minimization (%s; %s)..." % (MaxStepsMsg, ToleranceMsg))
205
206 if OutputParams["MinimizationDataStdout"]:
207 HeaderLine = SetupMinimizationDataOutHeaderLine()
208 print("\n%s" % HeaderLine)
209
210 Simulation.minimizeEnergy(tolerance = SimulationParams["MinimizationTolerance"], maxIterations = MaxSteps, reporter = MinimizeReporter)
211
212 if OutputParams["MinimizationDataLog"]:
213 WriteMinimizationDataLogFile(MinimizeReporter.DataOutValues)
214
215 if OutputParams["PDBOutMinimized"]:
216 MiscUtil.PrintInfo("\nWriting PDB file %s..." % OptionsInfo["MinimizedPDBOutfile"])
217 OpenMMUtil.WriteSimulationStatePDBFile(Simulation, OptionsInfo["MinimizedPDBOutfile"], OutputParams["PDBOutKeepIDs"])
218
219 def PerformSimulatedAnnealingWorkflow(Simulation, Integrator, Barostat):
220 """Perform simulated annealing workflow."""
221
222 MiscUtil.PrintInfo("\nPerforming simulated annealing...")
223
224 AnnealingParams = OptionsInfo["AnnealingParams"]
225 TotalSimulationSteps = 0
226
227 # Perform intial heating along with equilibration...
228 InitialStart = AnnealingParams["InitialStart"]
229 InitialEnd = AnnealingParams["InitialEnd"]
230 InitialChange = AnnealingParams["InitialChange"]
231 InitialSteps = AnnealingParams["InitialSteps"]
232
233 MiscUtil.PrintInfo("\nPerforming initial heating (Start: %.1f K; End: %.1f K; Change: %.1f K)..." % (InitialStart, InitialEnd, InitialChange))
234 TotalInitialSimulationSteps = OpenMMUtil.PerformAnnealing(Simulation, Integrator, Barostat, InitialStart, InitialEnd, InitialChange, InitialSteps)
235 MiscUtil.PrintInfo("Finished initial heating (TotalSteps: %s; TotalTime: %s)..." % (TotalInitialSimulationSteps, GetTotalSimulationTime(TotalInitialSimulationSteps)))
236 TotalSimulationSteps += TotalInitialSimulationSteps
237
238 # Perform equilibration after intial heating...
239 InitialEquilibrationSteps = AnnealingParams["InitialEquilibrationSteps"]
240 MiscUtil.PrintInfo("\nPerforming equilibration after initial heating (Steps: %s; Time: %s)..." % (InitialEquilibrationSteps, GetTotalSimulationTime(InitialEquilibrationSteps)))
241 Simulation.step(InitialEquilibrationSteps)
242 TotalSimulationSteps += InitialEquilibrationSteps
243
244 # Peform heating and coolling annealing cycles along with equilibration...
245 Cycles = AnnealingParams["Cycles"]
246 CycleStart = AnnealingParams["CycleStart"]
247 CycleEnd = AnnealingParams["CycleEnd"]
248 CycleChange = AnnealingParams["CycleChange"]
249 CycleSteps = AnnealingParams["CycleSteps"]
250 CycleEquilibrationSteps = AnnealingParams["CycleEquilibrationSteps"]
251
252 MiscUtil.PrintInfo("\nPerforming heating and cooling cycles (NumCycles: %s)..." % (Cycles))
253 for Cycle in range(Cycles):
254 MiscUtil.PrintInfo("\nPerforming heating and cooling cycle %s..." % (Cycle + 1))
255
256 MiscUtil.PrintInfo("\nPerforming heating (Start: %.1f K; End: %.1f K; Change: %.1f K)..." % (CycleStart, CycleEnd, CycleChange))
257 TotalCycleSimulationSteps = OpenMMUtil.PerformAnnealing(Simulation, Integrator, Barostat, CycleStart, CycleEnd, CycleChange, CycleSteps)
258 MiscUtil.PrintInfo("Finished heating cycle (TotalSteps: %s; TotalTime: %s)..." % (TotalCycleSimulationSteps, GetTotalSimulationTime(TotalCycleSimulationSteps)))
259 TotalSimulationSteps += TotalCycleSimulationSteps
260
261 MiscUtil.PrintInfo("\nPerforming equilibration (Steps: %s; Time: %s)..." % (CycleEquilibrationSteps, GetTotalSimulationTime(CycleEquilibrationSteps)))
262 Simulation.step(CycleEquilibrationSteps)
263 TotalSimulationSteps += CycleEquilibrationSteps
264
265 MiscUtil.PrintInfo("\nPerforming cooling (Start: %.1f K; End: %.1f K; Change: %.1f K)..." % (CycleEnd, CycleStart, CycleChange))
266 TotalCycleSimulationSteps = OpenMMUtil.PerformAnnealing(Simulation, Integrator, Barostat, CycleEnd, CycleStart, CycleChange, CycleSteps)
267 MiscUtil.PrintInfo("Finished cooling cycle (TotalSteps: %s; TotalTime: %s)..." % (TotalCycleSimulationSteps, GetTotalSimulationTime(TotalCycleSimulationSteps)))
268 TotalSimulationSteps += TotalCycleSimulationSteps
269
270 MiscUtil.PrintInfo("\nPerforming equilibration (Steps: %s; Time: %s)..." % (CycleEquilibrationSteps, GetTotalSimulationTime(CycleEquilibrationSteps)))
271 Simulation.step(CycleEquilibrationSteps)
272 TotalSimulationSteps += CycleEquilibrationSteps
273
274 MiscUtil.PrintInfo("\nFinished heating and cooling cycle %s..." % (Cycle + 1))
275
276 # Perform final equilibration...
277 FinalEquilibrationSteps = AnnealingParams["FinalEquilibrationSteps"]
278 MiscUtil.PrintInfo("\nPerforming final equilibration after heating and cooling cycles (Steps: %s; Time: %s)..." % (FinalEquilibrationSteps, GetTotalSimulationTime(FinalEquilibrationSteps)))
279 Simulation.step(FinalEquilibrationSteps)
280
281 TotalSimulationSteps += FinalEquilibrationSteps
282 MiscUtil.PrintInfo("\nFinishing simulated annealing (TotalSteps: %s; TotalTime: %s)..." % (TotalSimulationSteps, GetTotalSimulationTime(TotalSimulationSteps)))
283
284 def GetTotalSimulationTime(SimulationSteps):
285 """Get total simulation time. """
286
287 IntegratorParamsInfo = OpenMMUtil.SetupIntegratorParameters(OptionsInfo["IntegratorParams"])
288 StepSize = IntegratorParamsInfo["StepSize"]
289
290 TotalTime = OpenMMUtil.GetFormattedTotalSimulationTime(StepSize, SimulationSteps)
291
292 return TotalTime
293
294 def SetupReporters(Simulation):
295 """Setup reporters. """
296
297 DataAppend = False
298 (TrajReporter, DataLogReporter, DataStdoutReporter, CheckpointReporter) = OpenMMUtil.InitializeReporters(OptionsInfo["OutputParams"], OptionsInfo["SimulationParams"]["Steps"], DataAppend)
299
300 if TrajReporter is None and DataLogReporter is None and DataStdoutReporter is None and CheckpointReporter is None:
301 MiscUtil.PrintInfo("\nSkip adding reporters...")
302 return
303
304 MiscUtil.PrintInfo("\nAdding reporters...")
305
306 OutputParams = OptionsInfo["OutputParams"]
307 AppendMsg = ""
308 if TrajReporter is not None:
309 MiscUtil.PrintInfo("Adding trajectory reporter (Steps: %s; File: %s%s)..." % (OutputParams["TrajSteps"], OutputParams["TrajFile"], AppendMsg))
310 Simulation.reporters.append(TrajReporter)
311
312 if CheckpointReporter is not None:
313 MiscUtil.PrintInfo("Adding checkpoint reporter (Steps: %s; File: %s)..." % (OutputParams["CheckpointSteps"], OutputParams["CheckpointFile"]))
314 Simulation.reporters.append(CheckpointReporter)
315
316 if DataLogReporter is not None:
317 MiscUtil.PrintInfo("Adding data log reporter (Steps: %s; File: %s%s)..." % (OutputParams["DataLogSteps"], OutputParams["DataLogFile"], AppendMsg))
318 Simulation.reporters.append(DataLogReporter)
319
320 if DataStdoutReporter is not None:
321 MiscUtil.PrintInfo("Adding data stdout reporter (Steps: %s)..." % (OutputParams["DataStdoutSteps"]))
322 Simulation.reporters.append(DataStdoutReporter)
323
324 class LocalMinimizationReporter(mm.MinimizationReporter):
325 """Setup a local minimization reporter. """
326
327 (DataSteps, DataOutTypeList, DataOutDelimiter, StdoutStatus) = [None] * 4
328
329 DataOutValues = []
330 First = True
331
332 def report(self, Iteration, PositonsList, GradientsList, DataStatisticsMap):
333 """Report and track minimization."""
334
335 if self.First:
336 # Initialize...
337 self.DataSteps = OptionsInfo["OutputParams"]["MinimizationDataSteps"]
338 self.DataOutTypeList = OptionsInfo["OutputParams"]["MinimizationDataOutTypeOpenMMNameList"]
339 self.DataOutDelimiter = OptionsInfo["OutputParams"]["DataOutDelimiter"]
340 self.StdoutStatus = True if OptionsInfo["OutputParams"]["MinimizationDataStdout"] else False
341
342 self.First = False
343
344 if Iteration % self.DataSteps == 0:
345 # Setup data values...
346 DataValues = []
347 DataValues.append("%s" % Iteration)
348 for DataType in self.DataOutTypeList:
349 DataValue = "%.4f" % DataStatisticsMap[DataType]
350 DataValues.append(DataValue)
351
352 # Track data...
353 self.DataOutValues.append(DataValues)
354
355 # Print data values...
356 if self.StdoutStatus:
357 print("%s" % self.DataOutDelimiter.join(DataValues))
358
359 # This method must return a bool. You may return true for early termination.
360 return False
361
362 def WriteMinimizationDataLogFile(DataOutValues):
363 """Write minimization data log file. """
364
365 OutputParams = OptionsInfo["OutputParams"]
366
367 Outfile = OutputParams["MinimizationDataLogFile"]
368 OutDelimiter = OutputParams["DataOutDelimiter"]
369
370 MiscUtil.PrintInfo("\nWriting minimization log file %s..." % Outfile)
371 OutFH = open(Outfile, "w")
372
373 HeaderLine = SetupMinimizationDataOutHeaderLine()
374 OutFH.write("%s\n" % HeaderLine)
375
376 for LineWords in DataOutValues:
377 Line = OutDelimiter.join(LineWords)
378 OutFH.write("%s\n" % Line)
379
380 OutFH.close()
381
382 def SetupMinimizationDataOutHeaderLine():
383 """Setup minimization data output header line. """
384
385 LineWords = ["Iteration"]
386 for Label in OptionsInfo["OutputParams"]["MinimizationDataOutTypeList"]:
387 if re.match("^(SystemEnergy|RestraintEnergy)$", Label, re.I):
388 LineWords.append("%s(kjoules/mol)" % Label)
389 elif re.match("^RestraintStrength$", Label, re.I):
390 LineWords.append("%s(kjoules/mol/nm^2)" % Label)
391 else:
392 LineWords.append(Label)
393
394 Line = OptionsInfo["OutputParams"]["DataOutDelimiter"].join(LineWords)
395
396 return Line
397
398 def FreezeRestraintAtoms(System, Topology, Positions):
399 """Handle freezing and restraining of atoms."""
400
401 # Get atoms for freezing...
402 FreezeAtomList = None
403 if OptionsInfo["FreezeAtoms"]:
404 FreezeAtomList = OpenMMUtil.GetAtoms(Topology, OptionsInfo["FreezeAtomsParams"]["CAlphaProteinStatus"], OptionsInfo["FreezeAtomsParams"]["ResidueNames"], OptionsInfo["FreezeAtomsParams"]["Negate"])
405 if FreezeAtomList is None:
406 MiscUtil.PrintError("The freeze atoms parameters specified, \"selection, %s, selectionSpec, %s, negate, %s\", - using \"--freezeAtomsParams\" option didn't match any atoms in the system. You must specify a valid set of parameters for freezing atoms or disable freezing using \"No\" value for \"--freezeAtoms\" option.." % (OptionsInfo["FreezeAtomsParams"]["Selection"], OptionsInfo["FreezeAtomsParams"]["SelectionSpec"], OptionsInfo["FreezeAtomsParams"]["Negate"]))
407
408 # Get atoms for restraining...
409 RestraintAtomList = None
410 if OptionsInfo["RestraintAtoms"]:
411 RestraintAtomList = OpenMMUtil.GetAtoms(Topology, OptionsInfo["RestraintAtomsParams"]["CAlphaProteinStatus"], OptionsInfo["RestraintAtomsParams"]["ResidueNames"], OptionsInfo["RestraintAtomsParams"]["Negate"])
412 if RestraintAtomList is None:
413 MiscUtil.PrintError("The restraint atoms parameters specified, \"selection, %s, selectionSpec, %s, negate, %s\", - using \"--restraintAtomsParams\" option didn't match any atoms in the system. You must specify a valid set of parameters for restraining atoms or disable restraining using \"No\" value for \"--restraintAtoms\" option." % (OptionsInfo["RestraintAtomsParams"]["Selection"], OptionsInfo["RestraintAtomsParams"]["SelectionSpec"], OptionsInfo["RestraintAtomsParams"]["Negate"]))
414
415 # Check for atoms to freeze or restraint...
416 if FreezeAtomList is None and RestraintAtomList is None:
417 return
418
419 # Check for overlap between freeze and restraint atoms...
420 if OpenMMUtil.DoAtomListsOverlap(FreezeAtomList, RestraintAtomList):
421 MiscUtil.PrintError("The atoms specified using \"--freezeAtomsParams\" and \"--restraintAtomsParams\" options appear to overlap. You must specify unique sets of atoms to freeze and restraint.")
422
423 # Check overlap of freeze atoms with system constraints...
424 if OpenMMUtil.DoesAtomListOverlapWithSystemConstraints(System, FreezeAtomList):
425 MiscUtil.PrintError("The atoms specified using \"--freezeAtomsParams\" appear to overlap with atoms being constrained corresponding to the value specified, %s, for paramater name \"constraints\" using \"--systemParams\" option.\n\nYou must specify a unique set of atoms to freeze or turn off system constaints by specifying value, None, for \"constraints\" parameter using option \"--systemsParams\".\n\nIn addtion, you may want specify, no, value for \"rigidWater\" option.\n\nThe atoms are frozen by setting their particle mass to zero. OpenMM doesn't allow to both constraint atoms and set their mass to zero." % (OptionsInfo["SystemParams"]["Constraints"]))
426
427 # Check overlap of restraint atoms with system constraints...
428 if OpenMMUtil.DoesAtomListOverlapWithSystemConstraints(System, RestraintAtomList):
429 MiscUtil.PrintInfo("")
430 MiscUtil.PrintWarning("The atoms specified using \"--restraintAtomsParams\" appear to overlap with atoms being constrained corresponding to the value specified, %s, for paramater name \"constraints\" using \"--systemParams\" option. You may want to specify a unique set of atoms to restraints or turn off system constaints by specifying value, None, for \"constraints\" parameter using option \"--systemsParams\"." % (OptionsInfo["SystemParams"]["Constraints"]))
431
432 # Check and adjust step size...
433 if FreezeAtomList is not None or RestraintAtomList is not None:
434 CheckAndAdjustStepSizeDuringFreezeRestraintAtoms()
435
436 # Freeze atoms...
437 if FreezeAtomList is None:
438 MiscUtil.PrintInfo("\nSkipping freezing of atoms...")
439 else:
440 MiscUtil.PrintInfo("\nFreezing atoms (Selection: %s)..." % OptionsInfo["FreezeAtomsParams"]["Selection"])
441 OpenMMUtil.FreezeAtoms(System, FreezeAtomList)
442
443 # Restraint atoms...
444 if RestraintAtomList is None:
445 MiscUtil.PrintInfo("\nSkipping restraining of atoms...")
446 else:
447 MiscUtil.PrintInfo("\nRestraining atoms (Selection: %s)..." % OptionsInfo["RestraintAtomsParams"]["Selection"])
448
449 SprintConstantInKcal = OptionsInfo["RestraintSpringConstant"]
450 OpenMMUtil.RestraintAtoms(System, Positions, RestraintAtomList, SprintConstantInKcal)
451
452 def CheckAndAdjustStepSizeDuringFreezeRestraintAtoms():
453 """Check and set step size during freezing or restraining of atoms """
454
455 if re.match("^auto$", OptionsInfo["IntegratorParams"]["StepSizeSpecified"], re.I):
456 # Automatically set stepSize to 2.0 fs..
457 OptionsInfo["IntegratorParams"]["StepSize"] = 2.0
458 MiscUtil.PrintInfo("")
459 MiscUtil.PrintWarning("The time step has been automatically set to %s fs during freezing or restraining of atoms. You may specify an explicit value for parameter name, stepSize, using \"--integratorParams\" option." % (OptionsInfo["IntegratorParams"]["StepSize"]))
460 elif OptionsInfo["IntegratorParams"]["StepSize"] > 2:
461 MiscUtil.PrintInfo("")
462 MiscUtil.PrintWarning("A word to the wise: The parameter value specified, %s, for parameter name, stepSize, using \"--integratorParams\" option may be too large. You may want to consider using a smaller time step. Othwerwise, your simulation may blow up." % (OptionsInfo["IntegratorParams"]["StepSize"] ))
463 MiscUtil.PrintInfo("")
464
465 def WriteSimulationSetupFiles(System, Integrator):
466 """Write simulation setup files for system and integrator."""
467
468 OutputParams = OptionsInfo["OutputParams"]
469
470 if OutputParams["XmlSystemOut"] or OutputParams["XmlIntegratorOut"]:
471 MiscUtil.PrintInfo("")
472
473 if OutputParams["XmlSystemOut"]:
474 Outfile = OutputParams["XmlSystemFile"]
475 MiscUtil.PrintInfo("Writing system setup XML file %s..." % Outfile)
476 with open(Outfile, mode = "w") as OutFH:
477 OutFH.write(mm.XmlSerializer.serialize(System))
478
479 if OutputParams["XmlIntegratorOut"]:
480 Outfile = OutputParams["XmlIntegratorFile"]
481 MiscUtil.PrintInfo("Writing integrator setup XML file %s..." % Outfile)
482 with open(Outfile, mode = "w") as OutFH:
483 OutFH.write(mm.XmlSerializer.serialize(Integrator))
484
485 def WriteFinalStateFiles(Simulation):
486 """Write final state files. """
487
488 OutputParams = OptionsInfo["OutputParams"]
489
490 if OutputParams["SaveFinalStateCheckpoint"] or OutputParams["SaveFinalStateXML"] or OutputParams["PDBOutFinal"]:
491 MiscUtil.PrintInfo("")
492
493 if OutputParams["SaveFinalStateCheckpoint"]:
494 Outfile = OutputParams["SaveFinalStateCheckpointFile"]
495 MiscUtil.PrintInfo("Writing final state checkpoint file %s..." % Outfile)
496 Simulation.saveCheckpoint(Outfile)
497
498 if OutputParams["SaveFinalStateXML"]:
499 Outfile = OutputParams["SaveFinalStateXMLFile"]
500 MiscUtil.PrintInfo("Writing final state XML file %s..." % Outfile)
501 Simulation.saveState(Outfile)
502
503 if OutputParams["PDBOutFinal"]:
504 MiscUtil.PrintInfo("\nWriting PDB file %s..." % OptionsInfo["FinalPDBOutfile"])
505 OpenMMUtil.WriteSimulationStatePDBFile(Simulation, OptionsInfo["FinalPDBOutfile"], OutputParams["PDBOutKeepIDs"])
506
507 def ProcessTrajectory(System, Topology):
508 """Reimage and realign trajectory for periodic systems. """
509
510 ReimageFrames = True if OpenMMUtil.DoesSystemUsesPeriodicBoundaryConditions(System) else False
511 if not ReimageFrames:
512 MiscUtil.PrintInfo("\nSkipping reimaging and realigning of trajectory for a system not using periodic boundary conditions...")
513 return
514
515 MiscUtil.PrintInfo("\nReimaging and realigning trajectory for a system using periodic boundary conditions...")
516
517 # Reimage and realign trajectory file...
518 RealignFrames = True
519 TrajTopologyFile = OptionsInfo["PDBOutfile"]
520 TrajFile = OptionsInfo["OutputParams"]["TrajFile"]
521 Traj, ReimagedStatus, RealignedStatus = OpenMMUtil.ReimageRealignTrajectory(Topology, TrajFile, ReimageFrames, RealignFrames)
522
523 if (Traj is None) or (not ReimagedStatus and not RealignedStatus):
524 MiscUtil.PrintInfo("Skipping writing first frame to PDB file %s..." % PDBOutfile)
525 MiscUtil.PrintInfo("Skippig writing trajectory file %s..." % TrajOutfile)
526 return
527
528 PDBOutFormat = OptionsInfo["OutputParams"]["PDBOutFormat"]
529 PDBOutfile = OptionsInfo["ReimagedPDBOutfile"]
530 TrajOutfile = OptionsInfo["ReimagedTrajOutfile"]
531
532 # Write out first frame...
533 MiscUtil.PrintInfo("Writing first frame to PDB file %s..." % PDBOutfile)
534 if re.match("^CIF$", PDBOutFormat, re.I):
535 # MDTraj doesn't appear to support CIF format. Write it out as a temporary
536 # PDB file and convert it using OpenMM...
537 FileDir, FileRoot, FileExt = MiscUtil.ParseFileName(PDBOutfile)
538 TmpPDBOutfile = "%s_PID%s.pdb" % (FileRoot, os.getpid())
539 Traj[0].save(TmpPDBOutfile)
540
541 PDBHandle = OpenMMUtil.ReadPDBFile(TmpPDBOutfile)
542 OpenMMUtil.WritePDBFile(PDBOutfile, PDBHandle.topology, PDBHandle.positions, OptionsInfo["OutputParams"]["PDBOutKeepIDs"])
543
544 os.remove(TmpPDBOutfile)
545 else:
546 Traj[0].save(PDBOutfile)
547
548 # Write out reimaged and realinged trajectory...
549 MiscUtil.PrintInfo("Writing trajectory file %s..." % TrajOutfile)
550 Traj.save(TrajOutfile)
551
552 def ProcessOutfilePrefixParameter():
553 """Process outfile prefix paramater."""
554
555 OutfilePrefix = Options["--outfilePrefix"]
556
557 if not re.match("^auto$", OutfilePrefix, re.I):
558 OptionsInfo["OutfilePrefix"] = OutfilePrefix
559 return
560
561 if OptionsInfo["SmallMolFileMode"]:
562 OutfilePrefix = "%s_%s_Complex" % (OptionsInfo["InfileRoot"], OptionsInfo["SmallMolFileRoot"])
563 else:
564 OutfilePrefix = "%s" % (OptionsInfo["InfileRoot"])
565
566 if re.match("^yes$", Options["--waterBox"], re.I):
567 OutfilePrefix = "%s_Solvated" % (OutfilePrefix)
568
569 OutfilePrefix = "%s_%s" % (OutfilePrefix, OptionsInfo["Mode"])
570
571 OptionsInfo["OutfilePrefix"] = OutfilePrefix
572
573 def ProcessOutfileNames():
574 """Process outfile names."""
575
576 OutputParams = OptionsInfo["OutputParams"]
577 OutfileParamNames = ["CheckpointFile", "DataLogFile", "MinimizationDataLogFile", "SaveFinalStateCheckpointFile", "SaveFinalStateXMLFile", "TrajFile", "XmlSystemFile", "XmlIntegratorFile"]
578 for OutfileParamName in OutfileParamNames:
579 OutfileParamValue = OutputParams[OutfileParamName]
580 if not Options["--overwrite"]:
581 if os.path.exists(OutfileParamValue):
582 MiscUtil.PrintError("The file specified, %s, for parameter name, %s, using option \"--outfileParams\" already exist. Use option \"--ov\" or \"--overwrite\" and try again. " % (OutfileParamValue, OutfileParamName))
583
584 PDBOutfile = "%s.%s" % (OptionsInfo["OutfilePrefix"], OutputParams["PDBOutfileExt"])
585 ReimagedPDBOutfile = "%s_Reimaged.%s" % (OptionsInfo["OutfilePrefix"], OutputParams["PDBOutfileExt"])
586 ReimagedTrajOutfile = "%s_Reimaged.%s" % (OptionsInfo["OutfilePrefix"], OutputParams["TrajFileExt"])
587
588 MinimizedPDBOutfile = "%s_Minimized.%s" % (OptionsInfo["OutfilePrefix"], OutputParams["PDBOutfileExt"])
589 FinalPDBOutfile = "%s_Final.%s" % (OptionsInfo["OutfilePrefix"], OutputParams["PDBOutfileExt"])
590
591 for Outfile in [PDBOutfile, ReimagedPDBOutfile, ReimagedTrajOutfile, MinimizedPDBOutfile, FinalPDBOutfile]:
592 if not Options["--overwrite"]:
593 if os.path.exists(Outfile):
594 MiscUtil.PrintError("The file name, %s, generated using option \"--outfilePrefix\" already exist. Use option \"--ov\" or \"--overwrite\" and try again. " % (Outfile))
595 OptionsInfo["PDBOutfile"] = PDBOutfile
596 OptionsInfo["ReimagedPDBOutfile"] = ReimagedPDBOutfile
597 OptionsInfo["ReimagedTrajOutfile"] = ReimagedTrajOutfile
598
599 OptionsInfo["MinimizedPDBOutfile"] = MinimizedPDBOutfile
600 OptionsInfo["FinalPDBOutfile"] = FinalPDBOutfile
601
602 def ProcessWaterBoxParameters():
603 """Process water box parameters. """
604
605 OptionsInfo["WaterBox"] = True if re.match("^yes$", Options["--waterBox"], re.I) else False
606 OptionsInfo["WaterBoxParams"] = OpenMMUtil.ProcessOptionOpenMMWaterBoxParameters("--waterBoxParams", Options["--waterBoxParams"])
607
608 if OptionsInfo["WaterBox"]:
609 if OptionsInfo["ForcefieldParams"]["ImplicitWater"]:
610 MiscUtil.PrintInfo("")
611 MiscUtil.PrintWarning("The value, %s, specified using option \"--waterBox\" may not be valid for the combination of biopolymer and water forcefields, %s and %s, specified using \"--forcefieldParams\". You may consider using a valid combination of biopolymer and water forcefields for explicit water during the addition of a water box." % (Options["--waterBox"], OptionsInfo["ForcefieldParams"]["Biopolymer"], OptionsInfo["ForcefieldParams"]["Water"]))
612
613 def ProcessOptions():
614 """Process and validate command line arguments and options."""
615
616 MiscUtil.PrintInfo("Processing options...")
617
618 ValidateOptions()
619
620 OptionsInfo["Infile"] = Options["--infile"]
621 FileDir, FileName, FileExt = MiscUtil.ParseFileName(OptionsInfo["Infile"])
622 OptionsInfo["InfileRoot"] = FileName
623
624 SmallMolFile = Options["--smallMolFile"]
625 SmallMolID = Options["--smallMolID"]
626 SmallMolFileMode = False
627 SmallMolFileRoot = None
628 if SmallMolFile is not None:
629 FileDir, FileName, FileExt = MiscUtil.ParseFileName(SmallMolFile)
630 SmallMolFileRoot = FileName
631 SmallMolFileMode = True
632
633 OptionsInfo["SmallMolFile"] = SmallMolFile
634 OptionsInfo["SmallMolFileRoot"] = SmallMolFileRoot
635 OptionsInfo["SmallMolFileMode"] = SmallMolFileMode
636 OptionsInfo["SmallMolID"] = SmallMolID.upper()
637
638 OptionsInfo["Mode"] = Options["--mode"].upper()
639 OptionsInfo["NPTMode"] = True if re.match("^NPT$", OptionsInfo["Mode"]) else False
640 OptionsInfo["NVTMode"] = True if re.match("^NVT$", OptionsInfo["Mode"]) else False
641
642 ProcessOutfilePrefixParameter()
643
644 if OptionsInfo["NVTMode"]:
645 ParamsDefaultInfoOverride = {"DataOutType": "Step Speed PotentialEnergy Temperature Time Volume"}
646 else:
647 ParamsDefaultInfoOverride = {"DataOutType": "Step Speed PotentialEnergy Temperature Time Density"}
648 OptionsInfo["OutputParams"] = OpenMMUtil.ProcessOptionOpenMMOutputParameters("--outputParams", Options["--outputParams"], OptionsInfo["OutfilePrefix"], ParamsDefaultInfoOverride)
649 ProcessOutfileNames()
650
651 OptionsInfo["AnnealingParams"] = OpenMMUtil.ProcessOptionOpenMMAnnealingParameters("--annealingParams", Options["--annealingParams"])
652
653 OptionsInfo["ForcefieldParams"] = OpenMMUtil.ProcessOptionOpenMMForcefieldParameters("--forcefieldParams", Options["--forcefieldParams"])
654
655 OptionsInfo["FreezeAtoms"] = True if re.match("^yes$", Options["--freezeAtoms"], re.I) else False
656 if OptionsInfo["FreezeAtoms"]:
657 OptionsInfo["FreezeAtomsParams"] = OpenMMUtil.ProcessOptionOpenMMAtomsSelectionParameters("--freezeAtomsParams", Options["--freezeAtomsParams"])
658 else:
659 OptionsInfo["FreezeAtomsParams"] = None
660
661 ParamsDefaultInfoOverride = {"Name": Options["--platform"], "Threads": 1}
662 OptionsInfo["PlatformParams"] = OpenMMUtil.ProcessOptionOpenMMPlatformParameters("--platformParams", Options["--platformParams"], ParamsDefaultInfoOverride)
663
664 OptionsInfo["RestraintAtoms"] = True if re.match("^yes$", Options["--restraintAtoms"], re.I) else False
665 if OptionsInfo["RestraintAtoms"]:
666 OptionsInfo["RestraintAtomsParams"] = OpenMMUtil.ProcessOptionOpenMMAtomsSelectionParameters("--restraintAtomsParams", Options["--restraintAtomsParams"])
667 else:
668 OptionsInfo["RestraintAtomsParams"] = None
669 OptionsInfo["RestraintSpringConstant"] = float(Options["--restraintSpringConstant"])
670
671 OptionsInfo["SystemParams"] = OpenMMUtil.ProcessOptionOpenMMSystemParameters("--systemParams", Options["--systemParams"])
672
673 OptionsInfo["IntegratorParams"] = OpenMMUtil.ProcessOptionOpenMMIntegratorParameters("--integratorParams", Options["--integratorParams"], HydrogenMassRepartioningStatus = OptionsInfo["SystemParams"]["HydrogenMassRepartioning"])
674 OptionsInfo["IntegratorParams"]["Temperature"] = OptionsInfo["AnnealingParams"]["InitialStart"]
675
676 OptionsInfo["SimulationParams"] = OpenMMUtil.ProcessOptionOpenMMSimulationParameters("--simulationParams", Options["--simulationParams"])
677
678 ProcessWaterBoxParameters()
679
680 OptionsInfo["Overwrite"] = Options["--overwrite"]
681
682 def RetrieveOptions():
683 """Retrieve command line arguments and options."""
684
685 # Get options...
686 global Options
687 Options = docopt(_docoptUsage_)
688
689 # Set current working directory to the specified directory...
690 WorkingDir = Options["--workingdir"]
691 if WorkingDir:
692 os.chdir(WorkingDir)
693
694 # Handle examples option...
695 if "--examples" in Options and Options["--examples"]:
696 MiscUtil.PrintInfo(MiscUtil.GetExamplesTextFromDocOptText(_docoptUsage_))
697 sys.exit(0)
698
699 def ValidateOptions():
700 """Validate option values."""
701
702 MiscUtil.ValidateOptionFilePath("-i, --infile", Options["--infile"])
703 MiscUtil.ValidateOptionFileExt("-i, --infile", Options["--infile"], "pdb cif")
704
705 FileDir, FileName, FileExt = MiscUtil.ParseFileName(Options["--infile"])
706 OutfilePrefix = Options["--outfilePrefix"]
707 if not re.match("^auto$", OutfilePrefix, re.I):
708 if re.match("^(%s)$" % OutfilePrefix, FileName, re.I):
709 MiscUtil.PrintError("The value specified, %s, for option \"--outfilePrefix\" is not valid. You must specify a value different from, %s, the root of infile name." % (OutfilePrefix, FileName))
710
711 if Options["--smallMolFile"] is not None:
712 MiscUtil.ValidateOptionFilePath("-l, --smallMolFile", Options["--smallMolFile"])
713 MiscUtil.ValidateOptionFileExt("-l, --smallMolFile", Options["--smallMolFile"], "sd sdf")
714
715 SmallMolID = Options["--smallMolID"]
716 if len(SmallMolID) != 3:
717 MiscUtil.PrintError("The value specified, %s, for option \"--smallMolID\" is not valid. You must specify a three letter small molecule ID." % (SmallMolID))
718
719 MiscUtil.ValidateOptionTextValue("--freezeAtoms", Options["--freezeAtoms"], "yes no")
720 if re.match("^yes$", Options["--freezeAtoms"], re.I):
721 if Options["--freezeAtomsParams"] is None:
722 MiscUtil.PrintError("No value specified for option \"--freezeAtomsParams\". You must specify valid values during, yes, value for \"--freezeAtoms\" option.")
723
724 MiscUtil.ValidateOptionTextValue("-m, --mode", Options["--mode"], "NPT NVT")
725
726 MiscUtil.ValidateOptionTextValue("-p, --platform", Options["--platform"], "CPU CUDA OpenCL Reference")
727
728 MiscUtil.ValidateOptionTextValue("--restraintAtoms", Options["--restraintAtoms"], "yes no")
729 if re.match("^yes$", Options["--restraintAtoms"], re.I):
730 if Options["--restraintAtomsParams"] is None:
731 MiscUtil.PrintError("No value specified for option \"--restraintAtomsParams\". You must specify valid values during, yes, value for \"--restraintAtoms\" option.")
732
733 MiscUtil.ValidateOptionFloatValue("--restraintSpringConstant", Options["--restraintSpringConstant"], {">": 0})
734
735 MiscUtil.ValidateOptionTextValue("--waterBox", Options["--waterBox"], "yes no")
736
737 # Setup a usage string for docopt...
738 _docoptUsage_ = """
739 OpenMMPerformSimulatedAnnealing.py - Perform simulated annealing.
740
741 Usage:
742 OpenMMPerformSimulatedAnnealing.py [--annealingParams <Name,Value,..> ] [--forcefieldParams <Name,Value,..>]
743 [--freezeAtoms <yes or no>] [--freezeAtomsParams <Name,Value,..>] [--integratorParams <Name,Value,..>]
744 [--mode <NVT or NPT>] [--outputParams <Name,Value,..>] [--outfilePrefix <text>]
745 [--overwrite] [--platform <text>] [--platformParams <Name,Value,..>] [--restraintAtoms <yes or no>]
746 [--restraintAtomsParams <Name,Value,..>] [--restraintSpringConstant <number>]
747 [--simulationParams <Name,Value,..>] [--smallMolFile <SmallMolFile>] [--smallMolID <text>]
748 [--systemParams <Name,Value,..>] [--waterBox <yes or no>]
749 [--waterBoxParams <Name,Value,..>] [-w <dir>] -i <infile>
750 OpenMMPerformSimulatedAnnealing.py -h | --help | -e | --examples
751
752 Description:
753 Perform simulated annealing using an NPT or NVT statistical ensemble. You may
754 run a simulation using a macromolecule or a macromolecule in a complex with
755 small molecule. By default, the system is minimized and equilibrated before the
756 simulated annealing.
757
758 The input file must contain a macromolecule already prepared for simulation.
759 The preparation of the macromolecule for a simulation generally involves the
760 following: tasks: Identification and replacement non-standard residues;
761 Addition of missing residues; Addition of missing heavy atoms; Addition of
762 missing hydrogens; Addition of a water box which is optional.
763
764 In addition, the small molecule input file must contain a molecule already
765 prepared for simulation. It must contain appropriate 3D coordinates relative
766 to the macromolecule along with no missing hydrogens.
767
768 You may optionally add a water box and freeze/restraint atoms for the
769 simulation.
770
771 The simulated annealing workflow involves the following steps: Initial heating
772 and equilibration after each change in temperature; Heating and cooling cycles
773 along with equilibration after each change in temperature; Final equilibration.
774 By default, the simulated annealing is performed for a total of 3.35 ns as shown
775 below:
776
777 ... ... ...
778 Simulated annealing (StepSize: 4 fs)
779
780 Initial heating (Start: 0.0 K; End: 300.0 K; Change: 5.0 K)
781 TotalSteps: 305,000; TotalTime: 1.22 ns
782
783 Equilibration after initial heating (Steps: 100,000; Time: 400.00 ps)
784
785 Heating and cooling cycles (NumCycles: 1)
786
787 Heating and cooling cycle 1
788 Heating (Start: 300.0 K; End: 315.0 K; Change: 1.0 K)
789 TotalSteps: 16,000; TotalTime: 64.00 ps
790
791 Equilibration after heating (Steps: 100,000; Time: 400.00 ps)
792
793 Cooling (Start: 315.0 K; End: 300.0 K; Step: 1.0 K)
794 TotalSteps: 16,000; TotalTime: 64.00 ps
795
796 Equilibration after cooling (Steps: 100,000; Time: 400.00 ps)
797
798 Final equilibration after heating and cooling cycles (Steps: 200,000; Time: 800.00 ps)
799
800 Simulated annealing summary (TotalSteps: 837,000; TotalTime: 3.35 ns)...
801 ... ... ...
802
803 The supported macromolecule input file formats are: PDB (.pdb) and
804 CIF (.cif)
805
806 The supported small molecule input file format are : SD (.sdf, .sd)
807
808 Possible outfile prefixes:
809
810 <InfieRoot>_<Mode>
811 <InfieRoot>_Solvated_<Mode>
812 <InfieRoot>_<SmallMolFileRoot>_Complex_<Mode>,
813 <InfieRoot>_<SmallMolFileRoot>_Complex_Solvated_<Mode>
814
815 Possible output files:
816
817 <OutfilePrefix>.<pdb or cif> [ Initial sytem ]
818 <OutfilePrefix>.<dcd or xtc>
819
820 <OutfilePrefix>_Reimaged.<pdb or cif> [ First frame ]
821 <OutfilePrefix>_Reimaged.<dcd or xtc>
822
823 <OutfilePrefix>_Minimized.<pdb or cif>
824 <OutfilePrefix>_Final.<pdb or cif> [ Final system ]
825
826 <OutfilePrefix>.chk
827 <OutfilePrefix>.csv
828 <OutfilePrefix>_Minimization.csv
829 <OutfilePrefix>_FinalState.chk
830 <OutfilePrefix>_FinalState.xml
831
832 <OutfilePrefix>_System.xml
833 <OutfilePrefix>_Integrator.xml
834
835 The reimaged PDB file, <OutfilePrefix>_Reimaged.pdb, corresponds to the first
836 frame in the trajectory. The reimaged trajectory file contains all the frames
837 aligned to the first frame after reimaging of the frames for periodic systems.
838
839 Options:
840 -a, --annealingParams <Name,Value,..> [default: auto]
841 A comma delimited list of parameter name and value pairs for simulated
842 annealing.
843
844 The supported parameter names along with their default values are
845 are shown below:
846
847 Initial heating parameters:
848
849 initialStart, 0.0 [ Units: kelvin ]
850 initialEnd, 300.0 [ Units: kelvin ]
851 initialChange, 5.0 [ Units: kelvin ]
852 initialSteps, 5000
853
854 initialEquilibrationSteps, 100000
855
856 Heating and cooling cycle parameters:
857
858 cycles, 1
859
860 cycleStart, auto [ Units: kelvin. The default value is set to
861 initialEnd ]
862 cycleEnd, 315.0 [ Units: kelvin ]
863 cycleChange, 1.0 [ Units: kelvin ]
864 cycleSteps, 1000
865
866 cycleEquilibrationSteps, 100000
867
868 Final equilibration parameters:
869
870 finalEquilibrationSteps, 200000
871
872 A brief description of parameters is provided below:
873
874 Initial heating parameters:
875
876 initialStart: Start temperature for initial heating.
877 initialEnd: End temperature for initial heating.
878 initialChange: Temperature change for increasing temperature
879 during initial heating.
880 initialSteps: Number of simulation steps after each
881 heating step during initial heating
882
883 initialEquilibrationSteps: Number of equilibration steps
884 after the completion of initial heating.
885
886 Heating and cooling cycles parameters:
887
888 cycles: Number of annealing cycles to perform. Each cycle
889 consists of a heating and a cooling phase. The heating phase
890 consists of the following steps: Heat system from start to
891 end temperature using change size and perform simulation for a
892 number of steps after each increase in temperature; Perform
893 equilibration after the completion of heating. The cooling
894 phase is reverse of the heating phase and cools the system
895 from end to start temperature.
896
897 cycleStart: Start temperature for annealing cycle.
898 cycleEnd: End temperature for annealing cycle.
899 cycleChange: Temperature change for increasing or decreasing
900 temperature during annealing cycle.
901 cycleSteps: Number of simulation steps after each heating and
902 cooling step during annealing cycle.
903
904 cycleEquilibrationSteps: Number of equilibration steps
905 after the completion of heating and cooling phase during a
906 annealing cycle.
907
908 Final equilibration parameters:
909
910 finalEquilibrationSteps: Number of final equilibration
911 steps after the completion of annealing cycles.
912
913 -e, --examples
914 Print examples.
915 -f, --forcefieldParams <Name,Value,..> [default: auto]
916 A comma delimited list of parameter name and value pairs for biopolymer,
917 water, and small molecule forcefields.
918
919 The supported parameter names along with their default values are
920 are shown below:
921
922 biopolymer, amber14-all.xml [ Possible values: Any Valid value ]
923 smallMolecule, openff-2.2.1 [ Possible values: Any Valid value ]
924 water, auto [ Possible values: Any Valid value ]
925 additional, none [ Possible values: Space delimited list of any
926 valid value ]
927
928 Possible biopolymer forcefield values:
929
930 amber14-all.xml, amber99sb.xml, amber99sbildn.xml, amber03.xml,
931 amber10.xml
932 charmm36.xml, charmm_polar_2019.xml
933 amoeba2018.xml
934
935 Possible small molecule forcefield values:
936
937 openff-2.2.1, openff-2.0.0, openff-1.3.1, openff-1.2.1,
938 openff-1.1.1, openff-1.1.0,...
939 smirnoff99Frosst-1.1.0, smirnoff99Frosst-1.0.9,...
940 gaff-2.11, gaff-2.1, gaff-1.81, gaff-1.8, gaff-1.4,...
941
942 The default water forcefield valus is dependent on the type of the
943 biopolymer forcefield as shown below:
944
945 Amber: amber14/tip3pfb.xml
946 CHARMM: charmm36/water.xml or None for charmm_polar_2019.xml
947 Amoeba: None (Explicit)
948
949 Possible water forcefield values:
950
951 amber14/tip3p.xml, amber14/tip3pfb.xml, amber14/spce.xml,
952 amber14/tip4pew.xml, amber14/tip4pfb.xml,
953 charmm36/water.xml, charmm36/tip3p-pme-b.xml,
954 charmm36/tip3p-pme-f.xml, charmm36/spce.xml,
955 charmm36/tip4pew.xml, charmm36/tip4p2005.xml,
956 charmm36/tip5p.xml, charmm36/tip5pew.xml,
957 implicit/obc2.xml, implicit/GBn.xml, implicit/GBn2.xml,
958 amoeba2018_gk.xml (Implict water)
959 None (Explicit water for amoeba)
960
961 The additional forcefield value is a space delimited list of any valid
962 forcefield values and is passed on to the OpenMMForcefields
963 SystemGenerator along with the specified forcefield values for
964 biopolymer, water, and mall molecule. Possible additional forcefield
965 values are:
966
967 amber14/DNA.OL15.xml amber14/RNA.OL3.xml
968 amber14/lipid17.xml amber14/GLYCAM_06j-1.xml
969 ... ... ...
970
971 You may specify any valid forcefield names supported by OpenMM. No
972 explicit validation is performed.
973 --freezeAtoms <yes or no> [default: no]
974 Freeze atoms during a simulation. The specified atoms are kept completely
975 fixed by setting their masses to zero. Their positions do not change during
976 local energy minimization and MD simulation, and they do not contribute
977 to the kinetic energy of the system.
978 --freezeAtomsParams <Name,Value,..>
979 A comma delimited list of parameter name and value pairs for freezing
980 atoms during a simulation. You must specify these parameters for 'yes'
981 value of '--freezeAtoms' option.
982
983 The supported parameter names along with their default values are
984 are shown below:
985
986 selection, none [ Possible values: CAlphaProtein, Ions, Ligand,
987 Protein, Residues, or Water ]
988 selectionSpec, auto [ Possible values: A space delimited list of
989 residue names ]
990 negate, no [ Possible values: yes or no ]
991
992 A brief description of parameters is provided below:
993
994 selection: Atom selection to freeze.
995 selectionSpec: A space delimited list of residue names for
996 selecting atoms to freeze. You must specify its value during
997 'Ligand' and 'Protein' value for 'selection'. The default values
998 are automatically set for 'CAlphaProtein', 'Ions', 'Protein',
999 and 'Water' values of 'selection' as shown below:
1000
1001 CAlphaProtein: List of stadard protein residues from pdbfixer
1002 for selecting CAlpha atoms.
1003 Ions: Li Na K Rb Cs Cl Br F I
1004 Water: HOH
1005 Protein: List of standard protein residues from pdbfixer.
1006
1007 negate: Negate atom selection match to select atoms for freezing.
1008
1009 In addition, you may specify an explicit space delimited list of residue
1010 names using 'selectionSpec' for any 'selection". The specified residue
1011 names are appended to the appropriate default values during the
1012 selection of atoms for freezing.
1013 -h, --help
1014 Print this help message.
1015 -i, --infile <infile>
1016 Input file name containing a macromolecule.
1017 --integratorParams <Name,Value,..> [default: auto]
1018 A comma delimited list of parameter name and value pairs for integrator
1019 during a simulation.
1020
1021 The supported parameter names along with their default values are
1022 are shown below:
1023
1024 integrator, LangevinMiddle [ Possible values: LangevinMiddle,
1025 Langevin, NoseHoover, Brownian ]
1026
1027 randomSeed, auto [ Possible values: > 0 ]
1028
1029 frictionCoefficient, 1.0 [ Units: 1/ps ]
1030 stepSize, auto [ Units: fs; Default value: 4 fs during yes value of
1031 hydrogen mass repartioning with no freezing/restraining of atoms;
1032 otherwsie, 2 fs ]
1033
1034 barostat, MonteCarlo [ Possible values: MonteCarlo or
1035 MonteCarloMembrane ]
1036 barostatInterval, 25
1037 pressure, 1.0 [ Units: atm ]
1038
1039 Parameters used only for MonteCarloMembraneBarostat with default
1040 values corresponding to Amber forcefields:
1041
1042 surfaceTension, 0.0 [ Units: atm*A. It is automatically converted
1043 into OpenMM default units of atm*nm before its usage. ]
1044 xymode, Isotropic [ Possible values: Anisotropic or Isotropic ]
1045 zmode, Free [ Possible values: Free or Fixed ]
1046
1047 A brief description of parameters is provided below:
1048
1049 integrator: Type of integrator
1050
1051 randomSeed: Random number seed for barostat and integrator. Not
1052 supported for NoseHoover integrator.
1053
1054 frictionCoefficient: Friction coefficient for coupling the system to
1055 the heat bath..
1056 stepSize: Simulation time step size.
1057
1058 barostat: Barostat type.
1059 barostatInterval: Barostat interval step size during NPT
1060 simulation for applying Monte Carlo pressure changes.
1061 pressure: Pressure during NPT simulation.
1062
1063 Parameters used only for MonteCarloMembraneBarostat:
1064
1065 surfaceTension: Surface tension acting on the system.
1066 xymode: Behavior along X and Y axes. You may allow the X and Y axes
1067 to vary independently of each other or always scale them by the same
1068 amount to keep the ratio of their lengths constant.
1069 zmode: Beahvior along Z axis. You may allow the Z axis to vary
1070 independently of the other axes or keep it fixed.
1071
1072 -m, --mode <NPT or NVT> [default: NPT]
1073 Type of statistical ensemble to use for simulation. Possible values:
1074 NPT (constant Number of particles, Pressure, and Temperature) or
1075 NVT ((constant Number of particles, Volume and Temperature)
1076 --outfilePrefix <text> [default: auto]
1077 File prefix for generating the names of output files. The default value
1078 depends on the names of input files for macromolecule and small molecule
1079 along with the type of statistical ensemble and the nature of the solvation.
1080
1081 The possible values for outfile prefix are shown below:
1082
1083 <InfieRoot>_<Mode>
1084 <InfieRoot>_Solvated_<Mode>
1085 <InfieRoot>_<SmallMolFileRoot>_Complex_<Mode>,
1086 <InfieRoot>_<SmallMolFileRoot>_Complex_Solvated_<Mode>
1087
1088 --outputParams <Name,Value,..> [default: auto]
1089 A comma delimited list of parameter name and value pairs for generating
1090 output during a simulation.
1091
1092 The supported parameter names along with their default values are
1093 are shown below:
1094
1095 checkpoint, no [ Possible values: yes or no ]
1096 checkpointFile, auto [ Default: <OutfilePrefix>.chk ]
1097 checkpointSteps, 10000
1098
1099 dataOutType, auto [ Possible values: A space delimited list of valid
1100 parameter names.
1101 NPT simulation default: Density Step Speed PotentialEnergy
1102 Temperature Time.
1103 NVT simulation default: Step Speed PotentialEnergy Temperature
1104 Time Volumne
1105 Other valid names: ElapsedTime Progress RemainingTime
1106 KineticEnergy TotalEnergy ]
1107
1108 dataLog, yes [ Possible values: yes or no ]
1109 dataLogFile, auto [ Default: <OutfilePrefix>.csv ]
1110 dataLogSteps, 1000
1111
1112 dataStdout, no [ Possible values: yes or no ]
1113 dataStdoutSteps, 1000
1114
1115 minimizationDataSteps, 100
1116 minimizationDataStdout, no [ Possible values: yes or no ]
1117 minimizationDataLog, no [ Possible values: yes or no ]
1118 minimizationDataLogFile, auto [ Default:
1119 <OutfilePrefix>_MinimizationOut.csv ]
1120 minimizationDataOutType, auto [ Possible values: A space delimited
1121 list of valid parameter names. Default: SystemEnergy
1122 RestraintEnergy MaxConstraintError.
1123 Other valid names: RestraintStrength ]
1124
1125 pdbOutFormat, PDB [ Possible values: PDB or CIF ]
1126 pdbOutKeepIDs, yes [ Possible values: yes or no ]
1127
1128 pdbOutMinimized, no [ Possible values: yes or no ]
1129 pdbOutEquilibrated, no [ Possible values: yes or no ]
1130 pdbOutFinal, no [ Possible values: yes or no ]
1131
1132 saveFinalStateCheckpoint, yes [ Possible values: yes or no ]
1133 saveFinalStateCheckpointFile, auto [ Default:
1134 <OutfilePrefix>_FinalState.chk ]
1135 saveFinalStateXML, no [ Possible values: yes or no ]
1136 saveFinalStateXMLFile, auto [ Default:
1137 <OutfilePrefix>_FinalState.xml]
1138
1139 traj, yes [ Possible values: yes or no ]
1140 trajFile, auto [ Default: <OutfilePrefix>.<TrajFormat> ]
1141 trajFormat, DCD [ Possible values: DCD or XTC ]
1142 trajSteps, 10000 [ The default value corresponds to 40 ps for step
1143 size of 4 fs. ]
1144
1145 xmlSystemOut, no [ Possible values: yes or no ]
1146 xmlSystemFile, auto [ Default: <OutfilePrefix>_System.xml ]
1147 xmlIntegratorOut, no [ Possible values: yes or no ]
1148 xmlIntegratorFile, auto [ Default: <OutfilePrefix>_Integrator.xml ]
1149
1150 A brief description of parameters is provided below:
1151
1152 checkpoint: Write intermediate checkpoint file.
1153 checkpointFile: Intermediate checkpoint file name.
1154 checkpointSteps: Frequency of writing intermediate checkpoint file.
1155
1156 dataOutType: Type of data to write to stdout and log file.
1157
1158 dataLog: Write data to log file.
1159 dataLogFile: Data log file name.
1160 dataLogSteps: Frequency of writing data to log file.
1161
1162 dataStdout: Write data to stdout.
1163 dataStdoutSteps: Frequency of writing data to stdout.
1164
1165 minimizationDataSteps: Frequency of writing data to stdout
1166 and log file.
1167 minimizationDataStdout: Write data to stdout.
1168 minimizationDataLog: Write data to log file.
1169 minimizationDataLogFile: Data log fie name.
1170 minimizationDataOutType: Type of data to write to stdout
1171 and log file.
1172
1173 saveFinalStateCheckpoint: Save final state checkpoint file.
1174 saveFinalStateCheckpointFile: Name of final state checkpoint file.
1175 saveFinalStateXML: Save final state XML file.
1176 saveFinalStateXMLFile: Name of final state XML file.
1177
1178 pdbOutFormat: Format of output PDB files.
1179 pdbOutKeepIDs: Keep existing chain and residue IDs.
1180
1181 pdbOutMinimized: Write PDB file after minimization.
1182 pdbOutEquilibrated: Write PDB file after equilibration.
1183 pdbOutFinal: Write final PDB file after production run.
1184
1185 traj: Write out trajectory file.
1186 trajFile: Trajectory file name.
1187 trajFormat: Trajectory file format.
1188 trajSteps: Frequency of writing trajectory file.
1189
1190 xmlSystemOut: Write system XML file.
1191 xmlSystemFile: System XML file name.
1192 xmlIntegratorOut: Write integrator XML file.
1193 xmlIntegratorFile: Integrator XML file name.
1194
1195 --overwrite
1196 Overwrite existing files.
1197 -p, --platform <text> [default: CPU]
1198 Platform to use for running MD simulation. Possible values: CPU, CUDA,
1199 OpenCL, or Reference.
1200 --platformParams <Name,Value,..> [default: auto]
1201 A comma delimited list of parameter name and value pairs to configure
1202 platform for running MD simulation.
1203
1204 The supported parameter names along with their default values for
1205 different platforms are shown below:
1206
1207 CPU:
1208
1209 threads, 1 [ Possible value: >= 0 or auto. The value of 'auto'
1210 or zero implies the use of all available CPUs for threading. ]
1211
1212 CUDA:
1213
1214 deviceIndex, auto [ Possible values: 0, '0 1' etc. ]
1215 deterministicForces, auto [ Possible values: yes or no ]
1216 precision, single [ Possible values: single, double, or mix ]
1217 tempDirectory, auto [ Possible value: DirName ]
1218 useBlockingSync, auto [ Possible values: yes or no ]
1219 useCpuPme, auto [ Possible values: yes or no ]
1220
1221 OpenCL:
1222
1223 deviceIndex, auto [ Possible values: 0, '0 1' etc. ]
1224 openCLPlatformIndex, auto [ Possible value: Number]
1225 precision, single [ Possible values: single, double, or mix ]
1226 useCpuPme, auto [ Possible values: yes or no ]
1227
1228 A brief description of parameters is provided below:
1229
1230 CPU:
1231
1232 threads: Number of threads to use for simulation.
1233
1234 CUDA:
1235
1236 deviceIndex: Space delimited list of device indices to use for
1237 calculations.
1238 deterministicForces: Generate reproducible results at the cost of a
1239 small decrease in performance.
1240 precision: Number precision to use for calculations.
1241 tempDirectory: Directory name for storing temporary files.
1242 useBlockingSync: Control run-time synchronization between CPU and
1243 GPU.
1244 useCpuPme: Use CPU-based PME implementation.
1245
1246 OpenCL:
1247
1248 deviceIndex: Space delimited list of device indices to use for
1249 simulation.
1250 openCLPlatformIndex: Platform index to use for calculations.
1251 precision: Number precision to use for calculations.
1252 useCpuPme: Use CPU-based PME implementation.
1253
1254 --restraintAtoms <yes or no> [default: no]
1255 Restraint atoms during a simulation. The motion of specified atoms is
1256 restricted by adding a harmonic force that binds them to their starting
1257 positions. The atoms are not completely fixed unlike freezing of atoms.
1258 Their motion, however, is restricted and they are not able to move far away
1259 from their starting positions during local energy minimization and MD
1260 simulation.
1261 --restraintAtomsParams <Name,Value,..>
1262 A comma delimited list of parameter name and value pairs for restraining
1263 atoms during a simulation. You must specify these parameters for 'yes'
1264 value of '--restraintAtoms' option.
1265
1266 The supported parameter names along with their default values are
1267 are shown below:
1268
1269 selection, none [ Possible values: CAlphaProtein, Ions, Ligand,
1270 Protein, Residues, or Water ]
1271 selectionSpec, auto [ Possible values: A space delimited list of
1272 residue names ]
1273 negate, no [ Possible values: yes or no ]
1274
1275 A brief description of parameters is provided below:
1276
1277 selection: Atom selection to restraint.
1278 selectionSpec: A space delimited list of residue names for
1279 selecting atoms to restraint. You must specify its value during
1280 'Ligand' and 'Protein' value for 'selection'. The default values
1281 are automatically set for 'CAlphaProtein', 'Ions', 'Protein',
1282 and 'Water' values of 'selection' as shown below:
1283
1284 CAlphaProtein: List of stadard protein residues from pdbfixer
1285 for selecting CAlpha atoms.
1286 Ions: Li Na K Rb Cs Cl Br F I
1287 Water: HOH
1288 Protein: List of standard protein residues from pdbfixer.
1289
1290 negate: Negate atom selection match to select atoms for freezing.
1291
1292 In addition, you may specify an explicit space delimited list of residue
1293 names using 'selectionSpec' for any 'selection". The specified residue
1294 names are appended to the appropriate default values during the
1295 selection of atoms for restraining.
1296 --restraintSpringConstant <number> [default: 2.5]
1297 Restraint spring constant for applying external restraint force to restraint
1298 atoms relative to their initial positions during 'yes' value of '--restraintAtoms'
1299 option. Default units: kcal/mol/A**2. The default value, 2.5, corresponds to
1300 1046.0 kjoules/mol/nm**2. The default value is automatically converted into
1301 units of kjoules/mol/nm**2 before its usage.
1302 --simulationParams <Name,Value,..> [default: auto]
1303 A comma delimited list of parameter name and value pairs for simulation.
1304
1305 The supported parameter names along with their default values are
1306 are shown below:
1307
1308 minimization, yes [ Possible values: yes or no ]
1309 minimizationMaxSteps, auto [ Possible values: >= 0. The value of
1310 zero implies until the minimization is converged. ]
1311 minimizationTolerance, 0.24 [ Units: kcal/mol/A. The default value
1312 0.25, corresponds to OpenMM default of value of 10.04
1313 kjoules/mol/nm. It is automatically converted into OpenMM
1314 default units before its usage. ]
1315
1316 A brief description of parameters is provided below:
1317
1318 minimization: Perform minimization before equilibration and
1319 production run.
1320 minimizationMaxSteps: Maximum number of minimization steps. The
1321 value of zero implies until the minimization is converged.
1322 minimizationTolerance: Energy convergence tolerance during
1323 minimization.
1324
1325 -s, --smallMolFile <SmallMolFile>
1326 Small molecule input file name. The macromolecue and small molecule are
1327 merged for simulation and the complex is written out to a PDB file.
1328 --smallMolID <text> [default: LIG]
1329 Three letter small molecule residue ID. The small molecule ID corresponds
1330 to the residue name of the small molecule and is written out to a PDB file
1331 containing the complex.
1332 --systemParams <Name,Value,..> [default: auto]
1333 A comma delimited list of parameter name and value pairs to configure
1334 a system for simulation.
1335
1336 The supported parameter names along with their default values are
1337 are shown below:
1338
1339 constraints, BondsInvolvingHydrogens [ Possible values: None,
1340 WaterOnly, BondsInvolvingHydrogens, AllBonds, or
1341 AnglesInvolvingHydrogens ]
1342 constraintErrorTolerance, 0.000001
1343 ewaldErrorTolerance, 0.0005
1344
1345 nonbondedMethodPeriodic, PME [ Possible values: NoCutoff,
1346 CutoffNonPeriodic, or PME ]
1347 nonbondedMethodNonPeriodic, NoCutoff [ Possible values:
1348 NoCutoff or CutoffNonPeriodic]
1349 nonbondedCutoff, 1.0 [ Units: nm ]
1350
1351 hydrogenMassRepartioning, yes [ Possible values: yes or no ]
1352 hydrogenMass, 1.5 [ Units: amu]
1353
1354 removeCMMotion, yes [ Possible values: yes or no ]
1355 rigidWater, auto [ Possible values: yes or no. Default: 'No' for
1356 'None' value of constraints; Otherwise, yes ]
1357
1358 A brief description of parameters is provided below:
1359
1360 constraints: Type of system constraints to use for simulation. These
1361 constraints are different from freezing and restraining of any
1362 atoms in the system.
1363 constraintErrorTolerance: Distance tolerance for constraints as a
1364 fraction of the constrained distance.
1365 ewaldErrorTolerance: Ewald error tolerance for a periodic system.
1366
1367 nonbondedMethodPeriodic: Nonbonded method to use during the
1368 calculation of long range interactions for a periodic system.
1369 nonbondedMethodNonPeriodic: Nonbonded method to use during the
1370 calculation of long range interactions for a non-periodic system.
1371 nonbondedCutoff: Cutoff distance to use for long range interactions
1372 in both perioidic non-periodic systems.
1373
1374 hydrogenMassRepartioning: Use hydrogen mass repartioning. It
1375 increases the mass of the hydrogen atoms attached to the heavy
1376 atoms and decreasing the mass of the bonded heavy atom to
1377 maintain constant system mass. This allows the use of larger
1378 integration step size (4 fs) during a simulation.
1379 hydrogenMass: Hydrogen mass to use during repartioning.
1380
1381 removeCMMotion: Remove all center of mass motion at every time step.
1382 rigidWater: Keep water rigid during a simulation. This is determined
1383 automatically based on the value of 'constraints' parameter.
1384
1385 --waterBox <yes or no> [default: no]
1386 Add water box.
1387 --waterBoxParams <Name,Value,..> [default: auto]
1388 A comma delimited list of parameter name and value pairs for adding
1389 a water box.
1390
1391 The supported parameter names along with their default values are
1392 are shown below:
1393
1394 model, tip3p [ Possible values: tip3p, spce, tip4pew, tip5p or
1395 swm4ndp ]
1396 mode, Padding [ Possible values: Size or Padding ]
1397 padding, 1.0
1398 size, None [ Possible value: xsize ysize zsize ]
1399 shape, cube [ Possible values: cube, dodecahedron, or octahedron ]
1400 ionPositive, Na+ [ Possible values: Li+, Na+, K+, Rb+, or Cs+ ]
1401 ionNegative, Cl- [ Possible values: Cl-, Br-, F-, or I- ]
1402 ionicStrength, 0.0
1403
1404 A brief description of parameters is provided below:
1405
1406 model: Water model to use for adding water box. The van der
1407 Waals radii and atomic charges are determined using the
1408 specified water forcefield. You must specify an appropriate
1409 water forcefield. No validation is performed.
1410 mode: Specify the size of the waterbox explicitly or calculate it
1411 automatically for a macromolecule along with adding padding
1412 around ther macromolecule.
1413 padding: Padding around a macromolecule in nanometers for filling
1414 box with water. It must be specified during 'Padding' value of
1415 'mode' parameter.
1416 size: A space delimited triplet of values corresponding to water
1417 size in nanometers. It must be specified during 'Size' value of
1418 'mode' parameter.
1419 ionPositive: Type of positive ion to add during the addition of a
1420 water box.
1421 ionNegative: Type of negative ion to add during the addition of a
1422 water box.
1423 ionicStrength: Total concentration of both positive and negative
1424 ions to add excluding the ions added to neutralize the system
1425 during the addition of a water box.
1426
1427 -w, --workingdir <dir>
1428 Location of working directory which defaults to the current directory.
1429
1430 Examples:
1431 To perform simulated annealing for a macromolecule in a PDB file by using an NPT
1432 ensemble, applying system constraints for bonds involving hydrogens along with
1433 hydrogen mass repartioning, using a step size of 4 fs, performing minimization
1434 until it's converged, performing initial heating along with equilibration for
1435 305,000 steps (1.22 ns), performing one heating and cooling cycle along with
1436 equilibration 232,000 steps (928 ps), performing final equilibration for
1437 100,000 steps (400 ps), writing trajectory file every 10,000 steps (40 ps),
1438 writing data log file every 1,000 steps (4 ps), generating a checkpoint file
1439 file after the completion of the, and generating a PDB for the final system,
1440 type:
1441
1442 % OpenMMPerformSimulatedAnnealing.py -i Sample13.pdb --waterBox yes
1443
1444 To run the first example for performing OpenMM simulation using multi-
1445 threading employing all available CPUs on your machine and generate various
1446 output files, type:
1447
1448 % OpenMMPerformSimulatedAnnealing.py -i Sample13.pdb --waterBox yes
1449 --platformParams "threads,0"
1450
1451 To run the first example for performing OpenMM simulation using CUDA platform
1452 on your machine and generate various output files, type:
1453
1454 % OpenMMPerformSimulatedAnnealing.py -i Sample13.pdb --waterBox yes
1455 -p CUDA
1456
1457 To run the second example for a marcomolecule in a complex with a small
1458 molecule and generate various output files, type:
1459
1460 % OpenMMPerformSimulatedAnnealing.py -i Sample13.pdb
1461 -s Sample13Ligand.sdf --waterBox yes --platformParams "threads,0"
1462
1463 To run the second example for performing NVT simulation and generate various
1464 output files, type:
1465
1466 % OpenMMPerformSimulatedAnnealing.py -i Sample13.pdb
1467 -s Sample13Ligand.sdf --mode NVT --platformParams "threads,0"
1468
1469 % OpenMMPerformSimulatedAnnealing.py -i Sample13.pdb
1470 -s Sample13Ligand.sdf --mode NVT --waterBox yes
1471 --platformParams "threads,0"
1472
1473 To run the second example by freezing CAlpha atoms in a macromolecule without
1474 using any system constraints to avoid any issues with the freezing of the same atoms,
1475 using a step size of 2 fs, and generate various output files, type:
1476
1477 % OpenMMPerformSimulatedAnnealing.py -i Sample13.pdb --waterBox yes
1478 --freezeAtoms yes --freezeAtomsParams "selection,CAlphaProtein"
1479 --systemParams "constraints, None"
1480 --platformParams "threads,0" --integratorParams "stepSize,2"
1481
1482 To run the second example by restrainting CAlpha atoms in a macromolecule and
1483 and generate various output files, type:
1484
1485 % OpenMMPerformSimulatedAnnealing.py -i Sample13.pdb --waterBox yes
1486 --restraintAtomsParams "selection,CAlphaProtein"
1487 --platformParams "threads,0" --integratorParams "stepSize,2"
1488
1489 To run the second example by specifying explict values for various parametres
1490 and generate various output files, type:
1491
1492 % OpenMMPerformSimulatedAnnealing.py -m NPT -i Sample13.pdb
1493 --annealingParams "initialStart, 0.0, initialEnd, 300.0, initialChange, 5.0,
1494 initialSteps, 5000, initialEquilibrationSteps, 100000, cycles, 1,
1495 cycleStart, 300.0, cycleEnd, 315.0, cycleChange, 1.0,
1496 cycleSteps, 1000, finalEquilibrationSteps, 200000"
1497 -f " biopolymer,amber14-all.xml,smallMolecule, openff-2.2.1,
1498 water,amber14/tip3pfb.xml" --waterBox yes
1499 --integratorParams "integrator,LangevinMiddle,randomSeed,42,
1500 stepSize,2,,pressure, 1.0"
1501 --outputParams "checkpoint,yes,dataLog,yes,dataStdout,yes,
1502 minimizationDataStdout,yes,minimizationDataLog,yes,
1503 pdbOutFormat,CIF,pdbOutKeepIDs,yes,saveFinalStateCheckpoint, yes,
1504 traj,yes,xmlSystemOut,yes,xmlIntegratorOut,yes"
1505 -p CPU --platformParams "threads,0"
1506 --simulationParams "sminimization,yes, minimizationMaxSteps,
1507 500,equilibration,yes,"
1508 --systemParams "constraints,BondsInvolvingHydrogens,
1509 nonbondedMethodPeriodic,PME,nonbondedMethodNonPeriodic,NoCutoff,
1510 hydrogenMassRepartioning, yes"
1511
1512 Author:
1513 Manish Sud(msud@san.rr.com)
1514
1515 See also:
1516 OpenMMPrepareMacromolecule.py, OpenMMPerformMDSimulation.py,
1517 OpenMMPerformMinimization.py
1518
1519 Copyright:
1520 Copyright (C) 2025 Manish Sud. All rights reserved.
1521
1522 The functionality available in this script is implemented using OpenMM, an
1523 open source molecuar simulation package.
1524
1525 This file is part of MayaChemTools.
1526
1527 MayaChemTools is free software; you can redistribute it and/or modify it under
1528 the terms of the GNU Lesser General Public License as published by the Free
1529 Software Foundation; either version 3 of the License, or (at your option) any
1530 later version.
1531
1532 """
1533
1534 if __name__ == "__main__":
1535 main()
