BioSimSpace._SireWrappers.System#

class BioSimSpace._SireWrappers.System(system)[source]#

A container class for storing molecular systems.

__init__(system)[source]#

Constructor.

Parameters:

system (Sire.System.System, System, Sire.Mol._Mol.Molecule, Molecule, Molecules [Molecule]) – A Sire or BioSimSpace System object, a Sire or BioSimSpace Molecule object, a BioSimSpace Molecules object, or a list of BioSimSpace molecule objects.

Methods

__init__(system)

Constructor.

addMolecules(molecules)

Add a molecule, or list of molecules to the system.

charge([property_map, is_lambda1])

Return the total molecular charge.

copy()

Return a copy of this object.

fileFormat([property_map])

Return the file formats associated with the system.

getAminoAcids([property_map])

Return a list containing all of the amino acid residues in the system.

getAtom(index)

Return the atom specified by the absolute index.

getAtoms()

Return a list containing the atoms in the system.

getAxisAlignedBoundingBox([property_map])

Get the axis-aligned bounding box enclosing the object.

getBox([property_map])

Get the size of the periodic simulation box.

getIndex(item)

Convert indices of atoms and residues to their absolute values in this system.

getMolecule(index)

Return the molecule at the given index.

getMolecules([group])

Return a list containing all of the molecules in the specified group.

getNucleotides([property_map])

Return a list containing all of the nucleotide residues in the system.

getPerturbableMolecules()

Return a list containing all of the perturbable molecules in the system.

getResidue(index)

Return the residue specified by the absolute index.

getResidues()

Return a list containing the residues in the system.

getRestraintAtoms(restraint[, mol_index, ...])

Get the indices of atoms involved in a restraint.

getWaterMolecules([property_map])

Return a list containing all of the water molecules in the system.

isSame(other[, excluded_properties, ...])

Check whether "other" is the same as this system.

makeWhole([property_map])

Make all molecules in the system "whole", i.e. unwrap any molecules that have been split across the periodic boundary.

nAminoAcids()

Return the number of amino acid residues in the system.

nAtoms()

Return the number of atoms in the system.

nChains()

Return the number of chains in the system.

nMolecules()

Return the number of molecules in the system.

nNucleotides()

Return the number of nucleotide residues in the system.

nPerturbableMolecules()

Return the number of perturbable molecules in the system.

nResidues()

Return the number of residues in the system.

nWaterMolecules()

Return the number of water molecules in the system.

reduceBoxVectors([bias, property_map])

Reduce the box vectors.

removeMolecules(molecules)

Remove a molecule, or list of molecules from the system.

removeWaterMolecules()

Remove all of the water molecules from the system.

repartitionHydrogenMass([factor, water, ...])

Redistrubute mass of heavy atoms connected to bonded hydrogens into the hydrogen atoms.

rotateBoxVectors([origin, precision, ...])

Rotate the box vectors of the system so that the first vector is aligned with the x-axis, the second vector lies in the x-y plane, and the third vector has a positive z-component.

save(filebase)

Stream a wrapped Sire object to file.

search(query[, property_map])

Search the system for atoms, residues, and molecules.

setBox(box[, angles, property_map])

Set the size of the periodic simulation box.

smiles()

Return the SMILES string representation of this object.

translate(vector[, property_map])

Translate the system.

updateMolecule(index, molecule)

Update the molecule at the given index.

updateMolecules(molecules)

Update a molecule, or list of molecules in the system.

addMolecules(molecules)[source]#

Add a molecule, or list of molecules to the system.

Parameters:

molecules (Molecule, Molecules, [Molecule], System, SearchResult) – A Molecule, Molecules object, a list of Molecule objects, a System, or a SearchResult containing molecules.

charge(property_map={}, is_lambda1=False)[source]#

Return the total molecular charge.

Parameters:
  • property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

  • is_lambda1 (bool) – Whether to use the charge at lambda = 1 if the molecule is merged.

Returns:

charge – The molecular charge.

Return type:

Charge

copy()#

Return a copy of this object. The return type is same as the object on which copy is called.

Returns:

object – A copy of the object.

Return type:

Atom, Residue, Molecule, Molecules, System

fileFormat(property_map={})[source]#

Return the file formats associated with the system.

Parameters:

property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

Returns:

format – The file formats associated with the system.

Return type:

str

getAminoAcids(property_map={})[source]#

Return a list containing all of the amino acid residues in the system.

Parameters:

property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

Returns:

residues – The list of all amino acid residues in the system.

Return type:

[Residue]

getAtom(index)[source]#

Return the atom specified by the absolute index.

Parameters:

index (int) – The absolute index of the atom in the system.

Returns:

atom – The atom at the specified index.

Return type:

Atom

getAtoms()[source]#

Return a list containing the atoms in the system.

Returns:

atoms – The list of atoms in the system.

Return type:

[Atom]

getAxisAlignedBoundingBox(property_map={})#

Get the axis-aligned bounding box enclosing the object.

Parameters:

property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

Returns:

  • box_min ([Length]) – The minimum coordinates of the axis-aligned bounding box in each dimension.

  • box_max ([Length]) – The minimum coordinates of the axis-aligned bounding box in each dimension.

getBox(property_map={})[source]#

Get the size of the periodic simulation box.

Parameters:

property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

Returns:

  • box_size ([Length]) – The box vector magnitudes in each dimension.

  • angles ([Angle]) – The box vector angles: yz, xz, and xy.

getIndex(item)[source]#

Convert indices of atoms and residues to their absolute values in this system.

Parameters:

item (Atom, Residue Molecule) – An Atom, Residue, or Molecule object from the System, or a list containing objects of these types.

Returns:

index – The absolute index of the atom/residue/molecule in the system.

Return type:

int, [int]

getMolecule(index)[source]#

Return the molecule at the given index.

Parameters:

index (int) – The index of the molecule.

Returns:

molecule – The requested molecule.

Return type:

Molecule

getMolecules(group='all')[source]#

Return a list containing all of the molecules in the specified group.

Parameters:

group (str) – The name of the molecule group.

Returns:

molecules – The list of molecules in the group.

Return type:

[Molecule]

getNucleotides(property_map={})[source]#

Return a list containing all of the nucleotide residues in the system.

Parameters:

property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

Returns:

residues – The list of all nucleotide residues in the system.

Return type:

[Residue]

getPerturbableMolecules()[source]#

Return a list containing all of the perturbable molecules in the system.

Returns:

molecules – A container of perturbable molecule objects. The container will be empty if no perturbable molecules are present.

Return type:

Molecules

getResidue(index)[source]#

Return the residue specified by the absolute index.

Parameters:

index (int) – The absolute index of the residue in the system.

Returns:

residue – The residue at the specified index.

Return type:

Residue

getResidues()[source]#

Return a list containing the residues in the system.

Returns:

residues – The list of residues in the system.

Return type:

[Residue]

getRestraintAtoms(restraint, mol_index=None, is_absolute=True, allow_zero_matches=False, property_map={})[source]#

Get the indices of atoms involved in a restraint.

Parameters:
  • restraint (str) – The type of restraint.

  • mol_index (int) – The index of the molecule of interest. If None, then the entire system is searched.

  • is_absolute (bool) – Whether the indices are absolute, i.e. indexed within the entire system. If False, then indices are relative to the molecule in which the atom is found.

  • allow_zero_matches (bool) – Whether to raise an exception when no atoms match the restraint. Setting to false is useful if you need to determine restraints on a per-molecule basis, i.e. so molecules won’t contain atoms matching the restraint, but the entire system will.

  • property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their

Returns:

indices – A list of the backbone atom indices.

Return type:

[int]

getWaterMolecules(property_map={})[source]#

Return a list containing all of the water molecules in the system.

Parameters:

property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

Returns:

molecules – A container of water molecule objects. The container will be empty if no water molecules are present.

Return type:

Molecules

isSame(other, excluded_properties=[], property_map0={}, property_map1={}, skip_water=True)[source]#

Check whether “other” is the same as this system.

Parameters:
  • other (System) – Another BioSimSpace system.

  • excluded_properties ([str]) – A list of properties to exclude when comparing systems. This allows you to check whether a subset of the system is the same, e.g. checking whether the topology is the same, but allowing different coordinates.

  • property_map0 (dict) – A dictionary that maps “properties” in this system to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

  • property_map1 (dict) – A dictionary that maps “properties” in “other” to their user defined values.

  • skip_water (bool) – Whether to skip water molecules when comparing systems.

Returns:

is_same – Whether the systems are the same.

Return type:

bool

makeWhole(property_map={})[source]#

Make all molecules in the system “whole”, i.e. unwrap any molecules that have been split across the periodic boundary.

Parameters:

property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

nAminoAcids()[source]#

Return the number of amino acid residues in the system.

Returns:

num_aa – The number of amino acid residues in the system.

Return type:

int

nAtoms()[source]#

Return the number of atoms in the system.

Returns:

num_atoms – The number of atoms in the system.

Return type:

int

nChains()[source]#

Return the number of chains in the system.

Returns:

num_chains – The number of chains in the system.

Return type:

int

nMolecules()[source]#

Return the number of molecules in the system.

Returns:

num_molecules – The number of molecules in the system.

Return type:

int

nNucleotides()[source]#

Return the number of nucleotide residues in the system.

Returns:

num_nuc – The number of nucleotide residues in the system.

Return type:

int

nPerturbableMolecules()[source]#

Return the number of perturbable molecules in the system.

Returns:

num_perturbable – The number of perturbable molecules in the system.

Return type:

int

nResidues()[source]#

Return the number of residues in the system.

Returns:

num_residues – The number of residues in the system.

Return type:

int

nWaterMolecules()[source]#

Return the number of water molecules in the system.

Returns:

num_waters – The number of water molecules in the system.

Return type:

int

reduceBoxVectors(bias=0, property_map={})[source]#

Reduce the box vectors.

Parameters:
  • bias (float) – The bias to use when rounding during the lattice reduction. Negative values biases towards left-tilting boxes, whereas positive values biases towards right-tilting boxes. This can be used to ensure that rounding is performed in a consistent direction, avoiding oscillation when the box is instantiated from box vectors, or dimensions and angles, that have been read from fixed-precision input files.

  • property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

removeMolecules(molecules)[source]#

Remove a molecule, or list of molecules from the system.

Parameters:

molecules (Molecule, Molecules, [Molecule]) – A Molecule, Molecules object, or list of Molecule objects.

removeWaterMolecules()[source]#

Remove all of the water molecules from the system.

repartitionHydrogenMass(factor=4, water='no', use_coordinates=False, property_map={})[source]#

Redistrubute mass of heavy atoms connected to bonded hydrogens into the hydrogen atoms. This allows the use of larger simulation integration time steps without encountering instabilities related to high-frequency hydrogen motion.

Parameters:
  • factor (float) – The repartioning scale factor. Hydrogen masses are scaled by this amount.

  • water (str) – Whether to repartition masses for water molecules. Options are “yes”, “no”, and “exclusive”, which can be used to repartition masses for water molecules only.

  • use_coordinates (bool) – Whether to use the current molecular coordinates to work out the connectivity before repartitioning. If False, the information from the molecular topology, e.g. force field, will be used, if present.

  • property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

rotateBoxVectors(origin=(0.0000e+00 A, 0.0000e+00 A, 0.0000e+00 A), precision=0.0, property_map={})[source]#

Rotate the box vectors of the system so that the first vector is aligned with the x-axis, the second vector lies in the x-y plane, and the third vector has a positive z-component. This is a requirement of certain molecular dynamics engines and is used for simulation efficiency. All vector properties of the system, e.g. coordinates, velocities, etc., will be rotated accordingly.

Parameters:
  • origin (Coordinate) – The origin of the triclinic space. This is the point about which the lattice vectors are rotated.

  • precision (float) – The precision to use when sorting box vectors by their magnitude. This can be used to prevent unecessary box rotation when the vectors were read from a text file with limited precision.

  • property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

save(filebase)#

Stream a wrapped Sire object to file.

Parameters:
  • sire_object (System) – A wrapped Sire object.

  • filebase (str) – The base name of the binary output file.

search(query, property_map={})[source]#

Search the system for atoms, residues, and molecules. Search results will be reduced to their minimal representation, i.e. a molecule containing a single residue will be returned as a residue.

Parameters:
  • query (str) – The search query.

  • property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

Returns:

results – A list of objects matching the search query.

Return type:

[Atom, Residue, Molecule, …]

Examples

Search for residues names ALA.

>>> result = system.search("resname ALA")

Search for residues names ALA using a case insensitive search.

>>> result = system.search("resname /ala/i")

Search for all carbon atoms in residues named ALA.

>>> result = system.search("resname ALA and atomname C")

Search for all oxygen or hydrogen atoms.

>>> result = system.search("element oxygen or element hydrogen")

Search for atom index 23 in molecule index 10.

>>> result = system.search("molidx 10 and atomidx 23")
setBox(box, angles=[90.0000 degrees, 90.0000 degrees, 90.0000 degrees], property_map={})[source]#

Set the size of the periodic simulation box.

Parameters:
  • box ([Length]) – The box vector magnitudes in each dimension.

  • angles ([Angle]) – The box vector angles: yz, xz, and xy.

  • property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

smiles()#

Return the SMILES string representation of this object.

Returns:

smiles – The SMILES string representation of the object(s).

Return type:

str, [str]

translate(vector, property_map={})[source]#

Translate the system.

Parameters:
  • vector ([Length]) – The translation vector in Angstroms.

  • property_map (dict) – A dictionary that maps system “properties” to their user defined values. This allows the user to refer to properties with their own naming scheme, e.g. { “charge” : “my-charge” }

updateMolecule(index, molecule)[source]#

Update the molecule at the given index.

Parameters:
  • index (int) – The index of the molecule.

  • molecule (Molecule) – The updated (or replacement) molecule.

updateMolecules(molecules)[source]#

Update a molecule, or list of molecules in the system.

Parameters:

molecules (Molecule, [Molecule]) – A Molecule, or list of Molecule objects.