openforcefield.topology.Topology

class openforcefield.topology.Topology(other=None)[source]

A Topology is a chemical representation of a system containing one or more molecules appearing in a specified order.

Warning

This API is experimental and subject to change.

Examples

Import some utilities

>>> from simtk.openmm import app
>>> from openforcefield.tests.utils import get_data_filename, get_packmol_pdbfile
>>> pdb_filepath = get_packmol_pdbfile('cyclohexane_ethanol_0.4_0.6')
>>> monomer_names = ('cyclohexane', 'ethanol')

Create a Topology object from a PDB file and sdf files defining the molecular contents

>>> from openforcefield.topology import Molecule, Topology
>>> pdbfile = app.PDBFile(pdb_filepath)
>>> sdf_filepaths = [get_data_filename(f'systems/monomers/{name}.sdf') for name in monomer_names]
>>> unique_molecules = [Molecule.from_file(sdf_filepath) for sdf_filepath in sdf_filepaths]
>>> topology = Topology.from_openmm(pdbfile.topology, unique_molecules=unique_molecules)

Create a Topology object from a PDB file and IUPAC names of the molecular contents

>>> pdbfile = app.PDBFile(pdb_filepath)
>>> unique_molecules = [Molecule.from_iupac(name) for name in monomer_names]
>>> topology = Topology.from_openmm(pdbfile.topology, unique_molecules=unique_molecules)

Create an empty Topology object and add a few copies of a single benzene molecule

>>> topology = Topology()
>>> molecule = Molecule.from_iupac('benzene')
>>> molecule_topology_indices = [topology.add_molecule(molecule) for index in range(10)]
Attributes:
angles

Iterable of Tuple[TopologyAtom]: iterator over the angles in this Topology.

aromaticity_model

Get the aromaticity model applied to all molecules in the topology.

box_vectors

Return the box vectors of the topology, if specified

charge_model

Get the partial charge model applied to all molecules in the topology.

constrained_atom_pairs

Returns the constrained atom pairs of the Topology

fractional_bond_order_model

Get the fractional bond order model for the Topology.

impropers

Iterable of Tuple[TopologyAtom]: iterator over the improper torsions in this Topology.

n_angles

int: number of angles in this Topology.

n_impropers

int: number of improper torsions in this Topology.

n_propers

int: number of proper torsions in this Topology.

n_reference_molecules

Returns the number of reference (unique) molecules in in this Topology.

n_topology_atoms

Returns the number of topology atoms in in this Topology.

n_topology_bonds

Returns the number of TopologyBonds in in this Topology.

n_topology_molecules

Returns the number of topology molecules in in this Topology.

n_topology_particles

Returns the number of topology particles (TopologyAtoms and TopologyVirtualSites) in in this Topology.

n_topology_virtual_sites

Returns the number of TopologyVirtualSites in in this Topology.

propers

Iterable of Tuple[TopologyAtom]: iterator over the proper torsions in this Topology.

reference_molecules

Get an iterator of reference molecules in this Topology.

topology_atoms

Returns an iterator over the atoms in this Topology.

topology_bonds

Returns an iterator over the bonds in this Topology

topology_molecules

Returns an iterator over all the TopologyMolecules in this Topology

topology_particles

Returns an iterator over the particles (TopologyAtoms and TopologyVirtualSites) in this Topology.

topology_virtual_sites

Get an iterator over the virtual sites in this Topology

Methods

add_constraint(iatom, jatom[, distance]) Mark a pair of atoms as constrained.
add_molecule(molecule[, …]) Add a Molecule to the Topology.
add_particle(particle) Add a Particle to the Topology.
assert_bonded(atom1, atom2) Raise an exception if the specified atoms are not bonded in the topology.
atom(atom_topology_index) Get the TopologyAtom at a given Topology atom index.
bond(bond_topology_index) Get the TopologyBond at a given Topology bond index.
chemical_environment_matches(query[, …]) Retrieve all matches for a given chemical environment query.
from_bson(serialized) Instantiate an object from a BSON serialized representation.
from_dict(d) Static constructor from dictionary representation.
from_json(serialized) Instantiate an object from a JSON serialized representation.
from_mdtraj(mdtraj_topology[, unique_molecules]) Construct an openforcefield Topology object from an MDTraj Topology object.
from_messagepack(serialized) Instantiate an object from a MessagePack serialized representation.
from_molecules(molecules) Create a new Topology object containing one copy of each of the specified molecule(s).
from_openmm(openmm_topology[, unique_molecules]) Construct an openforcefield Topology object from an OpenMM Topology object.
from_parmed(parmed_structure[, unique_molecules])

Warning

This functionality will be implemented in a future toolkit release.

from_pickle(serialized) Instantiate an object from a pickle serialized representation.
from_toml(serialized) Instantiate an object from a TOML serialized representation.
from_xml(serialized) Instantiate an object from an XML serialized representation.
from_yaml(serialized) Instantiate from a YAML serialized representation.
get_bond_between(i, j) Returns the bond between two atoms
get_fractional_bond_order(iatom, jatom) Retrieve the fractional bond order for a bond.
is_bonded(i, j) Returns True if the two atoms are bonded
is_constrained(iatom, jatom) Check if a pair of atoms are marked as constrained.
to_bson() Return a BSON serialized representation.
to_dict() Convert to dictionary representation.
to_json([indent]) Return a JSON serialized representation.
to_messagepack() Return a MessagePack representation.
to_parmed()

Warning

This functionality will be implemented in a future toolkit release.

to_pickle() Return a pickle serialized representation.
to_toml() Return a TOML serialized representation.
to_xml([indent]) Return an XML representation.
to_yaml() Return a YAML serialized representation.
virtual_site(vsite_topology_index) Get the TopologyAtom at a given Topology atom index.
__init__(other=None)[source]

Create a new Topology.

Parameters:
other : optional, default=None

If specified, attempt to construct a copy of the Topology from the specified object. This might be a Topology object, or a file that can be used to construct a Topology object or serialized Topology object.

Methods

__init__([other]) Create a new Topology.
add_constraint(iatom, jatom[, distance]) Mark a pair of atoms as constrained.
add_molecule(molecule[, …]) Add a Molecule to the Topology.
add_particle(particle) Add a Particle to the Topology.
assert_bonded(atom1, atom2) Raise an exception if the specified atoms are not bonded in the topology.
atom(atom_topology_index) Get the TopologyAtom at a given Topology atom index.
bond(bond_topology_index) Get the TopologyBond at a given Topology bond index.
chemical_environment_matches(query[, …]) Retrieve all matches for a given chemical environment query.
from_bson(serialized) Instantiate an object from a BSON serialized representation.
from_dict(d) Static constructor from dictionary representation.
from_json(serialized) Instantiate an object from a JSON serialized representation.
from_mdtraj(mdtraj_topology[, unique_molecules]) Construct an openforcefield Topology object from an MDTraj Topology object.
from_messagepack(serialized) Instantiate an object from a MessagePack serialized representation.
from_molecules(molecules) Create a new Topology object containing one copy of each of the specified molecule(s).
from_openmm(openmm_topology[, unique_molecules]) Construct an openforcefield Topology object from an OpenMM Topology object.
from_parmed(parmed_structure[, unique_molecules])

Warning

This functionality will be implemented in a future toolkit release.

from_pickle(serialized) Instantiate an object from a pickle serialized representation.
from_toml(serialized) Instantiate an object from a TOML serialized representation.
from_xml(serialized) Instantiate an object from an XML serialized representation.
from_yaml(serialized) Instantiate from a YAML serialized representation.
get_bond_between(i, j) Returns the bond between two atoms
get_fractional_bond_order(iatom, jatom) Retrieve the fractional bond order for a bond.
is_bonded(i, j) Returns True if the two atoms are bonded
is_constrained(iatom, jatom) Check if a pair of atoms are marked as constrained.
to_bson() Return a BSON serialized representation.
to_dict() Convert to dictionary representation.
to_json([indent]) Return a JSON serialized representation.
to_messagepack() Return a MessagePack representation.
to_parmed()

Warning

This functionality will be implemented in a future toolkit release.

to_pickle() Return a pickle serialized representation.
to_toml() Return a TOML serialized representation.
to_xml([indent]) Return an XML representation.
to_yaml() Return a YAML serialized representation.
virtual_site(vsite_topology_index) Get the TopologyAtom at a given Topology atom index.

Attributes

angles Iterable of Tuple[TopologyAtom]: iterator over the angles in this Topology.
aromaticity_model Get the aromaticity model applied to all molecules in the topology.
box_vectors Return the box vectors of the topology, if specified Returns ——- box_vectors : simtk.unit.Quantity wrapped numpy array The unit-wrapped box vectors of this topology
charge_model Get the partial charge model applied to all molecules in the topology.
constrained_atom_pairs Returns the constrained atom pairs of the Topology
fractional_bond_order_model Get the fractional bond order model for the Topology.
impropers Iterable of Tuple[TopologyAtom]: iterator over the improper torsions in this Topology.
n_angles int: number of angles in this Topology.
n_impropers int: number of improper torsions in this Topology.
n_propers int: number of proper torsions in this Topology.
n_reference_molecules Returns the number of reference (unique) molecules in in this Topology.
n_topology_atoms Returns the number of topology atoms in in this Topology.
n_topology_bonds Returns the number of TopologyBonds in in this Topology.
n_topology_molecules Returns the number of topology molecules in in this Topology.
n_topology_particles Returns the number of topology particles (TopologyAtoms and TopologyVirtualSites) in in this Topology.
n_topology_virtual_sites Returns the number of TopologyVirtualSites in in this Topology.
propers Iterable of Tuple[TopologyAtom]: iterator over the proper torsions in this Topology.
reference_molecules Get an iterator of reference molecules in this Topology.
topology_atoms Returns an iterator over the atoms in this Topology.
topology_bonds Returns an iterator over the bonds in this Topology
topology_molecules Returns an iterator over all the TopologyMolecules in this Topology
topology_particles Returns an iterator over the particles (TopologyAtoms and TopologyVirtualSites) in this Topology.
topology_virtual_sites Get an iterator over the virtual sites in this Topology
reference_molecules

Get an iterator of reference molecules in this Topology.

Returns:
iterable of openforcefield.topology.Molecule
classmethod from_molecules(molecules)[source]

Create a new Topology object containing one copy of each of the specified molecule(s).

Parameters:
molecules : Molecule or iterable of Molecules

One or more molecules to be added to the Topology

Returns:
topology : Topology

The Topology created from the specified molecule(s)

assert_bonded(atom1, atom2)[source]

Raise an exception if the specified atoms are not bonded in the topology.

Parameters:
atom1, atom2 : openforcefield.topology.Atom or int

The atoms or atom topology indices to check to ensure they are bonded

aromaticity_model

Get the aromaticity model applied to all molecules in the topology.

Returns:
aromaticity_model : str

Aromaticity model in use.

box_vectors

Return the box vectors of the topology, if specified Returns ——- box_vectors : simtk.unit.Quantity wrapped numpy array

The unit-wrapped box vectors of this topology
charge_model

Get the partial charge model applied to all molecules in the topology.

Returns:
charge_model : str

Charge model used for all molecules in the Topology.

constrained_atom_pairs

Returns the constrained atom pairs of the Topology

Returns:
constrained_atom_pairs : dict

dictionary of the form d[(atom1_topology_index, atom2_topology_index)] = distance (float)

fractional_bond_order_model

Get the fractional bond order model for the Topology.

Returns:
fractional_bond_order_model : str

Fractional bond order model in use.

n_reference_molecules

Returns the number of reference (unique) molecules in in this Topology.

Returns:
n_reference_molecules : int
n_topology_molecules

Returns the number of topology molecules in in this Topology.

Returns:
n_topology_molecules : int
topology_molecules

Returns an iterator over all the TopologyMolecules in this Topology

Returns:
topology_molecules : Iterable of TopologyMolecule
n_topology_atoms

Returns the number of topology atoms in in this Topology.

Returns:
n_topology_atoms : int
topology_atoms

Returns an iterator over the atoms in this Topology. These will be in ascending order of topology index (Note that this is not necessarily the same as the reference molecule index)

Returns:
topology_atoms : Iterable of TopologyAtom
n_topology_bonds

Returns the number of TopologyBonds in in this Topology.

Returns:
n_bonds : int
topology_bonds

Returns an iterator over the bonds in this Topology

Returns:
topology_bonds : Iterable of TopologyBond
n_topology_particles

Returns the number of topology particles (TopologyAtoms and TopologyVirtualSites) in in this Topology.

Returns:
n_topology_particles : int
topology_particles

Returns an iterator over the particles (TopologyAtoms and TopologyVirtualSites) in this Topology. The TopologyAtoms will be in order of ascending Topology index (Note that this may differ from the order of atoms in the reference molecule index).

Returns:
topology_particles : Iterable of TopologyAtom and TopologyVirtualSite
n_topology_virtual_sites

Returns the number of TopologyVirtualSites in in this Topology.

Returns:
n_virtual_sites : iterable of TopologyVirtualSites
topology_virtual_sites

Get an iterator over the virtual sites in this Topology

Returns:
topology_virtual_sites : Iterable of TopologyVirtualSite
n_angles

int: number of angles in this Topology.

angles

Iterable of Tuple[TopologyAtom]: iterator over the angles in this Topology.

n_propers

int: number of proper torsions in this Topology.

propers

Iterable of Tuple[TopologyAtom]: iterator over the proper torsions in this Topology.

n_impropers

int: number of improper torsions in this Topology.

impropers

Iterable of Tuple[TopologyAtom]: iterator over the improper torsions in this Topology.

chemical_environment_matches(query, aromaticity_model='MDL', toolkit_registry=<openforcefield.utils.toolkits.ToolkitRegistry object>)[source]

Retrieve all matches for a given chemical environment query.

TODO: * Do we want to generalize this to other kinds of queries too, like mdtraj DSL, pymol selections, atom index slices, etc?

We could just call it topology.matches(query)
Parameters:
query : str or ChemicalEnvironment

SMARTS string (with one or more tagged atoms) or ChemicalEnvironment query Query will internally be resolved to SMARTS using query.as_smarts() if it has an .as_smarts method.

aromaticity_model : str

Override the default aromaticity model for this topology and use the specified aromaticity model instead. Allowed values: [‘MDL’]

Returns:
matches : list of TopologyAtom tuples

A list of all matching Atom tuples

to_dict()[source]

Convert to dictionary representation.

classmethod from_dict(d)[source]

Static constructor from dictionary representation.

classmethod from_openmm(openmm_topology, unique_molecules=None)[source]

Construct an openforcefield Topology object from an OpenMM Topology object.

Parameters:
openmm_topology : simtk.openmm.app.Topology

An OpenMM Topology object

unique_molecules : iterable of objects that can be used to construct unique Molecule objects

All unique molecules must be provided, in any order, though multiple copies of each molecule are allowed. The atomic elements and bond connectivity will be used to match the reference molecules to molecule graphs appearing in the OpenMM Topology. If bond orders are present in the OpenMM Topology, these will be used in matching as well. If all bonds have bond orders assigned in mdtraj_topology, these bond orders will be used to attempt to construct the list of unique Molecules if the unique_molecules argument is omitted.

Returns:
topology : openforcefield.topology.Topology

An openforcefield Topology object

static from_mdtraj(mdtraj_topology, unique_molecules=None)[source]

Construct an openforcefield Topology object from an MDTraj Topology object.

Parameters:
mdtraj_topology : mdtraj.Topology

An MDTraj Topology object

unique_molecules : iterable of objects that can be used to construct unique Molecule objects

All unique molecules mult be provided, in any order, though multiple copies of each molecule are allowed. The atomic elements and bond connectivity will be used to match the reference molecules to molecule graphs appearing in the MDTraj Topology. If bond orders are present in the MDTraj Topology, these will be used in matching as well. If all bonds have bond orders assigned in mdtraj_topology, these bond orders will be used to attempt to construct the list of unique Molecules if the unique_molecules argument is omitted.

Returns:
topology : openforcefield.Topology

An openforcefield Topology object

static from_parmed(parmed_structure, unique_molecules=None)[source]

Warning

This functionality will be implemented in a future toolkit release.

Construct an openforcefield Topology object from a ParmEd Structure object.

Parameters:
parmed_structure : parmed.Structure

A ParmEd structure object

unique_molecules : iterable of objects that can be used to construct unique Molecule objects

All unique molecules mult be provided, in any order, though multiple copies of each molecule are allowed. The atomic elements and bond connectivity will be used to match the reference molecules to molecule graphs appearing in the structure’s topology object. If bond orders are present in the structure’s topology object, these will be used in matching as well. If all bonds have bond orders assigned in the structure’s topology object, these bond orders will be used to attempt to construct the list of unique Molecules if the unique_molecules argument is omitted.

Returns:
topology : openforcefield.Topology

An openforcefield Topology object

to_parmed()[source]

Warning

This functionality will be implemented in a future toolkit release.

Create a ParmEd Structure object.

Returns:
parmed_structure : parmed.Structure

A ParmEd Structure objecft

get_bond_between(i, j)[source]

Returns the bond between two atoms

Parameters:
i, j : int or TopologyAtom

Atoms or atom indices to check

Returns:
bond : TopologyBond

The bond between i and j.

is_bonded(i, j)[source]

Returns True if the two atoms are bonded

Parameters:
i, j : int or TopologyAtom

Atoms or atom indices to check

Returns:
is_bonded : bool

True if atoms are bonded, False otherwise.

atom(atom_topology_index)[source]

Get the TopologyAtom at a given Topology atom index.

Parameters:
atom_topology_index : int

The index of the TopologyAtom in this Topology

Returns:
An openforcefield.topology.TopologyAtom
virtual_site(vsite_topology_index)[source]

Get the TopologyAtom at a given Topology atom index.

Parameters:
vsite_topology_index : int

The index of the TopologyVirtualSite in this Topology

Returns:
An openforcefield.topology.TopologyVirtualSite
bond(bond_topology_index)[source]

Get the TopologyBond at a given Topology bond index.

Parameters:
bond_topology_index : int

The index of the TopologyBond in this Topology

Returns:
An openforcefield.topology.TopologyBond
add_particle(particle)[source]

Add a Particle to the Topology.

Parameters:
particle : Particle

The Particle to be added. The Topology will take ownership of the Particle.

add_molecule(molecule, local_topology_to_reference_index=None)[source]

Add a Molecule to the Topology.

Parameters:
molecule : Molecule

The Molecule to be added.

local_topology_to_reference_index: dict, optional, default = None

Dictionary of {TopologyMolecule_atom_index : Molecule_atom_index} for the TopologyMolecule that will be built

Returns:
index : int

The index of this molecule in the topology

classmethod from_bson(serialized)

Instantiate an object from a BSON serialized representation.

Specification: http://bsonspec.org/

Parameters:
serialized : bytes

A BSON serialized representation of the object

Returns:
instance : cls

An instantiated object

classmethod from_json(serialized)

Instantiate an object from a JSON serialized representation.

Specification: https://www.json.org/

Parameters:
serialized : str

A JSON serialized representation of the object

Returns:
instance : cls

An instantiated object

classmethod from_messagepack(serialized)

Instantiate an object from a MessagePack serialized representation.

Specification: https://msgpack.org/index.html

Parameters:
serialized : bytes

A MessagePack-encoded bytes serialized representation

Returns:
instance : cls

Instantiated object.

classmethod from_pickle(serialized)

Instantiate an object from a pickle serialized representation.

Warning

This is not recommended for safe, stable storage since the pickle specification may change between Python versions.

Parameters:
serialized : str

A pickled representation of the object

Returns:
instance : cls

An instantiated object

classmethod from_toml(serialized)

Instantiate an object from a TOML serialized representation.

Specification: https://github.com/toml-lang/toml

Parameters:
serlialized : str

A TOML serialized representation of the object

Returns:
instance : cls

An instantiated object

classmethod from_xml(serialized)

Instantiate an object from an XML serialized representation.

Specification: https://www.w3.org/XML/

Parameters:
serialized : bytes

An XML serialized representation

Returns:
instance : cls

Instantiated object.

classmethod from_yaml(serialized)

Instantiate from a YAML serialized representation.

Specification: http://yaml.org/

Parameters:
serialized : str

A YAML serialized representation of the object

Returns:
instance : cls

Instantiated object

to_bson()

Return a BSON serialized representation.

Specification: http://bsonspec.org/

Returns:
serialized : bytes

A BSON serialized representation of the objecft

to_json(indent=None)

Return a JSON serialized representation.

Specification: https://www.json.org/

Parameters:
indent : int, optional, default=None

If not None, will pretty-print with specified number of spaces for indentation

Returns:
serialized : str

A JSON serialized representation of the object

to_messagepack()

Return a MessagePack representation.

Specification: https://msgpack.org/index.html

Returns:
serialized : bytes

A MessagePack-encoded bytes serialized representation of the object

to_pickle()

Return a pickle serialized representation.

Warning

This is not recommended for safe, stable storage since the pickle specification may change between Python versions.

Returns:
serialized : str

A pickled representation of the object

to_toml()

Return a TOML serialized representation.

Specification: https://github.com/toml-lang/toml

Returns:
serialized : str

A TOML serialized representation of the object

to_xml(indent=2)

Return an XML representation.

Specification: https://www.w3.org/XML/

Parameters:
indent : int, optional, default=2

If not None, will pretty-print with specified number of spaces for indentation

Returns:
serialized : bytes

A MessagePack-encoded bytes serialized representation.

to_yaml()

Return a YAML serialized representation.

Specification: http://yaml.org/

Returns:
serialized : str

A YAML serialized representation of the object

add_constraint(iatom, jatom, distance=True)[source]

Mark a pair of atoms as constrained.

Constraints between atoms that are not bonded (e.g., rigid waters) are permissible.

Parameters:
iatom, jatom : Atom

Atoms to mark as constrained These atoms may be bonded or not in the Topology

distance : simtk.unit.Quantity, optional, default=True

Constraint distance True if distance has yet to be determined False if constraint is to be removed

is_constrained(iatom, jatom)[source]

Check if a pair of atoms are marked as constrained.

Parameters:
iatom, jatom : int

Indices of atoms to mark as constrained.

Returns:
distance : simtk.unit.Quantity or bool

True if constrained but constraints have not yet been applied Distance if constraint has already been added to System

get_fractional_bond_order(iatom, jatom)[source]

Retrieve the fractional bond order for a bond.

An Exception is raised if it cannot be determined.

Parameters:
iatom, jatom : Atom

Atoms for which a fractional bond order is to be retrieved.

Returns:
order : float

Fractional bond order between the two specified atoms.