parmed.openmm.parameters module

This module contains the class for storing and creating/converting/writing OpenMM-style ffxml files defining a force field

Author(s): Jason Swails

class parmed.openmm.parameters.OpenMMParameterSet(*filenames)[source]

Bases: parmed.parameters.ParameterSet, parmed.charmm.parameters.CharmmImproperMatchingMixin

Class storing parameters from an OpenMM parameter set

Parameters
filenamesstr, list of str, file-like, or list of file-like; optional

Either the name of a file or a list of filenames from which parameters should be parsed.

Notes

Order is important in the list of files provided. The parameters are loaded in the order they are provided, and any parameters that are specified in multiple places are overwritten (that is, the last occurrence is the parameter type that is used)

Attributes
combining_rule

Methods

condense([do_dihedrals])

This function goes through each of the parameter type dicts and eliminates duplicate types.

from_parameterset(params[, copy, …])

Instantiates an OpenMMParameterSet from another ParameterSet (or subclass).

from_structure(struct[, …])

Extracts known parameters from a Structure instance

id_format(filename)

Identifies the file type as either an Amber-style frcmod or parm.dat file.

match_improper_type(a1, a2, a3, a4)

Matches an improper type based on atom type names

typeify_templates()

Assign atom types to atom names in templates

write(dest[, provenance, write_unused, …])

Write the parameter set to an XML file for use with OpenMM

classmethod from_parameterset(params, copy=False, remediate_residues=True, unique_atom_types=False)[source]

Instantiates an OpenMMParameterSet from another ParameterSet (or subclass). The main thing this feature is responsible for is converting lower-case atom type names into all upper-case and decorating the name to ensure each atom type name is unique.

Parameters
paramsparmed.parameters.ParameterSet

ParameterSet containing the list of parameters to be converted to as OpenMM-compatible parameter set

copybool, optional, default=False

If True, the returned parameter set is a deep copy of params. If False, the returned parameter set is a shallow copy. Default False.

remediate_residuesbool, optional, default=True

If True, will remove non-chemical bonds and drop Residue definitions that are missing parameters

unique_atom_typesbool

If True, a unique OpenMM atom type will be created for every atom of every residue. In this case, the AtomType objects correspond to atom classes rather than atom types (in the OpenMM terminology).

Returns
new_paramsOpenMMParameterSet

OpenMMParameterSet with the same parameters as that defined in the input parameter set

static id_format(filename)[source]

Identifies the file type as either an Amber-style frcmod or parm.dat file.

Parameters
filenamestr

Name of the file to check format for

Returns
is_fmtbool

True if it is an Amber-style parameter file. False otherwise.

write(dest, provenance=None, write_unused=True, separate_ljforce=False, improper_dihedrals_ordering='default', charmm_imp=False, skip_duplicates=True)[source]

Write the parameter set to an XML file for use with OpenMM

Parameters
deststr or file-like

The name of the file or the file-like object (with a write attribute) to which the XML file will be written

provenancedict, optional

If present, the XML file will be tagged with the available fields. Keys of the dictionary become XML etree.Element tags, the values of the dictionary must be instances of any of: - str / unicode (Py2) or str (Py3) - one XML element with this content is written - list - one XML element per each item of the list is written, all these XML elements use the same tag (key in provenance dict) - dict - one of the keys of this dict must be the same as the key of of the provenance dict under which this dict is nested. The value under this key becomes the content of the XML element. Remaining keys and their values are used to construct attributes of the XML element. Note that OrderedDict’s should be used to ensure appropriate order of the XML elements and their attributes. Default is no provenance. Example (unordered): provenance = {‘Reference’ : [‘Nature’, ‘Cell’],

‘Source’ : {‘Source’: ‘leaprc.ff14SB’, sourcePackage : ‘AmberTools’, sourcePackageVersion : ‘15’}, ‘User’ : ‘Mark’}

write_unusedbool

If False: a) residue templates using unavailable atom types will not be written, b) atom types that are not used in any of the residue templates remaining and parameters including those atom types will not be written. A ParameterWarning is issued if any such residues are found in a).

separate_ljforcebool

If True will use a separate LennardJonesForce to create a CostumNonbondedForce to compute L-J interactions. It will set sigma to 1 and epsilon to 0 in the NonbondedForce so that the NonbondedForce only calculates the electrostatic contribution. It should be set to True when converting a CHARMM force field file that doesn’t have pair-specific L-J modifications (NBFIX in CHARMM) so that the ffxml conversion is compatible with the main charmm36.xml file. Note: —- When pair-specific L-J modifications are present (NBFIX in CHARMM), this behavior is always present and this flag is ignored.

improper_dihedrals_orderingstr

The ordering to use when assigning improper torsions in OpenMM. Default is ‘default’, other option is ‘amber’

charmm_imp: bool

If True, will check for existence of IMPR in each residue and patch template, and write out the explicit improper definition without wildcards in the ffxml file.

skip_duplicatesbool

If True: residues which appear identical to an existing residue will not be written. This is usually the best choice. The most common reason for setting skip_duplicates to false is if you have different parametrizations for stereoisomers. Note that since OpenMM’s residue hashing is not aware of chirality, if you wish to use the results in simulations you will need to explicitly provide the template names for affected residues.

Notes

The generated XML file will have the XML tag DateGenerated added to the provenance information set to the current date. Therefore, you should not provide this information in provenance (it will be removed if it is provided).

parmed.openmm.parameters.needs_lxml(func)[source]

Decorator to raise an ImportError if a function requires lxml but it is not present

parmed.openmm.parameters.pretty_print(tree)