OpenSMOG

The OpenSMOG classes perform molecular dynamics using Structure-Based Models (SBM) for biomolecular simulations. OpenSMOG uses force fields generated by SMOG 2, and it allows the simulations of a wide variety of potential forms, including commonly employed C-alpha and all-atom variants. Details about the default models in SMOG 2 can be found in the following resources:

  • SMOG server: https://smog-server.org/smog2/

  • C-alpha: Clementi, C., Nymeyer, H. and Onuchic, J.N., 2000. Topological and energetic factors: what determines the structural details of the transition state ensemble and “en-route” intermediates for protein folding? An investigation for small globular proteins. Journal of molecular biology, 298(5), pp.937-953.

  • All-Atom: Whitford, P.C., Noel, J.K., Gosavi, S., Schug, A., Sanbonmatsu, K.Y. and Onuchic, J.N., 2009. An all‐atom structure‐based potential for proteins: bridging minimal models with all‐atom empirical forcefields. Proteins: Structure, Function, and Bioinformatics, 75(2), pp.430-441.

class OpenSMOG.SBM(time_step, collision_rate, r_cutoff, temperature=- 1, cmm=True, pbc=False, name='OpenSMOG', warn=True)[source]

Bases: object

createReporters(trajectory=True, trajectoryName=None, trajectoryFormat='dcd', energies=True, energiesName=None, energy_components=False, energy_componentsName=None, logFileName='OpenSMOG.log', interval=1000)[source]

Creates the reporters to provide the output data.

Parameters
  • logFileName (str, optional) – Name of log file. (Default value: OpenSMOG.log).

  • trajectory (bool, optional) – Whether to save the trajectory .dcd file containing the position of the atoms as a function of time. (Default value: True).

  • trajectoryName (str, optional) – Name of the trajectory file.

  • trajectoryFormat (str, optional) – File format of the trajectory file. Options are dcd, pdb, pdbx, hdf5, netcdf, and xtc. Saving the trajectory in the file formats hdf5, netcdf, and xtc require the MDTraj package. (Default value: dcd).

  • energies (bool, optional) – Whether to save the energies in a .txt file containing five columns, comma-delimited. The header of the files shows the information of each collum: #”Step”,”Potential Energy (kJ/mole)”,”Kinetic Energy (kJ/mole)”,”Total Energy (kJ/mole)”,”Temperature (K)”. (Default value: True).

  • energiesName (str, optional) – Name of the energy file.

  • energy_components (bool, optional) – Whether to save the potential energy for each applied force in a .txt file containing several columns, comma-delimited. The header of the files shows the information of each column. An example of the header is: #”Step”,”electrostatic”,”Non-Contacts”,”Bonds”,”Angles”,”Dihedrals”,”contact_1-10-12”. (Default value: False).

  • energy_componentsName (str, optional) – Name of the energy_components file.

  • logFileName – Name of the log file.

  • interval (int, required) – Frequency to write the data to the output files. (Default value: 10**3)

createSimulation()[source]

Creates the simulation context and loads into the OpenMM platform.

help()[source]

Prints information about using OpenSMOG.

loadGro(Grofile)[source]

Loads the .gro file format in the OpenMM system platform. The input files are generated using SMOG2 software with the flag -OpenSMOG. Details on how to create the files can be found in the SMOG2 User Manual. Tutorial on how to generate the inputs files can be found on the smog-server page, or here.

Parameters

Grofile (file, required) – Initial structure for the MD simulations in .gro file format. This is typically generated by SMOG2 with the flag -OpenSMOG. (Default value: None).

loadSystem(Grofile, Topfile, Xmlfile)[source]

Loads the input files in the OpenMM system platform. The input files are generated using SMOG2 software with the flag -OpenSMOG. Details on how to create the files can be found in the SMOG2 User Manual. Tutorials for how to generate the inputs files can be found on the smog-server page, or here.

Parameters
  • Grofile (file, required) – Initial structure for the MD simulations in .gro file format generated by SMOG2 software with the flag -OpenSMOG. (Default value: None).

  • Topfile (file, required) – Topology .top file format generated by SMOG2 software with the flag -OpenSMOG. The topology file lists the interactions between the system atoms except for the “Native Contacts” potential that is provided to OpenSMOG in a .xml file. (Default value: None).

  • Xmlfile (file, required) – The .xml file can contain the information that defines the “Contact” and nonbonded potentials. The .xml file is generated by SMOG2 software with the flag -OpenSMOG, which support custom potential functions. (Default value: None).

loadTop(Topfile)[source]

Loads the .top file format in the OpenMM system platform. The input files are generated using SMOG2 software with the flag -OpenSMOG. Details on how to create the files can be found in the SMOG2 User Manual. Tutorial on how to generate the inputs files can be found on the smog-server page, or here.

Parameters

Topfile (file, required) – Topology .top file format generated by SMOG2 software with the flag -OpenSMOG. The topology file defines the interactions between atoms, but the “Contacts” and nonbonded potentials can be provided with the .xml file. (Default value: None).

loadXml(Xmlfile)[source]

Loads the .xml file format in the OpenMM system platform. The input files are generated using SMOG2 software with the flag -OpenSMOG. Details on how to create the files can be found in the SMOG2 User Manual. Tutorial on how to generate the inputs files can be found on the smog-server page, or here.

Parameters

Xmlfile (file, required) – The .xml file can contain information that defines the “Contact” and nonbonded potentials. The .xml file is generated by SMOG2 software with the flag -OpenSMOG, which support custom potential energy functions. (Default value: None).

minimize(tolerance=1.0, maxIterations=0)[source]

Wrapper for L-BFGS energy minimization.

Parameters
  • tolerance (float, required) – Stopping criteria value between iterations. When the error between iteration is below this value, the minimization stops. (Default value: 1.0).

  • maxIteration (int, required) – Number of maximum steps to be performed in the minimization simulation. (Default value: 0).

opensmog_quit()[source]
opensmogcheck()[source]
run(nsteps, report=True, interval=10000)[source]
runForClockTime(time, report=True, interval=10000, checkpointFile=None, stateFile=None, checkpointInterval=None)[source]
saveFolder(folder)[source]

Sets the folder path to save data.

Parameters

folder (str, required) – Folder path to save the simulation data. If the folder path does not exist, the function will create the directory.

setup_openmm(platform='opencl', precision='single', GPUindex='default', integrator='langevin')[source]

Sets up the parameters of the simulation OpenMM platform.

Parameters
  • platform (str, optional) – Platform to use in the simulations. Opitions are CUDA, OpenCL, HIP, CPU, Reference. (Default value: OpenCL).

  • precision (str, optional) – Numerical precision type of the simulation. Options are single, mixed, double. (Default value: single). For details check the OpenMM Documentation.

  • GPUIndex (str, optional) – Set of Platform device index IDs. Ex: 0,1,2 for the system to use the devices 0, 1 and 2.

  • integrator (str, or integrator object) – Integrator to use in the simulations. Options are langevin, langevinMiddle, *variableLangevin and, brownian. You may also build your own integrator object and pass it, rather than use a named integrator. (Default value: langevin).

version = '1.1.1'

//smog-server.org

The SBM sets the environment to start the molecular dynamics simulations.

Parameters
  • time_step (float, required) – Simulation time step in units of \(\tau\).

  • collision_rate (float, required) – Friction/Damping constant in units of reciprocal time (\(1/\tau\)).

  • r_cutoff (float, required) – Cutoff distance to consider non-bonded interactions in units of nanometers.

  • pbc (boolean, optional) – Turn PBC on/off. Default value: False

  • cmm (boolean, optional) – Remove center of mass motion. Default value: True

  • temperature (float,required) – Temperature in reduced units

  • name (str) – Name used in the output files. (Default value: OpenSMOG)

  • warn (boolean, optional) – Give cautionary warnings… (Default value: True)

Type

The SBM class performs Molecular dynamics simulations using structure-based (SMOG) models to investigate a broad range of biomolecular dynamics, including domain rearrangements in proteins, folding and ligand binding in RNA, and large-scale rearrangements in ribonucleoprotein assemblies. In its simplest form, a structure-based model defines a particular structure (usually obtained from X-ray, cryo-EM, or NMR methods) as the energetic global minimum. Find more information about SMOG models and OpenSMOG at http