openff.toolkit.topology.Topology
- class openff.toolkit.topology.Topology(other=None)[source]
A Topology is a chemical representation of a system containing one or more molecules appearing in a specified order.
As of the 0.7.0 release, the Topology particle indexing system puts all atoms before all virtualsites. This ensures that atoms keep the same Topology particle index value, even if the Topology is modified during system creation by the addition of virtual sites.
Warning
This API is experimental and subject to change.
Examples
Import some utilities
>>> from openmm import app >>> from openff.toolkit.tests.utils import get_data_file_path, get_packmol_pdb_file_path >>> pdb_filepath = get_packmol_pdb_file_path('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 openff.toolkit.topology import Molecule, Topology >>> pdbfile = app.PDBFile(pdb_filepath) >>> sdf_filepaths = [get_data_file_path(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)]
- __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 OpenFF 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 OpenFF Topology object from an OpenMM Topology object.
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
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.
nth_degree_neighbors
(n_degrees)Return canonicalized pairs of atoms whose shortest separation is exactly n bonds.
to_bson
()Return a BSON serialized representation.
to_dict
()Convert to dictionary representation.
to_file
(filename, positions[, file_format, ...])Save coordinates and topology to a PDB file.
to_json
([indent])Return a JSON serialized representation.
Return a MessagePack representation.
to_openmm
([ensure_unique_atom_names])Create an OpenMM Topology object.
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
Iterate over improper torsions in the molecule, but only those with trivalent centers, reporting the central atom first in each improper.
iterator over the angles in this Topology.
Get the aromaticity model applied to all molecules in the topology.
Return the box vectors of the topology, if specified
Get the partial charge model applied to all molecules in the topology.
Returns the constrained atom pairs of the Topology
Get the fractional bond order model for the Topology.
iterator over the possible improper torsions in this Topology.
Return whether or not this Topology is intended to be described with periodic boundary conditions.
number of angles in this Topology.
number of possible improper torsions in this Topology.
number of proper torsions in this Topology.
Returns the number of reference (unique) molecules in in this Topology.
Returns the number of topology atoms in in this Topology.
Returns the number of TopologyBonds in in this Topology.
Returns the number of topology molecules in in this Topology.
Returns the number of topology particles (TopologyAtoms and TopologyVirtualSites) in in this Topology.
Returns the number of TopologyVirtualSites in in this Topology.
iterator over the proper torsions in this Topology.
Get an iterator of reference molecules in this Topology.
Iterate over improper torsions in the molecule, but only those with trivalent centers, reporting the central atom second in each improper.
Returns an iterator over the atoms in this Topology.
Returns an iterator over the bonds in this Topology
Returns an iterator over all the TopologyMolecules in this Topology
Returns an iterator over the particles (TopologyAtoms and TopologyVirtualSites) in this Topology.
Get an iterator over the virtual sites in this Topology
- property reference_molecules
Get an iterator of reference molecules in this Topology.
- Returns
iterable of openff.toolkit.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 (openff.toolkit.topology.Atom or int) – The atoms or atom topology indices to check to ensure they are bonded
atom2 (openff.toolkit.topology.Atom or int) – The atoms or atom topology indices to check to ensure they are bonded
- property aromaticity_model
Get the aromaticity model applied to all molecules in the topology.
- Returns
aromaticity_model (str) – Aromaticity model in use.
- property box_vectors
Return the box vectors of the topology, if specified
- Returns
box_vectors (openmm.unit.Quantity wrapped numpy array of shape (3, 3)) – The unit-wrapped box vectors of this topology
- property is_periodic
Return whether or not this Topology is intended to be described with periodic boundary conditions.
- property 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.
- property 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)
- property 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.
- property n_reference_molecules
Returns the number of reference (unique) molecules in in this Topology.
- Returns
n_reference_molecules (int)
- property n_topology_molecules
Returns the number of topology molecules in in this Topology.
- Returns
n_topology_molecules (int)
- property topology_molecules
Returns an iterator over all the TopologyMolecules in this Topology
- Returns
topology_molecules (Iterable of TopologyMolecule)
- property n_topology_atoms
Returns the number of topology atoms in in this Topology.
- Returns
n_topology_atoms (int)
- property 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)
- property n_topology_bonds
Returns the number of TopologyBonds in in this Topology.
- Returns
n_bonds (int)
- property topology_bonds
Returns an iterator over the bonds in this Topology
- Returns
topology_bonds (Iterable of TopologyBond)
- property n_topology_particles
Returns the number of topology particles (TopologyAtoms and TopologyVirtualSites) in in this Topology.
- Returns
n_topology_particles (int)
- property 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)
- property n_topology_virtual_sites
Returns the number of TopologyVirtualSites in in this Topology.
- Returns
n_virtual_sites (iterable of TopologyVirtualSites)
- property topology_virtual_sites
Get an iterator over the virtual sites in this Topology
- Returns
topology_virtual_sites (Iterable of TopologyVirtualSite)
- property angles
iterator over the angles in this Topology.
- Type
Iterable of Tuple[TopologyAtom]
- property propers
iterator over the proper torsions in this Topology.
- Type
Iterable of Tuple[TopologyAtom]
- property impropers
iterator over the possible improper torsions in this Topology.
- Type
Iterable of Tuple[TopologyAtom]
- property smirnoff_impropers
Iterate over improper torsions in the molecule, but only those with trivalent centers, reporting the central atom second in each improper.
Note that it’s possible that a trivalent center will not have an improper assigned. This will depend on the force field that is used.
Also note that this will return 6 possible atom orderings around each improper center. In current SMIRNOFF parameterization, three of these six orderings will be used for the actual assignment of the improper term and measurement of the angles. These three orderings capture the three unique angles that could be calculated around the improper center, therefore the sum of these three terms will always return a consistent energy.
For more details on the use of three-fold (‘trefoil’) impropers, see https://openforcefield.github.io/standards/standards/smirnoff/#impropertorsions
- Returns
impropers (set of tuple) – An iterator of tuples, each containing the indices of atoms making up a possible improper torsion. The central atom is listed second in each tuple.
See also
- property amber_impropers
Iterate over improper torsions in the molecule, but only those with trivalent centers, reporting the central atom first in each improper.
Note that it’s possible that a trivalent center will not have an improper assigned. This will depend on the force field that is used.
Also note that this will return 6 possible atom orderings around each improper center. In current AMBER parameterization, one of these six orderings will be used for the actual assignment of the improper term and measurement of the angle. This method does not encode the logic to determine which of the six orderings AMBER would use.
- Returns
impropers (set of tuple) – An iterator of tuples, each containing the indices of atoms making up a possible improper torsion. The central atom is listed first in each tuple.
See also
- nth_degree_neighbors(n_degrees: int)[source]
Return canonicalized pairs of atoms whose shortest separation is exactly n bonds. Only pairs with increasing atom indices are returned.
- Parameters
n (int) – The number of bonds separating atoms in each pair
- Returns
neighbors (iterator of tuple of TopologyAtom) – Tuples (len 2) of atom that are separated by
n
bonds.
Notes
The criteria used here relies on minimum distances; when there are multiple valid paths between atoms, such as atoms in rings, the shortest path is considered. For example, two atoms in “meta” positions with respect to each other in a benzene are separated by two paths, one length 2 bonds and the other length 4 bonds. This function would consider them to be 2 apart and would not include them if
n=4
was passed.
- chemical_environment_matches(query, aromaticity_model='MDL', unique=False, toolkit_registry=GLOBAL_TOOLKIT_REGISTRY)[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 usingquery.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 Topology._ChemicalEnvironmentMatch) – A list of tuples, containing the topology indices of the matching atoms.
- 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 OpenFF Topology object from an OpenMM Topology object.
- Parameters
openmm_topology (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 OpenMMTopology
, these will be used in matching as well.
- Returns
topology (openff.toolkit.topology.Topology) – An OpenFF Topology object
- to_openmm(ensure_unique_atom_names=True)[source]
Create an OpenMM Topology object.
The OpenMM
Topology
object will have one residue per topology molecule. Currently, the number of chains depends on how many copies of the same molecule are in theTopology
. Molecules with more than 5 copies are all assigned to a single chain, otherwise one chain is created for each molecule. This behavior may change in the future.- Parameters
ensure_unique_atom_names (bool, optional. Default=True) – Whether to check that the molecules in each molecule have unique atom names, and regenerate them if not. Note that this looks only at molecules, and does not guarantee uniqueness in the entire Topology.
- Returns
openmm_topology (openmm.app.Topology) – An OpenMM Topology object
- to_file(filename, positions, file_format='PDB', keepIds=False)[source]
Save coordinates and topology to a PDB file.
Reference: https://github.com/openforcefield/openff-toolkit/issues/502
Notes:
This doesn’t handle virtual sites (they’re ignored)
Atom numbering may not remain same, for example if the atoms in water are numbered as 1001, 1002, 1003, they would change to 1, 2, 3. This doesn’t affect the topology or coordinates or atom-ordering in any way.
Same issue with the amino acid names in the pdb file, they are not returned.
- Parameters
filename (str) – name of the pdb file to write to
positions (n_atoms x 3 numpy array or openmm.unit.Quantity-wrapped n_atoms x 3 iterable) – Can be an openmm ‘quantity’ object which has atomic positions as a list of Vec3s along with associated units, otherwise a 3D array of UNITLESS numbers are considered as “Angstroms” by default
file_format (str) – Output file format. Case insensitive. Currently only supported value is “pdb”.
- static from_mdtraj(mdtraj_topology, unique_molecules=None)[source]
Construct an OpenFF 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 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 MDTraj
Topology
. If bond orders are present in the MDTrajTopology
, these will be used in matching as well.
- Returns
topology (openff.toolkit.topology.Topology) – An OpenFF Topology object
- get_bond_between(i, j)[source]
Returns the bond between two atoms
- Parameters
i (int or TopologyAtom) – Atoms or atom indices to check
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 (int or TopologyAtom) – Atoms or atom indices to check
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 openff.toolkit.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 openff.toolkit.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 openff.toolkit.topology.TopologyBond
- 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_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. You can optionally request that the atoms be added to the Topology in a different order than they appear in the Molecule.
- 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. If None, this function will add the atoms to the Topology in the order that they appear in the reference molecule.
- Returns
index (int) – The index of this molecule in the topology
- 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 (Atom) – Atoms to mark as constrained These atoms may be bonded or not in the Topology
jatom (Atom) – Atoms to mark as constrained These atoms may be bonded or not in the Topology
distance (openmm.unit.Quantity, optional, default=True) – Constraint distance
True
if distance has yet to be determinedFalse
if constraint is to be removed
- is_constrained(iatom, jatom)[source]
Check if a pair of atoms are marked as constrained.