parmed package

Subpackages

Module contents

The parmed package manipulates molecular structures and provides a way to go between standard and amber file formats, manipulate structures, etc.

class parmed.AcceptorDonor(atom1, atom2)[source]

Bases: object

Just a holder for donors and acceptors in CHARMM speak

Parameters
atom1Atom

First atom in the donor/acceptor group

atom2Atom

Second atom in the donor/acceptor group

class parmed.AmoebaNonbondedExceptionType(vdw_weight, multipole_weight, direct_weight, polar_weight, mutual_weight, list=None)[source]

Bases: parmed.topologyobjects.NonbondedExceptionType

A parameter describing how the various nonbonded interactions between a particular pair of atoms is scaled in the AMOEBA force field

Parameters
vdw_weightfloat

The scaling factor by which van der Waals interactions are multiplied

multipole_weightfloat

The scaling factor by which multipole interactions are multiplied

direct_weightfloat

The scaling factor by which direct-space interactions are multiplied

polar_weightfloat

The scaling factor by which polarization interactions are multiplied

mutual_weightfloat

The scaling factor by which mutual interactions are multiplied

listTrackedList

The list containing this nonbonded exception

Attributes
idx
sigma
uepsilon
urmin
usigma
class parmed.Angle(atom1, atom2, atom3, type=None)[source]

Bases: object

A valence angle between 3 atoms separated by two covalent bonds.

Parameters
atom1Atom

An atom one end of the valence angle

atom2Atom

The atom in the middle of the valence angle bonded to both other atoms

atom3Atom

An atom on the other end of the valence angle to atom1

typeAngleType

The AngleType object containing the parameters for this angle

Notes

An Angle can contain bonds or atoms. A bond is contained if it exists between atoms 1 and 2 or atoms 2 and 3.

Examples

>>> a1, a2, a3 = Atom(), Atom(), Atom()
>>> angle = Angle(a1, a2, a3)
>>> Bond(a1, a2) in angle and Bond(a3, a2) in angle
True
>>> a1 in angle and a2 in angle and a3 in angle
True
>>> Bond(a1, a3) in angle # this is not part of the angle definition
False

Methods

delete()

Deletes this angle from the atoms that make it up.

energy()

Measures the current angle energy

measure()

Measures the current angle

uenergy()

Same as energy(), but with units

umeasure()

Same as measure(), but with units

delete()[source]

Deletes this angle from the atoms that make it up. This method removes this angle from the angles list for atom1, atom2, and atom3 as well as removing each atom form each others’ angle partner lists

energy()[source]

Measures the current angle energy

Returns
energyfloat or None

Angle strain energy in kcal/mol. Return value is None if either the coordinates of either atom is not set or the angle type is not set

measure()[source]

Measures the current angle

Returns
measurementfloat or None

If the atoms have coordinates, returns the angle between the three atoms. If any coordinates are missing, returns None

uenergy()[source]

Same as energy(), but with units

umeasure()[source]

Same as measure(), but with units

class parmed.AngleType(k, theteq, list=None)[source]

Bases: parmed.topologyobjects._ParameterType, parmed.topologyobjects._ListItem

An angle type with a set of angle parameters

Parameters
kfloat

Force constant in kcal/mol/radians^2

theteqfloat

Equilibrium angle in Degrees

listTrackedList

A list of `AngleType`s in which this is a member

Notes

Two AngleType`s are equal if their `k and theteq attributes are equal

Examples

>>> at1 = AngleType(10.0, 180.0)
>>> at2 = AngleType(10.0, 180.0)
>>> at1 is at2
False
>>> at1 == at2
True
>>> at1.idx # not part of any list or iterable
-1

As part of a list, they can be indexed

>>> angle_list = []
>>> angle_list.append(AngleType(10.0, 180.0, list=angle_list))
>>> angle_list.append(AngleType(10.0, 180.0, list=angle_list))
>>> angle_list[0].idx
0
>>> angle_list[1].idx
1
Attributes
idx
uk
utheteq
property uk
property utheteq
class parmed.Atom(list=None, atomic_number=0, name='', type='', charge=None, mass=0.0, nb_idx=0, solvent_radius=0.0, screen=0.0, tree='BLA', join=0.0, irotat=0.0, occupancy=0.0, bfactor=0.0, altloc='', number=- 1, rmin=None, epsilon=None, rmin14=None, epsilon14=None)[source]

Bases: parmed.topologyobjects._ListItem

An atom. Only use these as elements in AtomList instances, since AtomList will keep track of when indexes and other stuff needs to be updated. All parameters are optional.

Parameters
atomic_numberint

The atomic number of this atom

namestr

The name of this atom

typestr

The type name of this atom

chargefloat

The partial atomic charge of this atom in fractions of an electron

massfloat

The atomic mass of this atom in daltons

nb_idxint

The nonbonded index. This is a pointer that is relevant in the context of an Amber topology file and identifies its Lennard-Jones atom type

solvent_radiusfloat

The intrinsic solvation radius of this atom.

screenfloat

The Generalized Born screening factor for this atom.

occupancyfloat

The occupancy of the atom (see PDB file)

bfactorfloat

The B-factor of the atom (see PDB file)

altlocstr

Alternate location indicator (see PDB file)

Other Parameters
listAtomList

The AtomList that this atom belongs to. If None, this atom does not belong to any list. This can be any iterable, but should typically be an AtomList or None

treestr

The tree chain identifier assigned to this atom. Relevant in the context of an Amber topology file, and not used for very much.

joinint

The ‘join` property of atoms stored in the Amber topology file. At the time of writing this class, join is unused, but still oddly required. Add support for future-proofing

irotatint

The irotat property of atoms stored in the Amber topology file. Unused, but included for future-proofing.

numberint

The serial number given to the atom (see PDB file)

rminfloat

The Rmin/2 Lennard-Jones parameter for this atom. Default evaluates to 0

epsilonfloat

The epsilon (well depth) Lennard-Jones parameter for this atom. Default evaluates to 0

rmin14float

The Rmin/2 Lennard-Jones parameter for this atom in 1-4 interactions. Default evaluates to 0

epsilon14float

The epsilon (well depth) Lennard-Jones parameter for this atom in 1-4 interactions. Default evaluates to 0

Notes

The bond_partners, angle_partners, dihedral_partners, and exclusion_partners arrays are actually generated as properties by taking differences of sets and sorting them. As a result, accessing this attribute constantly can be less efficient than you would expect. Iterating over them in a loop requires minimal overhead. But if frequent access is needed and these sets are guaranteed not to change, you should save a reference to the object and use that instead.

Binary comparisons are done by atom index and are used primarily for sorting sublists of atoms. The == operator is not defined, so Atom equality should be done using the is operator. This allows Atom instances to be hashable (and so used as dict keys and put in `set`s)

Examples

>>> a1 = Atom(name='CO', type='C', charge=0.5, mass=12.01)
>>> a2 = Atom(name='OC', type='O', charge=-0.5, mass=12.01)
>>> a1.bond_to(a2)
>>> a1 in a2.bond_partners and a2 in a1.bond_partners
True
>>> a1.idx # Not part of a container
-1

This object also supports automatic indexing when it is part of a container

>>> atom_list = []
>>> atom_list.append(Atom(list=atom_list, name='CO', charge=0.5))
>>> atom_list.append(Atom(list=atom_list, name='OC', charge=-0.5))
>>> atom_list[0].idx
0
>>> atom_list[1].idx
1
Attributes
angle_partners

List of all angle partners that are NOT bond partners

bond_partners

Go through all bonded partners

charge
dihedral_partners

List of all dihedral partners that are NOT angle or bond partners

element
element_name
epsilon

Lennard-Jones epsilon parameter (the Lennard-Jones well depth)

epsilon_14

The 1-4 Lennard-Jones epsilon parameter

exclusion_partners

List of all exclusions not otherwise excluded by bonds/angles/torsions

idx
rmin

Lennard-Jones Rmin/2 parameter (the Lennard-Jones radius)

rmin_14

The 1-4 Lennard-Jones Rmin/2 parameter

sigma

Lennard-Jones sigma parameter – directly related to Rmin

sigma_14

Lennard-Jones sigma parameter – directly related to Rmin

tortor_partners

List of all 1-5 partners that are NOT in angle or bond partners.

ucharge

Charge with units

uepsilon

Lennard-Jones epsilon parameter with units

uepsilon_14

The 1-4 Lennard-Jones epsilon parameter

umass
urmin

Lennard-Jones Rmin/2 parameter with units

urmin_14

The 1-4 Lennard-Jones Rmin/2 parameter with units

usigma

Lennard-Jones sigma parameter with units

usigma_14

The 1-4 Lennard-Jones sigma parameter with units

usolvent_radius

Solvation radius with units attached

Methods

angle_to(other)

Log this atom as angled to another atom.

bond_to(other)

Log this atom as bonded to another atom.

dihedral_to(other)

Log this atom as dihedral-ed to another atom.

exclude(other)

Add one atom to my arbitrary exclusion list

nonbonded_exclusions([only_greater, index_from])

Returns the total number of nonbonded atom exclusions for this atom.

prune_exclusions()

For extra points, the exclusion partners may be filled before the bond, angle, dihedral, and tortor partners.

tortor_to(other)

Log this atom as 1-5 partners to another atom

property angle_partners

List of all angle partners that are NOT bond partners

angle_to(other)[source]

Log this atom as angled to another atom.

Parameters
otherAtom

An atom that will be added to angle_partners

Notes

This action adds self to other.angle_partners. Raises MoleculeError if other is self

property bond_partners

Go through all bonded partners

bond_to(other)[source]

Log this atom as bonded to another atom.

Parameters
otherAtom

An atom that will be added to bond_partners

Notes

This action adds self to other.bond_partners. Raises MoleculeError if other is self

property charge
property dihedral_partners

List of all dihedral partners that are NOT angle or bond partners

dihedral_to(other)[source]

Log this atom as dihedral-ed to another atom.

Parameters
otherAtom

An atom that will be added to dihedral_partners

Notes

This action adds self to other.dihedral_partners. Raises MoleculeError if other is self

property element
property element_name
property epsilon

Lennard-Jones epsilon parameter (the Lennard-Jones well depth)

property epsilon_14

The 1-4 Lennard-Jones epsilon parameter

exclude(other)[source]

Add one atom to my arbitrary exclusion list

Parameters
otherAtom

An atom that will be added to exclusion_partners.

Notes

This action adds self to other.exclusion_partners. Raises MoleculeError if other is self

property exclusion_partners

List of all exclusions not otherwise excluded by bonds/angles/torsions

nonbonded_exclusions(only_greater=True, index_from=0)[source]

Returns the total number of nonbonded atom exclusions for this atom. The general rules for building the exclusion list is to include both exceptions AND exclusions (i.e., the Amber scaling of 1-4 interactions means that the 1-4 terms are excluded and a special pairlist is built to handle those exceptions).

All atoms in the _partners arrays are nonbonded exclusions.

Parameters
only_greaterbool, optional

If True (default), only atoms whose idx value is greater than this Atom`s `idx will be counted as an exclusion (to avoid double- counting exclusions). If False, all exclusions will be counted.

index_fromint, optional

This is the index of the first atom, and is intended to be 0 (for C- and Python-style numbering, default) or 1 (for Fortran-style numbering, such as that used in the Amber and CHARMM topology files)

Returns
list of int

The returned list will be the atom indexes of the exclusion partners for this atom (indexing starts from index_from)

Notes

If this instance’s idx attribute evaluates to -1 – meaning it is not in an AtomList – a IndexError will be raised. If you have two extra points (i.e., those with atomic numbers of 0) bonded to each other, this routine may raise a RuntimeError if the recursion depth is exceeded.

prune_exclusions()[source]

For extra points, the exclusion partners may be filled before the bond, angle, dihedral, and tortor partners. Since we don’t want memory of these exclusions if any of those topological features were to break, we want to remove those from the exclusion list. This function makes sure that nothing in the bond, angle, dihedral, and tortor lists appears in the exclusion list.

property rmin

Lennard-Jones Rmin/2 parameter (the Lennard-Jones radius)

property rmin_14

The 1-4 Lennard-Jones Rmin/2 parameter

property sigma

Lennard-Jones sigma parameter – directly related to Rmin

property sigma_14

Lennard-Jones sigma parameter – directly related to Rmin

property tortor_partners

List of all 1-5 partners that are NOT in angle or bond partners. This is only used in the Amoeba force field

tortor_to(other)[source]

Log this atom as 1-5 partners to another atom

Parameters
otherAtom

An atom that will be added to tortor_partners

Notes

This action adds self to other.tortor_partners. Raises MoleculeError if other is self

property ucharge

Charge with units

property uepsilon

Lennard-Jones epsilon parameter with units

property uepsilon_14

The 1-4 Lennard-Jones epsilon parameter

property umass
property urmin

Lennard-Jones Rmin/2 parameter with units

property urmin_14

The 1-4 Lennard-Jones Rmin/2 parameter with units

property usigma

Lennard-Jones sigma parameter with units

property usigma_14

The 1-4 Lennard-Jones sigma parameter with units

property usolvent_radius

Solvation radius with units attached

class parmed.AtomList(*args)[source]

Bases: parmed.topologyobjects.TrackedList

Array of Atoms

Notes

Deleting an atom from the AtomList also deletes that atom from the residue it belongs to.

Methods

assign_nbidx_from_types()

Assigns the nb_idx attribute of every atom inside here from the atom_type definition.

claim()

This method causes this list to “claim” all of the items it contains and subsequently indexes all of its items.

clear(/)

Remove all items from list.

copy(/)

Return a shallow copy of the list.

count(value, /)

Return number of occurrences of value.

find_original_index(idx)

Finds an atom with the given original index.

index(value[, start, stop])

Return first index of value.

index_members()

Assigns the idx variable for every member of this list to its place in the list, if the members of this list permit

prune_unused()

This method inspects the used attribute of all of its members, if it has one, and deletes any item in which it is set to False

reverse(/)

Reverse IN PLACE.

sort(*[, key, reverse])

Stable sort IN PLACE.

unmark()

Unmark all atoms in this list

append

extend

insert

pop

remove

append(*args, **kwargs)[source]

Append object to the end of the list.

assign_nbidx_from_types()[source]

Assigns the nb_idx attribute of every atom inside here from the atom_type definition. If the atom_type is not assigned, RuntimeError is raised.

Returns
list of dict

Each element is a set of the nb_idx indices for which NBFIX alterations are defined for the type with that given that index (minus 1, to adjust for indexing from 0 and nb_idx starting from 1)

extend(*args, **kwargs)[source]

Extend list by appending elements from the iterable.

find_original_index(idx)[source]

Finds an atom with the given original index. Cannot assume that the original indexes are in order, since reordering may have been necessary. As a result, the complexity of this algorithm is O(N)

Parameters
idxint

The integer corresponding to the original index

Returns
atomAtom

The atom with the original index idx

Raises
IndexError

If no atom has the original index idx

insert(*args, **kwargs)[source]

Insert object before index.

pop(*args, **kwargs)[source]

Remove and return item at index (default last).

Raises IndexError if list is empty or index is out of range.

remove(*args, **kwargs)[source]

Remove first occurrence of value.

Raises ValueError if the value is not present.

unmark()[source]

Unmark all atoms in this list

class parmed.AtomType(name, number, mass, atomic_number=- 1, bond_type=None, charge=0.0)[source]

Bases: object

Atom types can either be compared by indexes or names. Can be assigned with a string, integer, (string is not automatically cast) or with both.

Parameters
namestr

The name of the atom type

numberint

The serial index of the atom type

massfloat

The mass of the atom type

atomic_numberint, optional

The atomic number of the element of the atom type. Default -1

bond_typestr, optional

If defined, this is the type name used to look up bonded parameters. Default is None (which falls back to name)

chargefloat, optional

If defined, this is the partial atomic charge in elementary charge units. Default is None

Notes

This object is primarily used to build parameter databases from parameter files. Also, sigma is related to Rmin, but rmin is Rmin/2, so there is an extra factor of 2 in the sigma for this reason.

Examples

>>> at = AtomType('HA', 1, 1.008, 1)
>>> at.name, at.number
('HA', 1)
>>> at2 = AtomType('CA', 2, 12.01, 6)
>>> at2.name, at2.number
('CA', 2)
>>> print("%s: %d" % (str(at), int(at)))
HA: 1
Attributes
bond_type
sigma

Sigma is Rmin / 2^(1/6)

sigma_14

Sigma is Rmin / 2^(1/6)

uepsilon
uepsilon_14
urmin
urmin_14
usigma
usigma_14

Methods

add_nbfix(typename, rmin, epsilon[, rmin14, …])

Adds a new NBFIX exclusion for this atom type

set_lj_params(eps, rmin[, eps14, rmin14])

Sets Lennard-Jones parameters on this atom type

add_nbfix(typename, rmin, epsilon, rmin14=None, epsilon14=None)[source]

Adds a new NBFIX exclusion for this atom type

Parameters
typenamestr

The name of the other type with which this NBFIX is defined

rminfloat

The combined Rmin value for this NBFIXed pair. If no units, assumed to be in Angstroms

epsilonfloat

The combined epsilon value for this NBFIXed pair. If no units, assumed to be in kcal/mol

rmin14float, optional

Same as rmin, but for 1-4 interactions. If None (default), it is given the same value as rmin

epsilon14float, optional

Same as epsilon, but for 1-4 interactions. If None (default), it is given the same value as rmin

property bond_type
set_lj_params(eps, rmin, eps14=None, rmin14=None)[source]

Sets Lennard-Jones parameters on this atom type

property sigma

Sigma is Rmin / 2^(1/6)

property sigma_14

Sigma is Rmin / 2^(1/6)

property uepsilon
property uepsilon_14
property urmin
property urmin_14
property usigma
property usigma_14
class parmed.Bond(atom1, atom2, type=None, order=1.0)[source]

Bases: object

A covalent bond connecting two atoms.

Parameters
atom1Atom

The first atom involved in the bond

atom2Atom

The other atom involved in the bond

typeBondType or None, optional

The bond type that defines the parameters for this bond. Default is None

orderfloat, optional
The bond order of this bond. Bonds are classified as follows:

1.0 – single bond 2.0 – double bond 3.0 – triple bond 1.5 – aromatic bond 1.25 – amide bond

Default is 1.0

Notes

You can test whether an Atom is contained within the bond using the in operator. A MoleculeError is raised if atom1 and atom2 are identical. This bond instance is append`ed to the `bonds list for both atom1 and atom2 and is automatically removed from those lists upon garbage collection

Examples

>>> a1, a2 = Atom(), Atom()
>>> bond = Bond(a1, a2)
>>> a1 in bond and a2 in bond
True
Attributes
order

Bond order.

Methods

delete()

Deletes this bond from the atoms that make it up.

energy()

Measures the current bond energy

measure()

Measures the current bond

uenergy()

Same as energy(), but with units

umeasure()

Same as “measure”, but with units

delete()[source]

Deletes this bond from the atoms that make it up. This method removes this bond from the bonds list for both atom1 and atom2 as well as removing atom1 from atom2.bond_partners and vice versa.

energy()[source]

Measures the current bond energy

Returns
energyfloat or None

Bond strain energy in kcal/mol. Return value is None if either the coordinates of either atom is not set or the bond type is not set

measure()[source]

Measures the current bond

Returns
measurementfloat or None

If the atoms have coordinates, returns the distance between the two atoms. If any coordinates are missing, returns None

property order

Bond order. See description in Bond argument list

uenergy()[source]

Same as energy(), but with units

umeasure()[source]

Same as “measure”, but with units

class parmed.BondType(k, req, list=None)[source]

Bases: parmed.topologyobjects._ParameterType, parmed.topologyobjects._ListItem

A bond type with a set of bond parameters

Parameters
kfloat

Force constant in kcal/mol/Angstrom^2

reqfloat

Equilibrium bond distance in Angstroms

listTrackedList

A list of :class:`BondType`s in which this is a member

Notes

Two BondType`s are equal if their `k and req attributes are equal

Examples

>>> bt1 = BondType(10.0, 1.0)
>>> bt2 = BondType(10.0, 1.0)
>>> bt1 is bt2
False
>>> bt1 == bt2
True
>>> bt1.idx # Not in a list or container
-1

As part of a list, they can be indexed

>>> bond_list = []
>>> bond_list.append(BondType(10.0, 1.0, list=bond_list))
>>> bond_list.append(BondType(10.0, 1.0, list=bond_list))
>>> bond_list[0].idx
0
>>> bond_list[1].idx
1
Attributes
idx
uk
ureq
property uk
property ureq
class parmed.ChiralFrame(atom1, atom2, chirality)[source]

Bases: object

A chiral frame as defined in the AMOEBA force field. It defines the frame of reference for a chiral center

Parameters
atom1Atom

The first atom defined in the chiral frame

atom2Atom

The second atom defined in the chiral frame

chiralityint

Either 1 or -1 to identify directionality. A ValueError is raised if a different value is provided

Notes

A chiral frame can only contain atoms.

class parmed.Cmap(atom1, atom2, atom3, atom4, atom5, type=None)[source]

Bases: object

A coupled-torsion correction map term defined between 5 atoms connected by four covalent bonds. This is a coupled-torsion potential in which the torsions are consecutive.

Parameters
atom1Atom

An atom on one end of the valence coupled-torsion bonded to atom2

atom2Atom

An atom in the middle of the CMAP bonded to atoms 1 and 3

atom3Atom

An atom in the middle of the CMAP bonded to atoms 2 and 4

atom4Atom

An atom in the middle of the CMAP bonded to atoms 3 and 5

atom5Atom

An atom in the middle of the CMAP bonded to atom 4

typeCmapType

The CmapType object containing the parameter map for this term

Notes

A CMAP can contain bonds or atoms. A bond is contained if it exists between atoms 1 and 2, between atoms 2 and 3, between atoms 3 and 4, or between atoms 4 and 5.

Examples

>>> a1, a2, a3, a4, a5 = Atom(), Atom(), Atom(), Atom(), Atom()
>>> cmap = Cmap(a1, a2, a3, a4, a5)
>>> Bond(a1, a2) in cmap and Bond(a2, a3) in cmap
True
>>> Bond(a1, a3) in cmap
False

Methods

delete()

Deletes this Cmap from the atoms that make it up.

extended(atom1, atom2, atom3, atom4, atom5, …)

Alternative constructor for correction maps defined with 8 atoms (each torsion being separately specified).

same_atoms(thing)

A coupled-torsion is equivalent if the 5 atoms are in the same or reverse order Allow comparison with another type of cmap or with a sequence of 5 indexes

delete()[source]

Deletes this Cmap from the atoms that make it up. This method removes the Cmap from the cmaps list for atom1, atom2, atom3, atom4, and atom5

classmethod extended(atom1, atom2, atom3, atom4, atom5, atom6, atom7, atom8, type=None)[source]

Alternative constructor for correction maps defined with 8 atoms (each torsion being separately specified). Correction maps are, to the best of my knowledge, used to parametrized “consecutive” torsions (i.e., atoms 2, 3, and 4 are common to both torsions). However, the CHARMM definition specifies the torsions separately.

If the torsions are _not_ consecutive (i.e., if atoms 2, 3, and 4 are not the same as atoms 5, 6, and 7), NotImplementedError is raised.

Parameters
atom1Atom

The first atom of the first torsion

atom2Atom

The second atom of the first torsion

atom3Atom

The third atom of the first torsion

atom4Atom

The fourth atom of the first torsion

atom5Atom

The first atom of the second torsion

atom6Atom

The second atom of the second torsion

atom7Atom

The third atom of the second torsion

atom8Atom

The fourth atom of the second torsion

typeCmapType

The CmapType object containing the parameter map for this term

same_atoms(thing)[source]

A coupled-torsion is equivalent if the 5 atoms are in the same or reverse order Allow comparison with another type of cmap or with a sequence of 5 indexes

class parmed.CmapType(resolution, grid, comments=None, list=None)[source]

Bases: parmed.topologyobjects._ParameterType, parmed.topologyobjects._ListItem

A CMAP type with a potential energy interpoloation grid mapping out the 2-D potential of coupled torsions.

Parameters
resolutionint

The number of grid points in the correction map potential in both torsion dimensions

griditerable of floats

This must be a 1-dimensional list of grid values. The dimension must be resolution*resolution, and must be row-major (i.e., the second dimension changes the fastest) when indexed with 2 indices.

listTrackedList

A list of `CmapType`s in which this is a member

Raises
TypeError is raised if the grid does not have the correct number of elements
for the given resolution

Notes

Two `CmapType`s are equal if their resolution is the same and each grid point is the same to within 1E-8

See the docs for _CmapGrid for information on how to access the interpolating grid data if necessary.

Examples

>>> ct = CmapType(2, [0, 1, 2, 3])
>>> ct.grid[0,0], ct.grid[0]
(0, 0)
>>> ct.grid[0,1], ct.grid[1]
(1, 1)
>>> ct.grid[1,0], ct.grid[2]
(2, 2)
>>> ct.grid[1,1], ct.grid[3]
(3, 3)
>>> ct.idx # not part of a list or iterable
-1
Attributes
resolutionint

Potential grid resolution (see description in Parameters)

grid_CmapGrid

A _CmapGrid object defining the interpolating potential energy grid, with each point having the units kcal/mol

commentslist(str)

List of strings that represent comments about this parameter type

listTrackedList

If not None, this is a list in which this instance _may_ be a member

idxint

The index of this CmapType inside its containing list

class parmed.Dihedral(atom1, atom2, atom3, atom4, improper=False, ignore_end=False, type=None)[source]

Bases: parmed.topologyobjects._FourAtomTerm

A valence dihedral between 4 atoms separated by three covalent bonds.

Parameters
atom1Atom

An atom on one end of the valence dihedral bonded to atom 2

atom2Atom

An atom in the middle of the valence dihedral bonded to atom1 and atom3

atom3Atom

An atom in the middle of the valence dihedral bonded to atom2 and atom4

atom4Atom

An atom on the other end of the valence dihedral bonded to atom 3

improperbool

If True, this is an Amber-style improper torsion, where atom3 is the “central” atom bonded to atoms 1, 2, and 4 (atoms 1, 2, and 4 are only bonded to atom 3 in this instance)

ignore_endbool

If True, the end-group interactions for this torsion are ignored, either because it is involved in a ring where the end-group atoms are excluded or because it is one term in a multi-term dihedral expansion

typeDihedralType

The DihedralType object containing the parameters for this dihedral

Notes

A Dihedral can contain bonds or atoms. A bond is contained if it exists between atoms 1 and 2, between atoms 2 and 3, or between atoms 3 and 4.

Examples

>>> a1, a2, a3, a4 = Atom(), Atom(), Atom(), Atom()
>>> dihed = Dihedral(a1, a2, a3, a4)
>>> Bond(a1,a2) in dihed and Bond(a3,a2) in dihed and Bond(a3,a4) in dihed
True
>>> a1 in dihed and a2 in dihed and a3 in dihed and a4 in dihed
True
>>> Bond(a1, a4) in dihed # this is not part of the angle definition
False

For improper torsions, the bond pattern is different

>>> a1, a2, a3, a4 = Atom(), Atom(), Atom(), Atom()
>>> dihed = Dihedral(a1, a2, a3, a4, improper=True)
>>> Bond(a1,a3) in dihed and Bond(a2,a3) in dihed and Bond(a3,a4) in dihed
True
>>> Bond(a1,a2) in dihed # Not like a normal dihedral!
False
>>> a1 in dihed and a2 in dihed and a3 in dihed and a4 in dihed
True
Attributes
funct

Methods

delete()

Deletes this dihedral from the atoms that make it up.

energy()

Measures the current dihedral angle energy

measure()

Measures the current dihedral angle

same_atoms(thing)

Determines if this dihedral has the same atoms (or atom indexes) as another object

uenergy()

Same as energy(), but with units

umeasure()

Same as measure(), but with units

delete()[source]

Deletes this dihedral from the atoms that make it up. This method removes this dihedral from the dihedrals list for atom1, atom2, atom3, and atom4 as well as removing each atom form each others’ dihedral partner lists

energy()[source]

Measures the current dihedral angle energy

Returns
energyfloat or None

Dihedral angle energy in kcal/mol. Return value is None if either the coordinates of either atom is not set or the angle type is not set

property funct
measure()[source]

Measures the current dihedral angle

Returns
measurementfloat or None

If the atoms have coordinates, returns the dihedral angle between the four atoms in degrees. If any coordinates are missing, returns None

same_atoms(thing)[source]

Determines if this dihedral has the same atoms (or atom indexes) as another object

Parameters
thingDihedral or other iterable

A Dihedral or an iterable with 4 indexes that will be used to identify whether this torsion has the same atoms

Returns
is_samebool

True if thing and self have the same atoms. False otherwise

Raises
TypeError

Notes

This raises a TypeError if thing is not a Dihedral and is not iterable

uenergy()[source]

Same as energy(), but with units

umeasure()[source]

Same as measure(), but with units

class parmed.DihedralType(phi_k, per, phase, scee=1.0, scnb=1.0, list=None)[source]

Bases: parmed.topologyobjects._ParameterType, parmed.topologyobjects._ListItem

A dihedral type with a set of dihedral parameters

Parameters
phi_kfloat

The force constant in kcal/mol

perint

The dihedral periodicity

phasefloat

The dihedral phase in degrees

sceefloat, optional

1-4 electrostatic scaling factor. Default is 1.0

scnbfloat, optional

1-4 Lennard-Jones scaling factor. Default is 1.0

listTrackedList

A list of `DihedralType`s in which this is a member

Notes

Two DihedralType`s are equal if their `phi_k, per, phase, scee, and scnb attributes are equal

Examples

>>> dt1 = DihedralType(10.0, 2, 180.0, 1.2, 2.0)
>>> dt2 = DihedralType(10.0, 2, 180.0, 1.2, 2.0)
>>> dt1 is dt2
False
>>> dt1 == dt2
True
>>> dt1.idx # not part of any list or iterable
-1

As part of a list, they can be indexed

>>> dihedral_list = []
>>> dihedral_list.append(DihedralType(10.0, 2, 180.0, 1.2, 2.0,
...                                   list=dihedral_list))
>>> dihedral_list.append(DihedralType(10.0, 2, 180.0, 1.2, 2.0,
...                                   list=dihedral_list))
>>> dihedral_list[0].idx
0
>>> dihedral_list[1].idx
1
Attributes
idx
uphase
uphi_k
property uphase
property uphi_k
class parmed.DihedralTypeList(*args, **kwargs)[source]

Bases: list, parmed.topologyobjects._ListItem

Dihedral types are a Fourier expansion of terms. In some cases, they are stored in a list like this one inside another TrackedList. In other cases, each term is a separate entry in the TrackedList.

In cases where DihedralType`s are stored with every term in the same container, this object supports list assignment and indexing like :class:`DihedralType.

Parameters
*argsobjects

Any arguments that list would take.

listTrackedList, optional

A list that “contains” this DihedralTypeList instance. This is a keyword-only argument. Default is None (i.e., belonging to no list)

**kwargskeyword argument list

All other keyword arguments passed directly to the list constructor

Attributes
idx
penalty

Methods

append(other[, override])

Adds a DihedralType to the DihedralTypeList

clear(/)

Remove all items from list.

copy(/)

Return a shallow copy of the list.

count(value, /)

Return number of occurrences of value.

extend(iterable, /)

Extend list by appending elements from the iterable.

from_rbtorsion(rbtorsion)

Creates a Fourier series of proper torsions from a Ryckaerts-Bellemans torsion.

index(value[, start, stop])

Return first index of value.

insert(index, object, /)

Insert object before index.

pop([index])

Remove and return item at index (default last).

remove(value, /)

Remove first occurrence of value.

reverse(/)

Reverse IN PLACE.

sort(*[, key, reverse])

Stable sort IN PLACE.

append(other, override=False)[source]

Adds a DihedralType to the DihedralTypeList

Parameters
otherDihedralType

The DihedralType instance to add to this list. It cannot have the same periodicity as any other DihedralType in the list

overridebool, optional, default=False

If True, this will override an existing torsion, but will raise a warning if it does so

Raises
TypeError if other is not an instance of DihedralTypeList
ParameterError if other has the same periodicity as another member in
this list and override is False
ParameterWarning if other has same periodicit as another member in this
list and override is True
classmethod from_rbtorsion(rbtorsion)[source]

Creates a Fourier series of proper torsions from a Ryckaerts-Bellemans torsion.

Parameters
rbtorsionRBTorsionType

The R-B torsion type to convert to a series of proper torsions

Raises
ValueError if all terms in rbtorsion are zero except for c0
property penalty
class parmed.DrudeAnisotropy(atom1, atom2, atom3, atom4, a11, a22)[source]

Bases: parmed.topologyobjects._FourAtomTerm

A description of an anisotropically polarizable atom.

Atom 1 is a DrudeAtom whose polarizability is anisotropic. The other three atoms define the coordinate frame.

Parameters
atom1DrudeAtom

the polarizable atom

atom2Atom

the second atom defining the coordinate frame

atom3Atom

the third atom defining the coordinate frame

atom4Atom

the fourth atom defining the coordinate frame

a11float

the scale factor for the polarizability along the direction defined by atom1 and atom2

a22float

the scale factor for the polarizability along the direction defined by atom3 and atom4

Methods

delete()

Sets all atoms in this term to None, and its type if it exists

class parmed.DrudeAtom(alpha=0.0, thole=1.3, drude_type='DRUD', **kwargs)[source]

Bases: parmed.topologyobjects.Atom

An Atom that has a Drude particle attached to it. This is a subclass of Atom, so it also has all the properties defined for regular Atoms.

Parameters
alphafloat

the atomic polarizability

tholefloat

the Thole damping facior

drude_typestr

the atom type to use for the Drude particle.

Attributes
angle_partners

List of all angle partners that are NOT bond partners

bond_partners

Go through all bonded partners

charge
dihedral_partners

List of all dihedral partners that are NOT angle or bond partners

drude_charge
element
element_name
epsilon

Lennard-Jones epsilon parameter (the Lennard-Jones well depth)

epsilon_14

The 1-4 Lennard-Jones epsilon parameter

exclusion_partners

List of all exclusions not otherwise excluded by bonds/angles/torsions

idx
rmin

Lennard-Jones Rmin/2 parameter (the Lennard-Jones radius)

rmin_14

The 1-4 Lennard-Jones Rmin/2 parameter

sigma

Lennard-Jones sigma parameter – directly related to Rmin

sigma_14

Lennard-Jones sigma parameter – directly related to Rmin

tortor_partners

List of all 1-5 partners that are NOT in angle or bond partners.

ucharge

Charge with units

uepsilon

Lennard-Jones epsilon parameter with units

uepsilon_14

The 1-4 Lennard-Jones epsilon parameter

umass
urmin

Lennard-Jones Rmin/2 parameter with units

urmin_14

The 1-4 Lennard-Jones Rmin/2 parameter with units

usigma

Lennard-Jones sigma parameter with units

usigma_14

The 1-4 Lennard-Jones sigma parameter with units

usolvent_radius

Solvation radius with units attached

Methods

angle_to(other)

Log this atom as angled to another atom.

bond_to(other)

Log this atom as bonded to another atom.

dihedral_to(other)

Log this atom as dihedral-ed to another atom.

exclude(other)

Add one atom to my arbitrary exclusion list

nonbonded_exclusions([only_greater, index_from])

Returns the total number of nonbonded atom exclusions for this atom.

prune_exclusions()

For extra points, the exclusion partners may be filled before the bond, angle, dihedral, and tortor partners.

tortor_to(other)

Log this atom as 1-5 partners to another atom

property drude_charge
class parmed.ExtraPoint(*args, **kwargs)[source]

Bases: parmed.topologyobjects.Atom

An extra point is a massless, virtual site that is used to extend the flexibility of partial atomic charge fitting. They can be used to, for instance, give atoms permanent multipoles in a fixed-charge format.

However, virtual sites are massless particles whose position is fixed with respect to a particular frame of reference. As a result, they must be treated specially when running dynamics. This class extends the Atom class with extra functionality specific to these “Extra points” or “virtual sites”. See the documentation for the Atom class for more information.

Parameters
weightslist of float, optional, keyword-only

This is the list of weights defining its frame.

See also

Atom

The ExtraPoint constructor also takes all Atom arguments as well, and shares all of the same properties

Attributes
angle_partners

List of all atoms angled to this atom and its parent

bond_partners

List of all atoms bonded to this atom and its parent

charge
dihedral_partners

List of all dihedral partners that are NOT angle or bond partners

element
element_name
epsilon

Lennard-Jones epsilon parameter (the Lennard-Jones well depth)

epsilon_14

The 1-4 Lennard-Jones epsilon parameter

exclusion_partners

List of all exclusions not otherwise excluded by bonds/angles/torsions

frame_type

The type of frame used for this extra point.

idx
parent

The parent atom of this extra point is the atom it is bonded to (there

radii
rmin

Lennard-Jones Rmin/2 parameter (the Lennard-Jones radius)

rmin_14

The 1-4 Lennard-Jones Rmin/2 parameter

sigma

Lennard-Jones sigma parameter – directly related to Rmin

sigma_14

Lennard-Jones sigma parameter – directly related to Rmin

tortor_partners

List of all 1-5 partners that are NOT in angle or bond partners.

ucharge

Charge with units

uepsilon

Lennard-Jones epsilon parameter with units

uepsilon_14

The 1-4 Lennard-Jones epsilon parameter

umass
urmin

Lennard-Jones Rmin/2 parameter with units

urmin_14

The 1-4 Lennard-Jones Rmin/2 parameter with units

usigma

Lennard-Jones sigma parameter with units

usigma_14

The 1-4 Lennard-Jones sigma parameter with units

usolvent_radius

Solvation radius with units attached

Methods

angle_to(other)

Log this atom as angled to another atom.

bond_to(other)

Log this atom as bonded to another atom.

dihedral_to(other)

Log this atom as dihedral-ed to another atom.

exclude(other)

Add one atom to my arbitrary exclusion list

nonbonded_exclusions([only_greater, index_from])

Returns the total number of nonbonded atom exclusions for this atom.

prune_exclusions()

For extra points, the exclusion partners may be filled before the bond, angle, dihedral, and tortor partners.

tortor_to(other)

Log this atom as 1-5 partners to another atom

property angle_partners

List of all atoms angled to this atom and its parent

property bond_partners

List of all atoms bonded to this atom and its parent

property dihedral_partners

List of all dihedral partners that are NOT angle or bond partners

property exclusion_partners

List of all exclusions not otherwise excluded by bonds/angles/torsions

property frame_type

The type of frame used for this extra point. Whether the EP lies beyond (outside) or between (inside) the is determined from the positions of each atom in the frame and the extra point. If those are not set, the EP is assumed to be between the two points in the frame

property parent

The parent atom of this extra point is the atom it is bonded to (there should only be 1 bond partner, but the parent is defined as its first)

property radii
property tortor_partners

List of all 1-5 partners that are NOT in angle or bond partners. This is only used in the Amoeba force field

class parmed.Group(atom, type, move)[source]

Bases: object

An ‘interacting’ group defined by CHARMM PSF files

Parameters
atomAtom

The first atom within a group

typeint

Flag for group information; 0 when all atoms have zero charge, 1 when group has a net zero charge but at least one atom has a non-zero partial charge, 2 when the net charge of the group is not zero

moveint

0 if the atoms are not fixed, 1 when they are

Notes

See the discussion on Github for the source of the meanings of these variables: https://github.com/ParmEd/ParmEd/pull/307#issuecomment-128244134

class parmed.Improper(atom1, atom2, atom3, atom4, type=None)[source]

Bases: parmed.topologyobjects._FourAtomTerm

A CHARMM-style improper torsion between 4 atoms. The first atom must be the central atom, as shown in the schematic below

A3 | |

A4 —– A1 —– A2

Parameters
atom1Atom

The central atom A1 in the schematic above

atom2Atom

An atom in the improper, A2 in the schematic above

atom3Atom

An atom in the improper, A3 in the schematic above

atom4Atom

An atom in the improper, A4 in the schematic above

typeImproperType

The ImproperType object containing the parameters for this improper torsion

Notes

An Improper torsion can contain bonds or atoms. A bond is contained if it exists between atom 1 and any other atom. Raises MoleculeError if any of the atoms are duplicates.

Examples

>>> a1, a2, a3, a4 = Atom(), Atom(), Atom(), Atom()
>>> imp = Improper(a1, a2, a3, a4)
>>> Bond(a1, a2) in imp and Bond(a1, a3) in imp and Bond(a1, a4) in imp
True
>>> Bond(a2, a3) in imp
False

Methods

delete()

Deletes this Improper from the atoms that make it up.

energy()

Measures the current dihedral angle energy

measure()

Measures the current torsional angle

same_atoms(thing)

An improper has the same 4 atoms if atom1 is the same for both and the other 3 can be in any order

uenergy()

Same as energy(), but with units

umeasure()

Same as measure(), but with units

delete()[source]

Deletes this Improper from the atoms that make it up. This method removes this Improper from the impropers list for atom1, atom2, atom3, and atom4

energy()[source]

Measures the current dihedral angle energy

Returns
energyfloat or None

Dihedral angle energy in kcal/mol. Return value is None if either the coordinates of either atom is not set or the angle type is not set

measure()[source]

Measures the current torsional angle

Returns
measurementfloat or None

If the atoms have coordinates, returns the torsional angle between the four atoms. If any coordinates are missing, returns None

same_atoms(thing)[source]

An improper has the same 4 atoms if atom1 is the same for both and the other 3 can be in any order

uenergy()[source]

Same as energy(), but with units

umeasure()[source]

Same as measure(), but with units

class parmed.ImproperType(psi_k, psi_eq, list=None)[source]

Bases: parmed.topologyobjects._ParameterType, parmed.topologyobjects._ListItem

An improper type with a set of improper torsion parameters

Parameters
psi_kfloat

Force constant in kcal/mol/radians^2

psi_eqfloat

Equilibrium torsion angle in Degrees

listTrackedList

A list of `ImproperType`s in which this is a member

Notes

Two ImproperType`s are equal if their `psi_k and psi_eq attributes are equal

Examples

>>> it1 = ImproperType(10.0, 180.0)
>>> it2 = ImproperType(10.0, 180.0)
>>> it1 is it2
False
>>> it1 == it2
True
>>> it1.idx # Not part of any list or iterable
-1

As part of a list, they can be indexed

>>> improper_list = []
>>> improper_list.append(ImproperType(10.0, 180.0, list=improper_list))
>>> improper_list.append(ImproperType(10.0, 180.0, list=improper_list))
>>> improper_list[0].idx
0
>>> improper_list[1].idx
1
Attributes
idx
upsi_eq
upsi_k
property upsi_eq
property upsi_k

Bases: object

An intra-residue “Link” as defined by the PDB standard:

See http://www.wwpdb.org/documentation/file-format-content/format33/sect6.html#LINK for more information

Parameters
atom1Atom

The first Atom involved in the Link

atom2Atom

The other atom to which atom1 is bonded in this link

lengthfloat

The length of the link

symmetry_op1str, optional

The first symmetry operator for the link

symmetry_op2str, optional

The second symmetry operator for the link

class parmed.MultipoleFrame(atom, frame_pt_num, vectail, vechead, nvec)[source]

Bases: object

This defines the frame of reference for computing multipole interactions in the AMOEBA force field.

Parameters
atomAtom

The atom for which the frame of reference is defined

frame_pt_numint

The frame point number

vectailint

The vector tail index

vecheadint

The vector head index

nvecint

The number of vectors

Examples

>>> atom = Atom()
>>> mf = MultipoleFrame(atom, 0, 1, 2, 3)
>>> atom in mf
True
>>> mf.frame_pt_num
0
>>> mf.vectail
1
>>> mf.vechead
2
>>> mf.nvec
3
class parmed.NonbondedException(atom1, atom2, type=None)[source]

Bases: object

The AMOEBA force field has complex exclusion and exception rules (referred to as “adjustments” in the Amber-converted files). This class stores per-particle exceptions.

Parameters
atom1Atom

One of the atoms in the exclusion pair

atom2Atom

The other atom in the exclusion pair

typeNonbondedExceptionType

The nonbonded exception type that describes how the various nonbonded interactions between these two atoms should work

Notes

NonbondedException objects “contain” two atoms and will return True when used with the binary in operator

class parmed.NonbondedExceptionType(rmin, epsilon, chgscale=1.0, list=None)[source]

Bases: parmed.topologyobjects._ParameterType, parmed.topologyobjects._ListItem

A parameter describing how the various nonbonded interactions between a particular pair of atoms behaves in a specified nonbonded exception (e.g., in 1-4 interacting terms)

Parameters
rminfloat

The combined Rmin value for this particular pair of atom types (dimension length, default units are Angstroms)

epsilonfloat

The combined well-depth value for this particular pair of atom types (dimension energy, default units are kcal/mol)

chgscalefloat, optional

The scaling factor by which to multiply the product of the charges for this pair. Default is 1.0.

listTrackedList

The list containing this nonbonded exception

Attributes
idx
sigma
uepsilon
urmin
usigma
property sigma
property uepsilon
property urmin
property usigma
class parmed.OutOfPlaneBend(atom1, atom2, atom3, atom4, type=None)[source]

Bases: parmed.topologyobjects._FourAtomTerm

Out-of-plane bending term in the AMOEBA force field. The bond pattern is the same as TrigonalAngle

Parameters
atom1Atom

The first atom involved in the trigonal angle

atom2Atom

The central atom involved in the trigonal angle

atom3Atom

The third atom involved in the trigonal angle

atom4Atom

The fourth atom involved in the trigonal angle

typeOutOfPlaneBendType

The angle type containing the parameters

Notes

Either `Atom`s or `Bond`s can be contained within this trigonal angle

Methods

delete()

Sets all atoms in this term to None, and its type if it exists

class parmed.OutOfPlaneBendType(k, list=None)[source]

Bases: parmed.topologyobjects._ParameterType, parmed.topologyobjects._ListItem

An angle type with a set of angle parameters

Parameters
kfloat

Force constant in kcal/mol/radians^2

listTrackedList

A list of `OutOfPlaneBendType`s in which this is a member

Notes

Two OutOfPlaneBendType`s are equal if their `k attribute is equal

Examples

>>> ot1 = OutOfPlaneBendType(10.0)
>>> ot2 = OutOfPlaneBendType(10.0)
>>> ot1 is ot2
False
>>> ot1 == ot2
True
>>> ot1.idx # not part of any list or iterable
-1

As part of a list, they can be indexed

>>> oopbend_list = []
>>> oopbend_list.append(OutOfPlaneBendType(10.0, list=oopbend_list))
>>> oopbend_list.append(OutOfPlaneBendType(10.0, list=oopbend_list))
>>> oopbend_list[0].idx
0
>>> oopbend_list[1].idx
1
Attributes
idx
uk
property uk
class parmed.OutOfPlaneExtraPointFrame(ep, angle=54.735)[source]

Bases: object

This class defines a frame of reference for a given extra point with a frame of reference defined by 3 particles, but with the virtual site out of the plane of those 3 particles. For example, TIP5P

Parameters
epExtraPoint

The extra point defined by this frame

anglefloat

The angle out-of-plane that the extra point is. By default, it is half of the angle of a tetrahedron. Given in degrees

Methods

get_atoms()

Returns the particles involved in the frame

get_weights()

Returns the weights for the three particles

get_atoms()[source]

Returns the particles involved in the frame

Returns
a1, a2, a3Atom, Atom, Atom

a1 is the parent atom of the extra point. a2 and a3 are the other atoms that are both bonded to the parent atom but are not EPs

get_weights()[source]

Returns the weights for the three particles

Returns
w1, w2, w3float, float, float

w1 and w2 are the weights with respect to the second two particles in the frame (i.e., NOT the parent atom). w3 is the weight of the cross-product

class parmed.ParameterSet[source]

Bases: object

Stores a parameter set defining a force field

Attributes
atom_typesdict(str:AtomType)

Dictionary mapping the names of the atom types to the corresponding AtomType instances

atom_types_intdict(int:AtomType)

Dictionary mapping the serial indexes of the atom types to the corresponding AtomType instances

atom_types_tupledict((str,int):AtomType)

Dictionary mapping the (name,number) tuple of the atom types to the corresponding AtomType instances

bond_typesdict((str,str):AtomType)

Dictionary mapping the 2-element tuple of the names of the two atom types involved in the bond to the BondType instances

angle_typesdict((str,str,str):AngleType)

Dictionary mapping the 3-element tuple of the names of the three atom types involved in the angle to the AngleType instances

urey_bradley_typesdict((str,str,str):BondType)

Dictionary mapping the 3-element tuple of the names of the three atom types involved in the angle to the BondType instances of the Urey-Bradley terms

dihedral_typesdict((str,str,str,str):list(DihedralType))

Dictionary mapping the 4-element tuple of the names of the four atom types involved in the dihedral to the DihedralType instances. Since each torsion term can be a multiterm expansion, each item corresponding to a key in this dict is a list of `DihedralType`s for each term in the expansion

improper_typesdict((str,str,str,str):ImproperType)

Dictionary mapping the 4-element tuple of the names of the four atom types involved in the improper torsion to the ImproperType instances

improper_periodic_typesdict((str,str,str,str):DihedralType)

Dictionary mapping the 4-element tuple of the names of the four atom types involved in the improper torsion (modeled as a Fourier series) to the DihedralType instances. Note, the central atom should always be put in the third position of the key

rb_torsion_typesdict((str,str,str,str):RBTorsionType)

Dictionary mapping the 4-element tuple of the names of the four atom types involved in the Ryckaert-Bellemans torsion to the RBTorsionType instances

cmap_typesdict((str,str,str,str,str,str,str,str):CmapType)

Dictionary mapping the 5-element tuple of the names of the five atom types involved in the correction map to the CmapType instances

nbfix_typesdict((str,str):(float,float))

Dictionary mapping the 2-element tuple of the names of the two atom types whose LJ terms are modified to the tuple of the (epsilon,rmin) terms for that off-diagonal term

pair_typesdict((str,str):NonbondedExceptionType)

Dictionary mapping the 2-element tuple of atom type names for which explicit exclusion rules should be applied

parametersetslist(str)

List of parameter set names processed in the current ParameterSet

residuesdict(str:ResidueTemplate|ResidueTemplateContainer)

A library of ResidueTemplate objects mapped to the residue name defined in the force field library files

Methods

condense([do_dihedrals])

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

from_structure(struct[, …])

Extracts known parameters from a Structure instance

typeify_templates()

Assign atom types to atom names in templates

property combining_rule
condense(do_dihedrals=True)[source]

This function goes through each of the parameter type dicts and eliminates duplicate types. After calling this function, every unique bond, angle, dihedral, improper, or cmap type will pair with EVERY key in the type mapping dictionaries that points to the equivalent type

Parameters
do_dihedralsbool=True

Dihedrals can take the longest time to compress since testing their equality takes the longest (this is complicated by the existence of multi-term torsions). This flag will allow you to skip condensing the dihedral parameter types (for large parameter sets, this can cut the compression time in half)

Returns
self

The instance that is being condensed

Notes

The return value allows you to condense the types at construction time.

classmethod from_structure(struct, allow_unequal_duplicates=True)[source]

Extracts known parameters from a Structure instance

Parameters
structparmed.structure.Structure

The parametrized Structure instance from which to extract parameters into a ParameterSet

allow_unequal_duplicatesbool, optional

If True, if two or more unequal parameter types are defined by the same atom types, the last one encountered will be assigned. If False, an exception will be raised. Default is True

Returns
paramsParameterSet

The parameter set with all parameters defined in the Structure

Raises
parmed.exceptions.ParameterError if allow_unequal_duplicates is
False and 2+ unequal parameters are defined between the same atom types.
NotImplementedError if any AMOEBA potential terms are defined in the
input structure

Notes

The parameters here are copies of the ones in the Structure, so modifying the generated ParameterSet will have no effect on struct. Furthermore, the first occurrence of each parameter will be used. If future ones differ, they will be silently ignored, since this is expected behavior in some instances (like with Gromacs topologies in the ff99sb-ildn force field) unless allow_unequal_duplicates is set to False

Dihedrals are a little trickier. They can be multi-term, which can be represented either as a single entry in dihedrals with a type of DihedralTypeList or multiple entries in dihedrals with a DihedralType parameter type. In this case, the parameter is constructed from either the first DihedralTypeList found or the first DihedralType of each periodicity found if no matching DihedralTypeList is found.

typeify_templates()[source]

Assign atom types to atom names in templates

class parmed.PiTorsion(atom1, atom2, atom3, atom4, atom5, atom6, type=None)[source]

Bases: object

Defines a pi-torsion term in the AMOEBA force field. The Pi-torsion is defined around a sp2-hybridized pi-delocalized orbital (like an amide) by 6 atoms, as shown in the schematic below.

A2 A5-AAA
/
/

A3—A4

/

/

AAA-A1 A6

In the above schematic, A3 and A4 are sp2-hybridized, and atoms A2 and A6 are bonded only to A3 and A4, respectively. Atoms A1 and A5 are each bonded to 3 other atoms.

Parameters
atom1Atom

atom A1 in the schematic above

atom2Atom

atom A2 in the schematic above

atom3Atom

atom A3 in the schematic above

atom4Atom

atom A4 in the schematic above

atom5Atom

atom A5 in the schematic above

atom6Atom

atom A6 in the schematic above

typeDihedralType

The parameters for this Pi-torsion

Notes

Both :class:`Bond`s and :class:`Atom`s can be contained in a pi-torsion

class parmed.RBTorsionType(c0, c1, c2, c3, c4, c5, scee=1.0, scnb=1.0, list=None)[source]

Bases: parmed.topologyobjects._ParameterType, parmed.topologyobjects._ListItem

A Ryckaert-Bellemans type with a set of dihedral parameters

Notes

Two `RBTorsionType`s are equal if their coefficients are all equal

Examples

>>> dt1 = RBTorsionType(10.0, 20.0, 30.0, 40.0, 50.0, 60.0)
>>> dt2 = RBTorsionType(10.0, 20.0, 30.0, 40.0, 50.0, 60.0)
>>> dt1 is dt2
False
>>> dt1 == dt2
True
>>> dt1.idx # not part of any list or iterable
-1

As part of a list, they can be indexed

>>> rb_torsion_list = []
>>> rb_torsion_list.append(RBTorsionType(10.0, 20.0, 30.0, 40.0, 50.0, 60.0,
...                                      list=rb_torsion_list))
>>> rb_torsion_list.append(RBTorsionType(10.0, 20.0, 30.0, 40.0, 50.0, 60.0,
...                                      list=rb_torsion_list))
>>> rb_torsion_list[0].idx
0
>>> rb_torsion_list[1].idx
1
Attributes
idx
uc0
uc1
uc2
uc3
uc4
uc5
property uc0
property uc1
property uc2
property uc3
property uc4
property uc5
class parmed.Residue(name, number=- 1, chain='', insertion_code='', segid='', list=None)[source]

Bases: parmed.topologyobjects._ListItem

A single residue that is composed of a small number of atoms

Parameters
namestr

Name of the residue. Typical convention is to choose a name that is 4 characters or shorter

numberint, optional

Residue number assigned in the input structure. Default is -1

chainstr, optional

The 1-letter chain identifier for this residue. Default is empty string

insertion_codestr, optional

The insertion code (used in PDB files) for this residue. Default is empty string

segidstr, optional

The segment identifier, used by CHARMM in a way similar to chain. Dfault is empty string

listTrackedList

List of residues in which this residue is a member

Notes

  • Iterating over a residue will iterate over the atoms. It is exactly equivalent to iterating over the atoms attribute

  • Supports testing if an Atom instance is contained in this residue

  • len() returns the number of atoms in this residue

Attributes
namestr

The name of this residue

numberint

The number of this residue in the input structure

idxint

The index of this residue inside the container. If this residue has no container, or it is not present in the container, idx is -1

chainstr

The 1-letter chain identifier for this residue

insertion_codestr

The insertion code (used in PDB files) for this residue

terbool

If True, there is a TER card directly after this residue (i.e., a molecule or chain ends). By default, it is False

listTrackedList

The container that _may_ have this residue contained inside

atomslist of :class`Atom instances

This is the list of `Atom`s that make up this residue

Methods

add_atom(atom)

Adds an atom to this residue

delete_atom(atom)

If an atom is present in this residue, delete it from the list of atoms.

is_empty()

Determines if there are any atoms in this residue

sort()

Sorts the atoms in this list by atom index

add_atom(atom)[source]

Adds an atom to this residue

Parameters
atomAtom

The atom to add to this residue

Notes

This action assigns the residue attribute to atom

delete_atom(atom)[source]

If an atom is present in this residue, delete it from the list of atoms. No change if an atom is not present in this residue.

is_empty()[source]

Determines if there are any atoms in this residue

Returns
empty: bool

True if there are no atoms left. False otherwise.

sort()[source]

Sorts the atoms in this list by atom index

class parmed.ResidueList(*args)[source]

Bases: parmed.topologyobjects.TrackedList

Array of Residue instances

Methods

add_atom(atom, resname, resnum[, chain, …])

Adds a new atom to the ResidueList, adding a new residue to this list if it has a different name or number as the last residue

claim()

This method causes this list to “claim” all of the items it contains and subsequently indexes all of its items.

clear(/)

Remove all items from list.

copy(/)

Return a shallow copy of the list.

count(value, /)

Return number of occurrences of value.

index(value[, start, stop])

Return first index of value.

index_members()

Assigns the idx variable for every member of this list to its place in the list, if the members of this list permit

prune()

This function goes through the residue list and removes all empty residues from the list.

prune_unused()

This method inspects the used attribute of all of its members, if it has one, and deletes any item in which it is set to False

reverse(/)

Reverse IN PLACE.

sort(*[, key, reverse])

Stable sort IN PLACE.

append

extend

insert

pop

remove

add_atom(atom, resname, resnum, chain='', inscode='', segid='')[source]

Adds a new atom to the ResidueList, adding a new residue to this list if it has a different name or number as the last residue

Parameters
atomAtom

The atom to add to this residue list

resnamestr

The name of the residue this atom belongs to

resnumint

The number of the residue this atom belongs to

chainstr

The chain ID character for this residue

inscodestr

The insertion code ID character for this residue (it is stripped)

segidstr

The segment identifier for this residue (it is stripped)

Notes

If the residue name and number differ from the last residue in this list, a new residue is added and the atom is added to that residue

prune()[source]

This function goes through the residue list and removes all empty residues from the list. This isn’t done automatically when atoms are deleted, since it will become very slow. You must remember to do this to avoid including empty residues

class parmed.StretchBend(atom1, atom2, atom3, type=None)[source]

Bases: object

This term models the stretching and bending of a standard valence angle, and is used in the AMOEBA force field

Parameters
atom1Atom

The first atom on one end of the angle

atom2Atom

The central atom in the angle

atom3Atom

The atom on the other end of the angle

typeStretchBendType

The type containing the stretch-bend parameters

Notes

Both :class:`Bond`s and :class:`Atom`s can be contained in a stretch-bend term

class parmed.StretchBendType(k1, k2, req1, req2, theteq, list=None)[source]

Bases: parmed.topologyobjects._ParameterType, parmed.topologyobjects._ListItem

A stretch-bend type with two distances and an angle in AMOEBA

Parameters
k1float

First force constant in kcal/mol/(radians*angstroms)

k2float

Second force constant in kcal/mol/(radians*angstroms)

req1float

Equilibrium bond distance for bond between the first and second atoms in Angstroms

req2float

Equilibrium bond distance for bond between the second and third atoms in Angstroms

theteqfloat

Equilibrium angle in degrees

listTrackedList

A list of `StretchBendType`s in which this is a member

Notes

Two StretchBendType`s are equal if their `req1, req2, theteq, and k attributes are equal

Examples

>>> sbt1 = StretchBendType(10.0, 10.0, 1.0, 1.0, 180.0)
>>> sbt2 = StretchBendType(10.0, 10.0, 1.0, 1.0, 180.0)
>>> sbt1 is sbt2
False
>>> sbt1 == sbt2
True
>>> sbt1.idx # Not part of any list or iterable
-1

As part of a list, they can be indexed

>>> strbnd_list = []
>>> strbnd_list.append(StretchBendType(10.0, 10.0, 1.0, 1.0, 180.0, strbnd_list))
>>> strbnd_list.append(StretchBendType(10.0, 10.0, 1.0, 1.0, 180.0, strbnd_list))
>>> strbnd_list[0].idx
0
>>> strbnd_list[1].idx
1
Attributes
idx
uk1
uk2
ureq1
ureq2
utheteq
property uk1
property uk2
property ureq1
property ureq2
property utheteq
class parmed.Structure[source]

Bases: object

A chemical structure composed of atoms, bonds, angles, torsions, and other topological features

Notes

This class also has a handful of type lists for each of the attributes above (excluding atoms, residues, chiral_frames, and multipole_frames). They are all TrackedList instances that are designed to hold the relevant parameter type. The list is:

bond_types, angle_types, dihedral_types, urey_bradley_types, improper_types, cmap_types, trigonal_angle_types, out_of_plane_bend_types, pi_torsion_types, stretch_bend_types, torsion_torsion_types, adjust_types

dihedral_types _may_ be a list of DihedralType instances, since torsion profiles are often represented by a Fourier series with multiple terms

Attributes
atomsAtomList

List of all atoms in the structure

residuesResidueList

List of all residues in the structure

bondsTrackedList (Bond)

List of all bonds in the structure

anglesTrackedList (Angle)

List of all angles in the structure

dihedralsTrackedList (Dihedral)

List of all dihedrals in the structure

rb_torsionsTrackedList (Dihedral)

List of all Ryckaert-Bellemans torsions in the structure

urey_bradleysTrackedList (UreyBradley)

List of all Urey-Bradley angle bends in the structure

impropersTrackedList (Improper)

List of all CHARMM-style improper torsions in the structure

cmapsTrackedList (Cmap)

List of all CMAP objects in the structure

trigonal_anglesTrackedList (TrigonalAngle)

List of all AMOEBA-style trigonal angles in the structure

out_of_plane_bendsTrackedList (OutOfPlaneBends)

List of all AMOEBA-style out-of-plane bending angles

pi_torsionsTrackedList (PiTorsion)

List of all AMOEBA-style pi-torsion angles

stretch_bendsTrackedList (StretchBend)

List of all AMOEBA-style stretch-bend compound bond/angle terms

torsion_torsionsTrackedList (TorsionTorsion)

List of all AMOEBA-style coupled torsion-torsion terms

chiral_framesTrackedList (ChiralFrame)

List of all AMOEBA-style chiral frames defined in the structure

multipole_framesTrackedList (MultipoleFrame)

List of all AMOEBA-style multipole frames defined in the structure

adjustsTrackedList (NonbondedException)

List of all nonbonded pair-exception rules

acceptorsTrackedList (AcceptorDonor)

List of all H-bond acceptors, if that information is present

donorsTrackedList (AcceptorDonor)

List of all H-bond donors, if that information is present

groupsTrackedList (Group)

List of all CHARMM-style GROUP objects (whatever those are used for)

boxlist of 6 floats

Box dimensions (a, b, c, alpha, beta, gamma) for the unit cell. If no box is defined, box is set to None

space_groupstr

The space group of the structure (default is “P 1”)

nrexclint

The number of bonds away that an atom can be in order to be excluded from the direct nonbonded summation

titlestr

Cosmetic only, it is an arbitrary title assigned to the system. Default value is an empty string

positionsu.Quantity(list(Vec3), u.angstroms)

A list of 3-element Quantity tuples of dimension length representing the atomic positions for every atom in the system.

coordinatesnp.ndarray of shape (nframes, natom, 3)

If no coordinates are set, this is set to None. The first frame will match the coordinates present on the atoms.

symmetrySymmetry

if no symmetry is set, this is set to None.

linksTrackedList (Link)

The list of Link definitions for this Structure

Methods

add_atom(atom, resname, resnum[, chain, …])

Adds a new atom to the Structure, adding a new residue to residues if it has a different name or number as the last residue added and adding it to the atoms list.

add_atom_to_residue(atom, residue)

Adds a new atom to the Structure at the end if the given residue

assign_bonds(*reslibs)

Assigns bonds to all atoms based on the provided residue template libraries.

copy(cls[, split_dihedrals])

Makes a copy of the current structure as an instance of a specified subclass

createSystem([nonbondedMethod, …])

Construct an OpenMM System representing the topology described by the prmtop file.

get_box([frame])

In some cases, multiple conformations may be stored in the Structure.

get_coordinates([frame])

In some cases, multiple conformations may be stored in the Structure.

has_NBFIX()

Returns whether or not any pairs of atom types have their LJ interactions modified by an NBFIX definition

is_changed()

Determines if any of the topology has changed for this structure

join_dihedrals()

Joins multi-term torsions into a single term and makes all of the parameters DihedralTypeList instances.

load_dataframe(df)

Loads atomic properties from an input DataFrame

omm_add_constraints(system, constraints, …)

Adds constraints to a given system

omm_angle_force([constraints, …])

Creates an OpenMM HarmonicAngleForce object (or AmoebaAngleForce if the angles are for an Amoeba-parametrized system)

omm_bond_force([constraints, rigidWater, …])

Creates an OpenMM Bond Force object (or AmoebaBondForce if the bonds are for an Amoeba-parametrized system)

omm_cmap_force()

Creates the OpenMM CMAP torsion force

omm_dihedral_force([split])

Creates the OpenMM PeriodicTorsionForce modeling dihedrals

omm_gbsa_force(implicitSolvent[, …])

Creates a Generalized Born force for running implicit solvent calculations

omm_improper_force()

Creates the OpenMM improper torsion force (quadratic bias)

omm_nonbonded_force([nonbondedMethod, …])

Creates the OpenMM NonbondedForce instance

omm_out_of_plane_bend_force()

Creates the Amoeba out-of-plane bend force

omm_pi_torsion_force()

Creates the Amoeba pi-torsion force

omm_rb_torsion_force()

Creates the OpenMM RBTorsionForce for Ryckaert-Bellemans torsions

omm_set_virtual_sites(system)

Sets the virtual sites in a given OpenMM System object from the extra points defined in this system

omm_stretch_bend_force()

Create the OpenMM Amoeba stretch-bend force for this system

omm_torsion_torsion_force()

Create the OpenMM Amoeba coupled-torsion (CMAP) force

omm_trigonal_angle_force()

Creates the Amoeba trigonal-angle force

omm_urey_bradley_force()

Creates the OpenMM Urey-Bradley force

prune_empty_terms()

Looks through all of the topological lists and gets rid of terms in which at least one of the atoms is None or has an idx attribute set to -1 (indicating that it has been removed from the atoms atom list)

save(fname[, format, overwrite])

Saves the current Structure in the requested file format.

split()

Split the current Structure into separate Structure instances for each unique molecule.

strip(selection)

Deletes a subset of the atoms corresponding to an atom-based selection.

to_dataframe()

Generates a DataFrame from the current Structure’s atomic properties

unchange()

Toggles all lists so that they do not indicate any changes

update_dihedral_exclusions()

Nonbonded exclusions and exceptions have the following priority:

visualize(*args, **kwargs)

Use nglview for visualization.

write_cif(dest[, renumber, coordinates, …])

Write a PDB file from the current Structure instance

write_pdb(dest[, renumber, coordinates, …])

Write a PDB file from a Structure instance

write_psf(dest[, vmd])

Writes a PSF file from the stored molecule

ANGLE_FORCE_GROUP = 1
BOND_FORCE_GROUP = 0
CMAP_FORCE_GROUP = 5
DIHEDRAL_FORCE_GROUP = 2
IMPROPER_FORCE_GROUP = 4
NONBONDED_FORCE_GROUP = 11
OUT_OF_PLANE_BEND_FORCE_GROUP = 7
PI_TORSION_FORCE_GROUP = 8
RB_TORSION_FORCE_GROUP = 12
STRETCH_BEND_FORCE_GROUP = 9
TORSION_TORSION_FORCE_GROUP = 10
TRIGONAL_ANGLE_FORCE_GROUP = 6
UREY_BRADLEY_FORCE_GROUP = 3
add_atom(atom, resname, resnum, chain='', inscode='', segid='')[source]

Adds a new atom to the Structure, adding a new residue to residues if it has a different name or number as the last residue added and adding it to the atoms list.

Parameters
atomAtom

The atom to add to this residue list

resnamestr

The name of the residue this atom belongs to

resnumint

The number of the residue this atom belongs to

chainstr

The chain ID character for this residue

inscodestr

The insertion code ID character for this residue (it is stripped)

segidstr

The segment identifier for this residue (it is stripped)

Notes

If the residue name and number differ from the last residue in this list, a new residue is added and the atom is added to that residue

add_atom_to_residue(atom, residue)[source]

Adds a new atom to the Structure at the end if the given residue

Parameters
atomAtom

The atom to add to the system

residueResidue

The residue to which to add this atom. It MUST be part of this Structure instance already or a ValueError is raised

Notes

This atom is added at the end of the residue and is inserted into the atoms list in such a way that all residues are composed of atoms contiguous in the atoms list. For large systems, this may be a relatively expensive operation

assign_bonds(*reslibs)[source]

Assigns bonds to all atoms based on the provided residue template libraries. Atoms whose names are not in the templates, as well as those residues for whom no template is found, is assigned to bonds based on distances.

Parameters
reslibsdict{str: ResidueTemplate}

Any number of residue template libraries. By default, assign_bonds knows about the standard amino acid, RNA, and DNA residues.

property box
property box_vectors

3, 3-element tuple of unit cell vectors that are Quantity objects of dimension length

property combining_rule
property coordinates
copy(cls, split_dihedrals=False)[source]

Makes a copy of the current structure as an instance of a specified subclass

Parameters
clsStructure subclass

The returned object is a copy of this structure as a cls instance

split_dihedralsbool

If True, then the Dihedral entries will be split up so that each one is paired with a single DihedralType (rather than a DihedralTypeList)

Returns
cls instance

The instance of the Structure subclass cls with a copy of the current Structure’s topology information

createSystem(nonbondedMethod=None, nonbondedCutoff=Quantity(value=8.0, unit=angstrom), switchDistance=Quantity(value=0.0, unit=angstrom), constraints=None, rigidWater=True, implicitSolvent=None, implicitSolventKappa=None, implicitSolventSaltConc=Quantity(value=0.0, unit=mole / liter), temperature=Quantity(value=298.15, unit=kelvin), soluteDielectric=1.0, solventDielectric=78.5, useSASA=False, removeCMMotion=True, hydrogenMass=None, ewaldErrorTolerance=0.0005, flexibleConstraints=True, verbose=False, splitDihedrals=False)[source]

Construct an OpenMM System representing the topology described by the prmtop file.

Parameters
nonbondedMethodcutoff method

This is the cutoff method. It can be either the NoCutoff, CutoffNonPeriodic, CutoffPeriodic, PME, or Ewald objects from the simtk.openmm.app namespace

nonbondedCutofffloat or distance Quantity

The nonbonded cutoff must be either a floating point number (interpreted as nanometers) or a Quantity with attached units. This is ignored if nonbondedMethod is NoCutoff.

switchDistancefloat or distance Quantity

The distance at which the switching function is turned on for van der Waals interactions. This is ignored when no cutoff is used, and no switch is used if switchDistance is 0, negative, or greater than the cutoff

constraintsNone, app.HBonds, app.HAngles, or app.AllBonds

Which type of constraints to add to the system (e.g., SHAKE). None means no bonds are constrained. HBonds means bonds with hydrogen are constrained

rigidWaterbool=True

If True, water is kept rigid regardless of the value of constraints. A value of False is ignored if constraints is not None.

implicitSolventNone, app.HCT, app.OBC1, app.OBC2, app.GBn, app.GBn2

The Generalized Born implicit solvent model to use.

implicitSolventKappafloat or 1/distance Quantity = None

This is the Debye kappa property related to modeling saltwater conditions in GB. It should have units of 1/distance (1/nanometers is assumed if no units present). A value of None means that kappa will be calculated from implicitSolventSaltConc (below)

implicitSolventSaltConcfloat or amount/volume Quantity=0 moles/liter

If implicitSolventKappa is None, the kappa will be computed from the salt concentration. It should have units compatible with mol/L

temperaturefloat or temperature Quantity = 298.15 kelvin

This is only used to compute kappa from implicitSolventSaltConc

soluteDielectricfloat=1.0

The dielectric constant of the protein interior used in GB

solventDielectricfloat=78.5

The dielectric constant of the water used in GB

useSASAbool=False

If True, use the ACE non-polar solvation model. Otherwise, use no SASA-based nonpolar solvation model.

removeCMMotionbool=True

If True, the center-of-mass motion will be removed periodically during the simulation. If False, it will not.

hydrogenMassfloat or mass quantity = None

If not None, hydrogen masses will be changed to this mass and the difference subtracted from the attached heavy atom (hydrogen mass repartitioning)

ewaldErrorTolerancefloat=0.0005

When using PME or Ewald, the Ewald parameters will be calculated from this value

flexibleConstraintsbool=True

If False, the energies and forces from the constrained degrees of freedom will NOT be computed. If True, they will (but those degrees of freedom will still be constrained).

verbosebool=False

If True, the progress of this subroutine will be printed to stdout

splitDihedralsbool=False

If True, the dihedrals will be split into two forces – proper and impropers. This is primarily useful for debugging torsion parameter assignments.

Notes

This function calls prune_empty_terms if any Topology lists have changed

get_box(frame='all')[source]

In some cases, multiple conformations may be stored in the Structure. This function retrieves a particular frame’s unit cell (box) dimensions

Parameters
frameint or ‘all’, optional

The frame number whose unit cell should be retrieved. Default is ‘all’

Returns
boxnp.ndarray, shape([#,] 6) or None

If frame is ‘all’, all unit cells are returned with shape (#, 6). Otherwise the requested frame is returned with shape (6,). If no unit cell exist and ‘all’ is requested, None is returned

Raises
IndexError if there are fewer than frame unit cell dimensions
get_coordinates(frame='all')[source]

In some cases, multiple conformations may be stored in the Structure. This function retrieves a particular frame’s coordinates

Parameters
frameint or ‘all’, optional

The frame number whose coordinates should be retrieved. Default is ‘all’

Returns
coordsnp.ndarray, shape([#,] natom, 3) or None

If frame is ‘all’, all coordinates are returned with shape (#, natom, 3). Otherwise the requested frame is returned with shape (natom, 3). If no coordinates exist and ‘all’ is requested, None is returned

Raises
IndexError if there are fewer than frame coordinates
has_NBFIX()[source]

Returns whether or not any pairs of atom types have their LJ interactions modified by an NBFIX definition

Returns
has_nbfixbool

If True, at least two atom types have NBFIX mod definitions

is_changed()[source]

Determines if any of the topology has changed for this structure

join_dihedrals()[source]

Joins multi-term torsions into a single term and makes all of the parameters DihedralTypeList instances. If any dihedrals are already DihedralTypeList instances, or any are not parametrized, or there are no dihedral_types, this method returns without doing anything

load_dataframe(df)[source]

Loads atomic properties from an input DataFrame

Parameters
dfpandas.DataFrame

A pandas DataFrame with atomic properties that will be used to set the properties on the current list of atoms

omm_add_constraints(system, constraints, rigidWater)[source]

Adds constraints to a given system

Parameters
systemmm.System

The OpenMM system for which constraints should be added

constraintsNone, app.HBonds, app.AllBonds, or app.HAngles

Which kind of constraints should be used

rigidWaterbool

If True, water bonds are constrained regardless of whether constrains is None

omm_angle_force(constraints=None, flexibleConstraints=True)[source]

Creates an OpenMM HarmonicAngleForce object (or AmoebaAngleForce if the angles are for an Amoeba-parametrized system)

Parameters
constraintsNone, app.HBonds, app.AllBonds, or app.HAngles

The types of constraints that are on the system. If flexibleConstraints is False, then the constrained bonds will not be added to the resulting Force

flexibleConstraintsbool=True

If True, all bonds are added to the force regardless of constraints

Returns
force

HarmonicAngleForce (or AmoebaAngleForce if this is an Amoeba system), or None if there are no angles to add

omm_bond_force(constraints=None, rigidWater=True, flexibleConstraints=True)[source]

Creates an OpenMM Bond Force object (or AmoebaBondForce if the bonds are for an Amoeba-parametrized system)

Parameters
constraintsNone, app.HBonds, app.AllBonds, or app.HAngles

The types of constraints that are on the system. If flexibleConstraints is False, then the constrained bonds will not be added to the resulting Force

rigidWaterbool=True

Should water-H bonds be constrained regardless of constraints?

flexibleConstraintsbool=True

If True, all bonds are added to the force regardless of constraints

Returns
force

HarmonicBondForce (or AmoebaBondForce if this is an Amoeba system), or None if there are no bonds to add

omm_cmap_force()[source]

Creates the OpenMM CMAP torsion force

Returns
CMAPTorsionForce

Or None, if no CMAP terms are present

omm_dihedral_force(split=False)[source]

Creates the OpenMM PeriodicTorsionForce modeling dihedrals

Parameters
splitbool, optional, default=False

If True, separate PeriodicTorsionForce instances with the propers in the first and impropers in the second return item. If no impropers or propers are present, the instances with zero terms are not returned.

Returns
PeriodicTorsionForce[, PeriodicTorsionForce]

Or None if no torsions are present in this system

omm_gbsa_force(implicitSolvent, nonbondedMethod=None, nonbondedCutoff=Quantity(value=30.0, unit=angstrom), soluteDielectric=1.0, solventDielectric=78.5, implicitSolventKappa=None, implicitSolventSaltConc=Quantity(value=0.0, unit=mole / liter), temperature=Quantity(value=298.15, unit=kelvin), useSASA=True)[source]

Creates a Generalized Born force for running implicit solvent calculations

Parameters
implicitSolventapp.HCT, app.OBC1, app.OBC2, app.GBn, app.GBn2

The Generalized Born implicit solvent model to use.

nonbondedMethodcutoff method

This is the cutoff method. It can be either the NoCutoff, CutoffNonPeriodic, CutoffPeriodic, PME, or Ewald objects from the simtk.openmm.app namespace. Default is NoCutoff

nonbondedCutofffloat or distance Quantity

The nonbonded cutoff must be either a floating opint number (interpreted as nanometers) or a Quantity with attached units. This is ignored if nonbondedMethod is NoCutoff

implicitSolventKappafloat or 1/distance Quantity = None

This is the Debye kappa property related to modeling saltwater conditions in GB. It should have units of 1/distance (1/nanometers is assumed if no units present). A value of None means that kappa will be calculated from implicitSolventSaltConc (below)

implicitSolventSaltConcfloat or amount/volume Quantity=0 moles/liter

If implicitSolventKappa is None, the kappa will be computed from the salt concentration. It should have units compatible with mol/L

temperaturefloat or temperature Quantity = 298.15 kelvin

This is only used to compute kappa from implicitSolventSaltConc

soluteDielectricfloat=1.0

The dielectric constant of the protein interior used in GB

solventDielectricfloat=78.5

The dielectric constant of the water used in GB

omm_improper_force()[source]

Creates the OpenMM improper torsion force (quadratic bias)

Returns
CustomTorsionForce

With the formula k*(phi-phi0)^2, or None if there are no impropers

omm_nonbonded_force(nonbondedMethod=None, nonbondedCutoff=Quantity(value=8, unit=angstrom), switchDistance=Quantity(value=0, unit=angstrom), ewaldErrorTolerance=0.0005, reactionFieldDielectric=78.5)[source]

Creates the OpenMM NonbondedForce instance

Parameters
nonbondedMethodcutoff method

This is the cutoff method. It can be either the NoCutoff, CutoffNonPeriodic, CutoffPeriodic, PME, or Ewald objects from the simtk.openmm.app namespace

nonbondedCutofffloat or distance Quantity

The nonbonded cutoff must be either a floating point number (interpreted as nanometers) or a Quantity with attached units. This is ignored if nonbondedMethod is NoCutoff.

switchDistancefloat or distance Quantity

The distance at which the switching function is turned on for van der Waals interactions. This is ignored when no cutoff is used, and no switch is used if switchDistance is 0, negative, or greater than the cutoff

ewaldErrorTolerancefloat=0.0005

When using PME or Ewald, the Ewald parameters will be calculated from this value

reactionFieldDielectricfloat=78.5

If the nonbondedMethod is CutoffPeriodic or CutoffNonPeriodic, the region beyond the cutoff is treated using a reaction field method with this dielectric constant. It should be set to 1 if another implicit solvent model is being used (e.g., GB)

Returns
NonbondedForce

This just implements the very basic NonbondedForce with the typical charge-charge and 12-6 Lennard-Jones interactions with the Lorentz-Berthelot combining rules.

Notes

Subclasses of Structure for which this nonbonded treatment is inadequate should override this method to implement what is needed.

If nrexcl is set to 3 and no exception parameters are stored in the adjusts list, the 1-4 interactions are determined from the list of dihedrals

omm_out_of_plane_bend_force()[source]

Creates the Amoeba out-of-plane bend force

Returns
AmoebaOutOfPlaneBendForce

The out-of-plane bend Angle force

omm_pi_torsion_force()[source]

Creates the Amoeba pi-torsion force

Returns
AmoebaPiTorsionForce

The pi-torsion force

omm_rb_torsion_force()[source]

Creates the OpenMM RBTorsionForce for Ryckaert-Bellemans torsions

Returns
RBTorsionForce

Or None if no torsions are present in this system

omm_set_virtual_sites(system)[source]

Sets the virtual sites in a given OpenMM System object from the extra points defined in this system

Parameters
systemmm.System

The system for which the virtual sites will be set. All particles must have already been added to this System before calling this method

omm_stretch_bend_force()[source]

Create the OpenMM Amoeba stretch-bend force for this system

Returns
AmoebaStretchBendForce

The stretch-bend force containing all terms in this system

omm_torsion_torsion_force()[source]

Create the OpenMM Amoeba coupled-torsion (CMAP) force

Returns
AmoebaTorsionTorsionForce

The torsion-torsion (CMAP) force with all coupled-torsion parameters for this system

omm_trigonal_angle_force()[source]

Creates the Amoeba trigonal-angle force

Returns
AmoebaInPlaneAngleForce

The trigonal in-plane Angle force

omm_urey_bradley_force()[source]

Creates the OpenMM Urey-Bradley force

Returns
HarmonicBondForce

Or None, if no urey-bradleys are present

property positions

A list of 3-element Quantity tuples of dimension length representing the atomic positions for every atom in the system. If set with unitless numbers, those numbers are assumed to be in angstroms. If any atoms do not have coordinates, this is simply None.

prune_empty_terms()[source]

Looks through all of the topological lists and gets rid of terms in which at least one of the atoms is None or has an idx attribute set to -1 (indicating that it has been removed from the atoms atom list)

save(fname, format=None, overwrite=False, **kwargs)[source]

Saves the current Structure in the requested file format. Supported formats can be specified explicitly or determined by file-name extension. The following formats are supported, with the recognized suffix and format keyword shown in parentheses:

  • PDB (.pdb, pdb)

  • PDBx/mmCIF (.cif, cif)

  • PQR (.pqr, pqr)

  • Amber topology file (.prmtop/.parm7, amber)

  • CHARMM PSF file (.psf, psf)

  • CHARMM coordinate file (.crd, charmmcrd)

  • Gromacs topology file (.top, gromacs)

  • Gromacs GRO file (.gro, gro)

  • Mol2 file (.mol2, mol2)

  • Mol3 file (.mol3, mol3)

  • Amber ASCII restart (.rst7/.inpcrd/.restrt, rst7)

  • Amber NetCDF restart (.ncrst, ncrst)

Parameters
fnamestr or file-like object

Name of the file or file-like object to save. If format is None (see below), the file type will be determined based on the filename extension. If fname is file-like object, format must be provided. If the type cannot be determined, a ValueError is raised.

formatstr, optional

The case-insensitive keyword specifying what type of file fname should be saved as. If None (default), the file type will be determined from filename extension of fname

overwritebool, optional

If True, allow the target file to be overwritten. Otherwise, an IOError is raised if the file exists. Default is False

kwargskeyword-arguments

Remaining arguments are passed on to the file writing routines that are called by this function

Raises
ValueError if either filename extension or format are not recognized
TypeError if the structure cannot be converted to the desired format for
whatever reason
IOError if the file cannot be written either because it exists and
overwrite is False, the filesystem is read-only, or write
permissions are not granted for the user
split()[source]

Split the current Structure into separate Structure instances for each unique molecule. A molecule is defined as all atoms connected by a graph of covalent bonds.

Returns
[structs, counts]list of (Structure, list) tuples

List of all molecules in the order that they appear first in the parent structure accompanied by the list of the molecule numbers in which that molecule appears in the Structure

strip(selection)[source]

Deletes a subset of the atoms corresponding to an atom-based selection.

Parameters
selectionAmberMask, str, or iterable

This is the selection of atoms that will be deleted from this structure. If it is a string, it will be interpreted as an AmberMask. If it is an AmberMask, it will be converted to a selection of atoms. If it is an iterable, it must be the same length as the atoms list.

to_dataframe()[source]

Generates a DataFrame from the current Structure’s atomic properties

Returns
dfDataFrame

DataFrame with all atomic properties

property topology

The OpenMM Topology object. Cached when possible, but any changes to the Structure instance results in the topology being deleted and rebuilt

Notes

This function calls prune_empty_terms if any topology lists have changed.

unchange()[source]

Toggles all lists so that they do not indicate any changes

update_dihedral_exclusions()[source]

Nonbonded exclusions and exceptions have the following priority:

bond -> angle -> dihedral

Since bonds and angles are completely excluded, any ring systems in which two atoms are attached by a bond or angle as well as a dihedral should be completely excluded as the bond and angle exclusion rules take precedence. If a Bond or Angle was _added_ to the structure between a pair of atoms previously connected only by a dihedral term, it’s possible that those two atoms have both an exclusion and an exception defined. The result of this scenario is that sander and pmemd will happily compute an energy, _including_ the 1-4 nonbonded terms between atoms now connected by a bond or an Angle. OpenMM, on the other hand, will complain about an exception specified multiple times. This method scans through all of the dihedrals in which ignore_end is False and turns it to True if the two end atoms are in the bond or angle partners arrays

property velocities
property view

Returns an indexable object that can be indexed like a standard Structure, but returns a view rather than a copy

See also

Structure.__getitem__
visualize(*args, **kwargs)[source]

Use nglview for visualization. This only works with Jupyter notebook and require to install nglview

Parameters
args and kwargspositional and keyword arguments given to nglview, optional

Examples

>>> import parmed as pmd
>>> parm = pmd.download_PDB('1tsu')
>>> parm.visualize()
write_cif(dest, renumber=True, coordinates=None, altlocs='all', write_anisou=False, standard_resnames=False)

Write a PDB file from the current Structure instance

Parameters
structStructure

The structure from which to write the PDBx/mmCIF file

deststr or file-like

Either a file name or a file-like object containing a write method to which to write the PDB file. If it is a filename that ends with .gz or .bz2, a compressed version will be written using either gzip or bzip2, respectively.

renumberbool

If True, renumber the atoms and residues sequentially as they are stored in the structure. If False, use the original numbering if it was assigned previously

coordinatesarray-like of float

If provided, these coordinates will be written to the PDB file instead of the coordinates stored in the structure. These coordinates should line up with the atom order in the structure (not necessarily the order of the “original” PDB file if they differ)

altlocsstr

Keyword controlling which alternate locations are printed to the resulting PDB file. Allowable options are:

  • ‘all’ : (default) print all alternate locations

  • ‘first’ : print only the first alternate locations

  • ‘occupancy’ : print the one with the largest occupancy. If two conformers have the same occupancy, the first one to occur is printed

Input is case-insensitive, and partial strings are permitted as long as it is a substring of one of the above options that uniquely identifies the choice.

write_anisoubool

If True, an ANISOU record is written for every atom that has one. If False, ANISOU records are not written

standard_resnamesbool, optional

If True, common aliases for various amino and nucleic acid residues will be converted into the PDB-standard values. Default is False

Notes

If multiple coordinate frames are present, these will be written as separate models (but only the unit cell from the first model will be written, as the PDBx standard dictates that only one set of unit cells shall be present).

write_pdb(dest, renumber=True, coordinates=None, altlocs='all', write_anisou=False, charmm=False, use_hetatoms=True, standard_resnames=False, increase_tercount=True, write_links=False)

Write a PDB file from a Structure instance

Parameters
structStructure

The structure from which to write the PDB file

deststr or file-like

Either a file name or a file-like object containing a write method to which to write the PDB file. If it is a filename that ends with .gz or .bz2, a compressed version will be written using either gzip or bzip2, respectively.

renumberbool, optional, default True

If True, renumber the atoms and residues sequentially as they are stored in the structure. If False, use the original numbering if it was assigned previously.

coordinatesarray-like of float, optional

If provided, these coordinates will be written to the PDB file instead of the coordinates stored in the structure. These coordinates should line up with the atom order in the structure (not necessarily the order of the “original” PDB file if they differ)

altlocsstr, optional, default ‘all’

Keyword controlling which alternate locations are printed to the resulting PDB file. Allowable options are:

  • ‘all’ : print all alternate locations

  • ‘first’ : print only the first alternate locations

  • ‘occupancy’ : print the one with the largest occupancy. If two conformers have the same occupancy, the first one to occur is printed

Input is case-insensitive, and partial strings are permitted as long as it is a substring of one of the above options that uniquely identifies the choice.

write_anisoubool, optional, default False

If True, an ANISOU record is written for every atom that has one. If False, ANISOU records are not written.

charmmbool, optional, default False

If True, SEGID will be written in columns 73 to 76 of the PDB file in the typical CHARMM-style PDB output. This will be omitted for any atom that does not contain a SEGID identifier.

use_hetatoms: bool, optional, default True

If True, certain atoms will have the HETATM tag instead of ATOM as per the PDB-standard.

standard_resnamesbool, optional, default False

If True, common aliases for various amino and nucleic acid residues will be converted into the PDB-standard values.

increase_tercountbool, optional, default True

If True, the TER atom number field increased by one compared to atom card preceding it; this conforms to PDB standard.

write_linksbool, optional, default False

If True, any LINK records stored in the Structure will be written to the LINK records near the top of the PDB file. If this is True, then renumber must be False or a ValueError will be thrown

Notes

If multiple coordinate frames are present, these will be written as separate models (but only the unit cell from the first model will be written, as the PDB standard dictates that only one set of unit cells shall be present).

write_psf(dest, vmd=False)

Writes a PSF file from the stored molecule

Parameters
structStructure

The Structure instance from which the PSF should be written

deststr or file-like

The place to write the output PSF file. If it has a “write” attribute, it will be used to print the PSF file. Otherwise, it will be treated like a string and a file will be opened, printed, then closed

vmdbool

If True, it will write out a PSF in the format that VMD prints it in (i.e., no NUMLP/NUMLPH or MOLNT sections)

Examples

>>> cs = CharmmPsfFile('testfiles/test.psf')
>>> cs.write_psf('testfiles/test2.psf')
class parmed.StructureView[source]

Bases: object

A view of a Structure. In many cases, this can serve as a duck-typed Structure and it has many of the same attributes. However, none of its lists own their objects, and the lists of atoms, residues, and parameters/topologies are regular lists, rather than TrackedList instances. Therefore, the indexes correspond to the indexes from the original Structure from which this view was taken. Furthermore, there are no “type” lists, as they would be exactly equivalent to the type lists of the parent Structure instance.

Attributes
atomsAtomList

List of all atoms in the structure

residuesResidueList

List of all residues in the structure

bondsTrackedList (Bond)

List of all bonds in the structure

anglesTrackedList (Angle)

List of all angles in the structure

dihedralsTrackedList (Dihedral)

List of all dihedrals in the structure

rb_torsionsTrackedList (RBTorsion)

List of all Ryckaert-Bellemans torsions in the structure

urey_bradleysTrackedList (UreyBradley)

List of all Urey-Bradley angle bends in the structure

impropersTrackedList (Improper)

List of all CHARMM-style improper torsions in the structure

cmapsTrackedList (Cmap)

List of all CMAP objects in the structure

trigonal_anglesTrackedList (TrigonalAngle)

List of all AMOEBA-style trigonal angles in the structure

out_of_plane_bendsTrackedList (OutOfPlaneBends)

List of all AMOEBA-style out-of-plane bending angles

pi_torsionsTrackedList (PiTorsion)

List of all AMOEBA-style pi-torsion angles

stretch_bendsTrackedList (StretchBend)

List of all AMOEBA-style stretch-bend compound bond/angle terms

torsion_torsionsTrackedList (TorsionTorsion)

List of all AMOEBA-style coupled torsion-torsion terms

chiral_framesTrackedList (ChiralFrame)

List of all AMOEBA-style chiral frames defined in the structure

multipole_framesTrackedList (MultipoleFrame)

List of all AMOEBA-style multipole frames defined in the structure

adjustsTrackedList (NonbondedException)

List of all nonbonded pair-exception rules

acceptorsTrackedList (AcceptorDonor)

List of all H-bond acceptors, if that information is present

donorsTrackedList (AcceptorDonor)

List of all H-bond donors, if that information is present

positionsu.Quantity(list(Vec3), u.angstroms)

Unit-bearing atomic coordinates. If not all atoms have coordinates, this property is None

coordinatesnp.ndarray of shape (nframes, natom, 3)

If no coordinates are set, this is set to None. The first frame will match the coordinates present on the atoms.

Methods

load_dataframe(df)

Loads atomic properties from an input DataFrame

to_dataframe()

Generates a DataFrame from the current Structure’s atomic properties

property coordinates
load_dataframe(df)[source]

Loads atomic properties from an input DataFrame

Parameters
dfpandas.DataFrame

A pandas DataFrame with atomic properties that will be used to set the properties on the current list of atoms

property positions
to_dataframe()[source]

Generates a DataFrame from the current Structure’s atomic properties

Returns
dfDataFrame

DataFrame with all atomic properties

class parmed.ThreeParticleExtraPointFrame(ep, inside=True)[source]

Bases: object

This class defines a frame of reference for a given extra point with a frame of reference defined by 3 particles

Parameters
epExtraPoint

The extra point defined by this frame

insidebool

If True, the extra point is contained inside the curve connecting all points in the frame of reference. If False, the point is outside that curve. See figures below for example

Methods

from_weights(parent, a1, a2, w1, w2[, dp1, …])

This function determines the necessary bond length between an ExtraPoint and its parent atom from the weights that are calculated in get_weights.

get_atoms()

Returns the particles involved in the frame

get_weights()

Returns the weights for the three particles

static from_weights(parent, a1, a2, w1, w2, dp1=None, dp2=None, theteq=None, d12=None)[source]

This function determines the necessary bond length between an ExtraPoint and its parent atom from the weights that are calculated in get_weights.

Parameters
parentAtom

The parent atom to the ExtraPoint

a1Atom

The first atom in the frame bonded to the parent

a2Atom

The second atom in the frame bonded to the parent

w1float

The first weight defining the ExtraPoint position wrt a1

w2float

The second weight defining the ExtraPoint position wrt a2

dp1float, optional

Equilibrium distance between parent and a1. If None, a bond with a bond type must be defined between parent and a1, and the distance will be taken from there. Units must be Angstroms. Default is None

dp2float, optional

Same as dp1 above, but between atoms parent and a2. Either both dp1 and dp2 should be None or neither should be. If one is None, the other will be ignored if not None.

theteqfloat, optional

Angle between bonds parent-a1 and parent-a2. If None, d12 (below) will be used to calculate the angle. If both are None, the angle or d12 distance will be calculated from parametrized bonds or angles between a1, parent, and a2. Units of degrees. Default None.

d12float, optional

Distance between a1 and a2 as an alternative way to specify theteq. Default is None.

Returns
distfloat

The distance between the ExtraPoint and its parent atom that satisfies the weights

Raises
ValueError if the necessary geometry requirements are not set or if the
two weights are different

Notes

parent must form a Bond with both a1 and a2. Then, if a1-parent-a2 forms an angle and it has an assigned type, that equilibrium value is used. Otherwise, the a1-a2 bond distance is used. If neither of those are defined, a ValueError is raised.

get_atoms()[source]

Returns the particles involved in the frame

Returns
a1, a2, a3Atom, Atom, Atom

a1 is the parent atom of the extra point. a2 and a3 are the other atoms that are both bonded to the parent atom

get_weights()[source]

Returns the weights for the three particles

Returns
w1, w2, w3float, float, float

w1 is the weight of particle 1 (the parent atom), whereas w2 and w3 are the weights of particles 2 and 3 (the atoms bonded to the parent atom)

class parmed.TorsionTorsion(atom1, atom2, atom3, atom4, atom5, type=None)[source]

Bases: parmed.topologyobjects.Cmap

This is a coupled-torsion map used in the AMOEBA force field similar to the correction-map (CMAP) potential used by the CHARMM force field

Parameters
atom1Atom

An atom on one end of the valence torsion-torsion bonded to atom2

atom2Atom

An atom in the middle of the torsion-torsion bonded to atoms 1 and 3

atom3Atom

An atom in the middle of the torsion-torsion bonded to atoms 2 and 4

atom4Atom

An atom in the middle of the torsion-torsion bonded to atoms 3 and 5

atom5Atom

An atom in the middle of the torsion-torsion bonded to atom 4

typeTorsionTorsionType

The TorsionTorsionType object containing the parameter map for this term

Notes

A TorsionTorsion can contain bonds or atoms. A bond is contained if it exists between atoms 1 and 2, between atoms 2 and 3, between atoms 3 and 4, or between atoms 4 and 5.

Examples

>>> a1, a2, a3, a4, a5 = Atom(), Atom(), Atom(), Atom(), Atom()
>>> tortor = TorsionTorsion(a1, a2, a3, a4, a5)
>>> Bond(a1, a2) in tortor and Bond(a2, a3) in tortor
True
>>> Bond(a1, a3) in tortor
False

Methods

delete()

Deletes this TorsionTorsion from the atoms that make it up.

extended(atom1, atom2, atom3, atom4, atom5, …)

Alternative constructor for correction maps defined with 8 atoms (each torsion being separately specified).

same_atoms(thing)

A coupled-torsion is equivalent if the 5 atoms are in the same or reverse order Allow comparison with another type of cmap or with a sequence of 5 indexes

delete()[source]

Deletes this TorsionTorsion from the atoms that make it up. This method removes the TorsionTorsion from the tortors list for atom1, atom2, atom3, atom4, and atom5, and removes each atom from the others’ tortor_partners list.

class parmed.TorsionTorsionType(dims, ang1, ang2, f, dfda1=None, dfda2=None, d2fda1da2=None, list=None)[source]

Bases: parmed.topologyobjects._ParameterType, parmed.topologyobjects._ListItem

The type containing the parameter maps for the Amoeba torsion-torsion potentials. It contains the original potential as well as interpolated first and second derivatives for the AMOEBA force field.

Parameters
dimstuple of 2 ints

The table dimensions

ang1list of floats

The list of angles in the first dimension

ang2list of floats

The list of angles in the second dimension

flist of floats

The interpolation table for the energy

dfda1list of floats

The interpolation table of the gradient w.r.t. angle 1

dfda2list of floats

The interpolation table of the gradient w.r.t. angle 2

d2fda1da2list of floats

The interpolation table of the 2nd derivative w.r.t. both angles

listTrackedList

The list containing this coupled torsion-torsion map

Attributes
dimstuple of 2 ints

The table dimensions

ang1list of floats

The list of angles in the first dimension

ang2list of floats

The list of angles in the second dimension

f_TorTorTable

The interpolation table for the energy as a _TorTorTable

dfda1_TorTorTable

The interpolation table for the first gradient as a _TorTorTable

dfda2_TorTorTable

The interpolation table for the second gradient as a _TorTorTable

d2fda1da2_TorTorTable

The interpolation table for the second derivative as a _TorTorTable

listTrackedList

The list that may, or may not, contain this TorsionTorsionType

idxint

The index of this item in the list or iterable defined by list

class parmed.TrackedList(*args)[source]

Bases: list

This creates a list type that allows you to see if anything has changed

Examples

>>> tl = TrackedList()
>>> tl.append(Atom())
>>> tl.append(Atom())
>>> tl.append(Atom())
>>> tl.needs_indexing, tl.changed
(True, True)
>>> tl.index_members()
>>> tl.needs_indexing, tl.changed
(False, True)
>>> tl.changed = False # Must do when changes have been incorporated
>>> tl.needs_indexing, tl.changed
(False, False)
Attributes
changedbool

Determines if something has been done to fundamentally change the underlying topology defined by this list such that the topology needs to be rebuilt

needs_indexingbool

A flag to determine whether or not the items in a tracked list need to be indexed or not.

Methods

claim()

This method causes this list to “claim” all of the items it contains and subsequently indexes all of its items.

clear(/)

Remove all items from list.

copy(/)

Return a shallow copy of the list.

count(value, /)

Return number of occurrences of value.

index(value[, start, stop])

Return first index of value.

index_members()

Assigns the idx variable for every member of this list to its place in the list, if the members of this list permit

prune_unused()

This method inspects the used attribute of all of its members, if it has one, and deletes any item in which it is set to False

reverse(/)

Reverse IN PLACE.

sort(*[, key, reverse])

Stable sort IN PLACE.

append

extend

insert

pop

remove

append(*args, **kwargs)

Append object to the end of the list.

claim()[source]

This method causes this list to “claim” all of the items it contains and subsequently indexes all of its items.

extend(*args, **kwargs)

Extend list by appending elements from the iterable.

index_members()[source]

Assigns the idx variable for every member of this list to its place in the list, if the members of this list permit

insert(*args, **kwargs)

Insert object before index.

pop(*args, **kwargs)[source]

Remove and return item at index (default last).

Raises IndexError if list is empty or index is out of range.

prune_unused()[source]

This method inspects the used attribute of all of its members, if it has one, and deletes any item in which it is set to False

remove(*args, **kwargs)[source]

Remove first occurrence of value.

Raises ValueError if the value is not present.

class parmed.TrigonalAngle(atom1, atom2, atom3, atom4, type=None)[source]

Bases: parmed.topologyobjects._FourAtomTerm

A trigonal-angle term in the AMOEBA force field. It exists in a pattern like the one shown below

A1 | |

A4—-A2—-A3

Parameters
atom1Atom

The first atom involved in the trigonal angle

atom2Atom

The central atom involved in the trigonal angle

atom3Atom

The third atom involved in the trigonal angle

atom4Atom

The fourth atom involved in the trigonal angle

typeAngleType

The angle type containing the parameters

Notes

Either `Atom`s or `Bond`s can be contained within this trigonal angle

Methods

delete()

Sets all atoms in this term to None, and its type if it exists

class parmed.TwoParticleExtraPointFrame(ep, inside=False)[source]

Bases: object

This class defines a frame of reference for a given extra point with a frame of reference defined by 2 particles

Parameters
epExtraPoint

The extra point defined by this frame

insidebool

If True, the extra point is contained inside the curve connecting all points in the frame of reference. If False, the point is outside that curve. See figures below for example

Methods

get_atoms()

Returns the particles involved in the frame

get_weights()

Returns the weights for the two particles

get_atoms()[source]

Returns the particles involved in the frame

Returns
a1, a2Atom, Atom

a1 is the parent atom of the extra point. a2 is the other atom bonded to the parent atom

get_weights()[source]

Returns the weights for the two particles

Returns
w1, w2float, float

w1 is the weight of particle 1 (the parent atom), whereas w2 is the weight of particle 2 (the atom bonded to the parent atom)

class parmed.UreyBradley(atom1, atom2, type=None)[source]

Bases: parmed.topologyobjects.Bond

A Urey-Bradley angle type with a set of parameters. It has the same functional form as a bond, but it is defined between two atoms forming a valence angle separated by two bonds.

Parameters
atom1Atom

The first atom involved in the Urey-Bradley bond

atom2Atom

The other atom involved in the Urey-Bradley bond

typeBondType

The Urey-Bradley bond type that defines the parameters for this bond

Notes

You can test whether an Atom is contained within the bond using the in operator. A MoleculeError is raised if atom1 and atom2 are identical. You can also test that a Bond is contained in this Urey-Bradley valence angle

Examples

>>> a1, a2, a3 = Atom(), Atom(), Atom()
>>> b1 = Bond(a1, a2)
>>> b2 = Bond(a2, a3)
>>> angle = Angle(a1, a2, a3)
>>> urey = UreyBradley(a1, a3)
>>> a1 in urey and a3 in urey
True
>>> b1 in urey and b2 in urey
True
Attributes
order

Bond order.

Methods

delete()

Deletes this Urey-Bradley from the atoms that make it up.

energy()

Measures the current bond energy

measure()

Measures the current bond

uenergy()

Same as energy(), but with units

umeasure()

Same as “measure”, but with units

delete()[source]

Deletes this Urey-Bradley from the atoms that make it up. This method removes this urey-bradley from the urey_bradleys list for atom1 and atom2.

class parmed.Vec3(x, y, z)[source]

Bases: parmed.vec3.Vec3

Vec3 is a 3-element tuple that supports many math operations.

Attributes
x

Alias for field number 0

y

Alias for field number 1

z

Alias for field number 2

Methods

count(value, /)

Return number of occurrences of value.

index(value[, start, stop])

Return first index of value.

parmed.download_CIF(pdb_id, timeout=10, saveto=None)

Goes to the wwPDB website and downloads the requested PDBx/mmCIF, loading it as a Structure instance

Parameters
pdb_idstr

The 4-letter PDB ID to try and download from the RCSB PDB database

timeoutfloat, optional

The number of seconds to wait before raising a timeout error. Default is 10 seconds

savetostr, optional

If provided, this will be treated as a file name to which the PDB file will be saved. If None (default), no CIF file will be written. This will be a verbatim copy of the downloaded CIF file, unlike the somewhat-stripped version you would get by using Structure.write_cif

Returns
structStructure

Structure instance populated by the requested PDBx/mmCIF

Raises
socket.timeout if the connection times out while trying to contact the
FTP server
IOError if there is a problem retrieving the requested PDB
ImportError if the gzip module is not available
TypeError if pdb_id is not a 4-character string
parmed.download_PDB(pdb_id, timeout=10, saveto=None)

Goes to the wwPDB website and downloads the requested PDB, loading it as a Structure instance

Parameters
pdb_idstr

The 4-letter PDB ID to try and download from the RCSB PDB database

timeoutfloat, optional

The number of seconds to wait before raising a timeout error. Default is 10 seconds

savetostr, optional

If provided, this will be treated as a file name to which the PDB file will be saved. If None (default), no PDB file will be written. This will be a verbatim copy of the downloaded PDB file, unlike the somewhat-stripped version you would get by using Structure.write_pdb

Returns
structStructure

Structure instance populated by the requested PDB

Raises
socket.timeout if the connection times out while trying to contact the
FTP server
IOError if there is a problem retrieving the requested PDB or writing a
requested saveto file
ImportError if the gzip module is not available
TypeError if pdb_id is not a 4-character string
parmed.load_file(filename, *args, **kwargs)[source]

Identifies the file format of the specified file and returns its parsed contents.

Parameters
filenamestr

The name of the file to try to parse. If the filename starts with http:// or https:// or ftp://, it is treated like a URL and the file will be loaded directly from its remote location on the web

structureobject, optional

For some classes, such as the Mol2 file class, the default return object is not a Structure, but can be made to return a Structure if the structure=True keyword argument is passed. To facilitate writing easy code, the structure keyword is always processed and only passed on to the correct file parser if that parser accepts the structure keyword. There is no default, as each parser has its own default.

natomint, optional

This is needed for some coordinate file classes, but not others. This is treated the same as structure, above. It is the # of atoms expected

hasboxbool, optional

Same as structure, but indicates whether the coordinate file has unit cell dimensions

skip_bondsbool, optional

Same as structure, but indicates whether or not bond searching will be skipped if the topology file format does not contain bond information (like PDB, GRO, and PQR files).

*argsother positional arguments

Some formats accept positional arguments. These will be passed along

**kwargsother options

Some formats can only be instantiated with other options besides just a file name.

Returns
object

The returned object is the result of the parsing function of the class associated with the file format being parsed

Raises
IOError

If filename does not exist

parmed.exceptions.FormatNotFound

If no suitable file format can be identified, a TypeError is raised

TypeError

If the identified format requires additional arguments that are not provided as keyword arguments in addition to the file name

Notes

Compressed files are supported and detected by filename extension. This applies both to local and remote files. The following names are supported:

  • .gz : gzip compressed file

  • .bz2 : bzip2 compressed file

SDF file is loaded via rdkit package.

Examples

Load a Mol2 file

>>> load_file('tripos1.mol2')
<ResidueTemplate DAN: 31 atoms; 33 bonds; head=None; tail=None>

Load a Mol2 file as a Structure

>>> load_file('tripos1.mol2', structure=True)
<Structure 31 atoms; 1 residues; 33 bonds; NOT parametrized>

Load an Amber topology file

>>> load_file('trx.prmtop', xyz='trx.inpcrd')
<AmberParm 1654 atoms; 108 residues; 1670 bonds; parametrized>

Load a CHARMM PSF file

>>> load_file('ala_ala_ala.psf')
<CharmmPsfFile 33 atoms; 3 residues; 32 bonds; NOT parametrized>

Load a PDB and CIF file

>>> load_file('4lzt.pdb')
<Structure 1164 atoms; 274 residues; 0 bonds; PBC (triclinic); NOT parametrized>
>>> load_file('4LZT.cif')
<Structure 1164 atoms; 274 residues; 0 bonds; PBC (triclinic); NOT parametrized>

Load a Gromacs topology file – only works with Gromacs installed

>>> load_file('1aki.ff99sbildn.top')
<GromacsTopologyFile 40560 atoms [9650 EPs]; 9779 residues; 30934 bonds; parametrized>

Load a SDF file – only works with rdkit installed

>>> load_file('mol.sdf', structure=True)
<Structure 34 atoms; 1 residues; 66 bonds; NOT parametrized>
parmed.load_rdkit(rmol)

Load a Mol object and return a populated Structure instance

Parameters
rmol: :class:`Mol`

RDKit Mol object to convert

Examples

>>> from rdkit import Chem
>>> import parmed as pmd
>>> mol = Chem.MolFromSmiles('Cc1ccccc1')
>>> struct = pmd.load_rdkit(mol)
parmed.load_rosetta(pose)

Load a Pose object and return a populated Structure instance

Parameters
posePose

PyRosetta Pose object to convert

parmed.read_CIF(filename, skip_bonds=False)

Read a PDBx or mmCIF file and return a populated Structure class

Parameters
filenamestr or file-like

Name of PDB file to read, or a file-like object that can iterate over the lines of a PDB. Compressed file names can be specified and are determined by file-name extension (e.g., file.pdb.gz, file.pdb.bz2)

skip_bondsbool, optional

If True, skip trying to assign bonds. This can save substantial time when parsing large files with non-standard residue names. However, no bonds are assigned. This is OK if, for instance, the CIF file is being parsed simply for its coordinates. Default is False.

Returns
structure1 [, structure2 [, structure3 [, …] ] ]
structure#Structure

The Structure object initialized with all of the information from the PDBx/mmCIF file. No bonds or other topological features are added by default. If multiple structures are defined in the CIF file, multiple Structure instances will be returned as a tuple.

Raises
ValueError if the file severely violates the PDB format specification.
If this occurs, check the formatting on each line and make sure it
matches the others.
parmed.read_PDB(filename, skip_bonds=False)

Read a PDB file and return a populated Structure class

Parameters
filenamestr or file-like

Name of the PDB file to read, or a file-like object that can iterate over the lines of a PDB. Compressed file names can be specified and are determined by file-name extension (e.g., file.pdb.gz, file.pdb.bz2)

skip_bondsbool, optional

If True, skip trying to assign bonds. This can save substantial time when parsing large files with non-standard residue names. However, no bonds are assigned. This is OK if, for instance, the PDB file is being parsed simply for its coordinates. This may also reduce element assignment if element information is not present in the PDB file already. Default is False.

Returns
structureStructure

The Structure object initialized with all of the information from the PDB file. No bonds or other topological features are added by default.

parmed.write_CIF(struct, dest, renumber=True, coordinates=None, altlocs='all', write_anisou=False, standard_resnames=False)

Write a PDB file from the current Structure instance

Parameters
structStructure

The structure from which to write the PDBx/mmCIF file

deststr or file-like

Either a file name or a file-like object containing a write method to which to write the PDB file. If it is a filename that ends with .gz or .bz2, a compressed version will be written using either gzip or bzip2, respectively.

renumberbool

If True, renumber the atoms and residues sequentially as they are stored in the structure. If False, use the original numbering if it was assigned previously

coordinatesarray-like of float

If provided, these coordinates will be written to the PDB file instead of the coordinates stored in the structure. These coordinates should line up with the atom order in the structure (not necessarily the order of the “original” PDB file if they differ)

altlocsstr

Keyword controlling which alternate locations are printed to the resulting PDB file. Allowable options are:

  • ‘all’ : (default) print all alternate locations

  • ‘first’ : print only the first alternate locations

  • ‘occupancy’ : print the one with the largest occupancy. If two conformers have the same occupancy, the first one to occur is printed

Input is case-insensitive, and partial strings are permitted as long as it is a substring of one of the above options that uniquely identifies the choice.

write_anisoubool

If True, an ANISOU record is written for every atom that has one. If False, ANISOU records are not written

standard_resnamesbool, optional

If True, common aliases for various amino and nucleic acid residues will be converted into the PDB-standard values. Default is False

Notes

If multiple coordinate frames are present, these will be written as separate models (but only the unit cell from the first model will be written, as the PDBx standard dictates that only one set of unit cells shall be present).

parmed.write_PDB(struct, dest, renumber=True, coordinates=None, altlocs='all', write_anisou=False, charmm=False, use_hetatoms=True, standard_resnames=False, increase_tercount=True, write_links=False)

Write a PDB file from a Structure instance

Parameters
structStructure

The structure from which to write the PDB file

deststr or file-like

Either a file name or a file-like object containing a write method to which to write the PDB file. If it is a filename that ends with .gz or .bz2, a compressed version will be written using either gzip or bzip2, respectively.

renumberbool, optional, default True

If True, renumber the atoms and residues sequentially as they are stored in the structure. If False, use the original numbering if it was assigned previously.

coordinatesarray-like of float, optional

If provided, these coordinates will be written to the PDB file instead of the coordinates stored in the structure. These coordinates should line up with the atom order in the structure (not necessarily the order of the “original” PDB file if they differ)

altlocsstr, optional, default ‘all’

Keyword controlling which alternate locations are printed to the resulting PDB file. Allowable options are:

  • ‘all’ : print all alternate locations

  • ‘first’ : print only the first alternate locations

  • ‘occupancy’ : print the one with the largest occupancy. If two conformers have the same occupancy, the first one to occur is printed

Input is case-insensitive, and partial strings are permitted as long as it is a substring of one of the above options that uniquely identifies the choice.

write_anisoubool, optional, default False

If True, an ANISOU record is written for every atom that has one. If False, ANISOU records are not written.

charmmbool, optional, default False

If True, SEGID will be written in columns 73 to 76 of the PDB file in the typical CHARMM-style PDB output. This will be omitted for any atom that does not contain a SEGID identifier.

use_hetatoms: bool, optional, default True

If True, certain atoms will have the HETATM tag instead of ATOM as per the PDB-standard.

standard_resnamesbool, optional, default False

If True, common aliases for various amino and nucleic acid residues will be converted into the PDB-standard values.

increase_tercountbool, optional, default True

If True, the TER atom number field increased by one compared to atom card preceding it; this conforms to PDB standard.

write_linksbool, optional, default False

If True, any LINK records stored in the Structure will be written to the LINK records near the top of the PDB file. If this is True, then renumber must be False or a ValueError will be thrown

Notes

If multiple coordinate frames are present, these will be written as separate models (but only the unit cell from the first model will be written, as the PDB standard dictates that only one set of unit cells shall be present).