Skip to content

Commit

Permalink
Switch to MoleculeGroup container to preserve ordering. [ref #193]
Browse files Browse the repository at this point in the history
  • Loading branch information
lohedges committed Apr 29, 2021
1 parent 2f79dc2 commit da1e68d
Show file tree
Hide file tree
Showing 4 changed files with 56 additions and 84 deletions.
4 changes: 2 additions & 2 deletions python/BioSimSpace/Solvent/_solvent.py
Original file line number Diff line number Diff line change
Expand Up @@ -789,7 +789,7 @@ def _solvate(molecule, box, angles, shell, model, num_point,
if molecule is not None:
if type(molecule) is _System:
# Extract the non-water molecules from the original system.
non_waters = _Molecules(molecule.search("not water")._sire_object.toMolecules())
non_waters = _Molecules(molecule.search("not water")._sire_object.toGroup())

# Create a system by adding these to the water molecules from
# gmx solvate, which will include the original waters.
Expand Down Expand Up @@ -987,7 +987,7 @@ def _solvate(molecule, box, angles, shell, model, num_point,
if molecule is not None:
if type(molecule) is _System:
# Extract the non-water molecules from the original system.
non_waters = _Molecules(molecule.search("not water")._sire_object.toMolecules())
non_waters = _Molecules(molecule.search("not water")._sire_object.toGroup())

# Create a system by adding these to the water and ion
# molecules from gmx solvate, which will include the
Expand Down
28 changes: 1 addition & 27 deletions python/BioSimSpace/_SireWrappers/_molecule.py
Original file line number Diff line number Diff line change
Expand Up @@ -131,11 +131,6 @@ def __add__(self, other):

# A single Molecule object.
elif type(other) is Molecule:
mol = self.copy()
mol._renumber(_SireMol.MolNum.getUniqueNumber())
molecules[0] = mol
other = other.copy()
other._renumber(_SireMol.MolNum.getUniqueNumber())
molecules.append(other)

# A Molecules object.
Expand All @@ -144,12 +139,7 @@ def __add__(self, other):

# A list of Molecule objects.
elif type(other) is list and all(isinstance(x, Molecule) for x in other):
mol = self.copy()
mol._renumber(_SireMol.MolNum.getUniqueNumber())
molecules[0] = mol
for mol in other:
mol = mol.copy()
mol._renumber(_SireMol.MolNum.getUniqueNumber())
molecules.append(mol)

# Unsupported.
Expand All @@ -170,7 +160,7 @@ def __add__(self, other):
for molecule in molecules:
molgrp.add(molecule._sire_object)

return _Molecules(molgrp.molecules())
return _Molecules(molgrp)

def number(self):
"""Return the number of the molecule. Each molecule has a unique
Expand Down Expand Up @@ -1122,22 +1112,6 @@ def makeCompatibleWith(self, molecule, property_map={}, overwrite=True,
# Finally, commit the changes to the internal object.
self._sire_object = edit_mol.commit()

def _renumber(self, mol_num):
"""Renumber the molecule with a unique MolNum.
Parameters
----------
mol_num : Sire.Mol.MolNum
The molecule number.
"""
if type(mol_num) is not _SireMol.MolNum:
raise TypeError("'mol_num' must be of type 'Sire.Mol.MolNum'")

edit_mol = self._sire_object.edit()
edit_mol.renumber(mol_num)
self._sire_object = edit_mol.commit()

def _getPropertyMap0(self):
"""Generate a property map for the lambda = 0 state of the merged molecule."""

Expand Down
32 changes: 12 additions & 20 deletions python/BioSimSpace/_SireWrappers/_molecules.py
Original file line number Diff line number Diff line change
Expand Up @@ -20,8 +20,8 @@
#####################################################################

"""
A thin wrapper around Sire.Mol.Molecules. This is an internal package and should
not be directly exposed to the user.
A thin wrapper around Sire.Mol.MoleculeGroup. This is an internal package and
should not be directly exposed to the user.
"""

__author__ = "Lester Hedges"
Expand All @@ -48,7 +48,7 @@ def __init__(self, molecules):
Parameters
----------
molecules : Sire.Mol.Molecules, Sire.System.System, \
molecules : Sire.Mol.MoleculeGroup, Sire.System.System, \
:class:`System <BioSimSpace._SireWrappers.System>` \
[:class:`Molecule <BioSimSpace._SireWrappers.Molecule>`]
A Sire Molecues object, a Sire or BioSimSpace System object,
Expand All @@ -62,7 +62,7 @@ def __init__(self, molecules):
molecules = list(molecules)

# A Sire Molecules object.
if type(molecules) is _SireMol.Molecules:
if type(molecules) is _SireMol.MoleculeGroup:
super().__init__(molecules)

# A Sire System object.
Expand All @@ -71,27 +71,21 @@ def __init__(self, molecules):

# A BioSimSpace System object.
elif type(molecules) is _System:
super().__init__(molecules._sire_object.molecules())
super().__init__(molecules._sire_object.group(_SireMol.MGName("all")))

# Another BioSimSpace Molecule object.
elif type(molecules) is list and all(isinstance(x, _Molecule) for x in molecules):
mols = _SireMol.Molecules()
molgrp = _SireMol.MoleculeGroup("all")
for molecule in molecules:
mols.add(molecule._sire_object)
super().__init__(mols)
molgrp.add(molecule._sire_object)
super().__init__(molgrp)

# Invalid type.
else:
raise TypeError("'molecules' must be of type 'Sire.Mol.Molecules' "
raise TypeError("'molecules' must be of type 'Sire.Mol.MoleculeGroup' "
"'Sire.System.System', 'BioSimSpace._SireWrappers.System', "
"or a list of 'BioSimSpace._SireWrappers.Molecule' types.")

# Store the list of MolNums.
self._mol_nums = self._sire_object.molNums()

# Sort the molecule numbers.
self._mol_nums.sort()

# Store the number of molecules.
self._num_mols = self._sire_object.nMolecules()

Expand Down Expand Up @@ -124,8 +118,6 @@ def __add__(self, other):

# A Molecule object.
elif type(other) is _Molecule:
mol = other.copy()
mol._renumber(_SireMol.MolNum.getUniqueNumber())
molecules.add(mol._sire_object)

# A Molecules object.
Expand All @@ -135,8 +127,8 @@ def __add__(self, other):
# A list of Molecule objects.
elif type(other) is list and all(isinstance(x, _Molecule) for x in other):
for molecule in other:
mol = molecule.copy()
mol._renumber(_SireMol.MolNum.getUniqueNumber())
#mol = molecule.copy()
#mol._renumber(_SireMol.MolNum.getUniqueNumber())
molecules.add(mol._sire_object)

# Unsupported.
Expand Down Expand Up @@ -178,7 +170,7 @@ def __getitem__(self, key):
key = key + self._num_mols

# Extract and return the corresponding molecule.
return _Molecule(self._sire_object.molecule(self._mol_nums[key]))
return _Molecule(self._sire_object[_SireMol.MolIdx(key)])

def __setitem__(self, key, value):
"""Set a molecule in the container."""
Expand Down
76 changes: 41 additions & 35 deletions python/BioSimSpace/_SireWrappers/_system.py
Original file line number Diff line number Diff line change
Expand Up @@ -42,10 +42,6 @@

from ._sire_wrapper import SireWrapper as _SireWrapper

class _MolWithResName(_SireMol.MolWithResID):
def __init__(self, resname):
super().__init__(_SireMol.ResName(resname))

class System(_SireWrapper):
"""A container class for storing molecular systems."""

Expand Down Expand Up @@ -121,14 +117,8 @@ def __init__(self, system):
# Initialise dictionary to map MolNum to MolIdx.
self._molecule_index = {}

# Store the sorted molecule numbers.
# Store the molecule numbers.
self._mol_nums = self._sire_object.molNums()
self._mol_nums.sort()

# Store a mapping from MolNum to index.
self._mol_num_to_idx = {}
for idx, num in enumerate(self._mol_nums):
self._mol_num_to_idx[num] = idx

# Intialise the iterator counter.
self._iter_count = 0
Expand Down Expand Up @@ -198,7 +188,7 @@ def __getitem__(self, key):
key = key + self.nMolecules()

# Extract and return the corresponding molecule.
return _Molecule(self._sire_object.molecule(self._mol_nums[key]))
return _Molecule(self._sire_object[_SireMol.MolIdx(key)])

def __setitem__(self, key, value):
"""Set a molecule in the container."""
Expand Down Expand Up @@ -417,7 +407,6 @@ def addMolecules(self, molecules):

# Update the molecule numbers.
self._mol_nums = self._sire_object.molNums()
self._mol_nums.sort()

def removeMolecules(self, molecules):
"""Remove a molecule, or list of molecules from the system.
Expand Down Expand Up @@ -467,7 +456,6 @@ def removeMolecules(self, molecules):

# Update the molecule numbers.
self._mol_nums = self._sire_object.molNums()
self._mol_nums.sort()

def removeWaterMolecules(self):
"""Remove all of the water molecules from the system."""
Expand All @@ -483,7 +471,6 @@ def removeWaterMolecules(self):

# Update the molecule numbers.
self._mol_nums = self._sire_object.molNums()
self._mol_nums.sort()

def updateMolecule(self, index, molecule):
"""Updated the molecule at the given index.
Expand All @@ -498,6 +485,9 @@ def updateMolecule(self, index, molecule):
The updated (or replacement) molecule.
"""

if type(index) is not int:
raise TypeError("'index' must be of type 'int'")

if index < -self.nMolecules() or index >= self.nMolecules():
raise IndexError("The molecule 'index' is out of range.")

Expand Down Expand Up @@ -533,7 +523,6 @@ def updateMolecule(self, index, molecule):

# Update the molecule numbers.
self._mol_nums = self._sire_object.molNums()
self._mol_nums.sort()

def updateMolecules(self, molecules):
"""Update a molecule, or list of molecules in the system.
Expand Down Expand Up @@ -579,7 +568,6 @@ def updateMolecules(self, molecules):

# Update the molecule numbers.
self._mol_nums = self._sire_object.molNums()
self._mol_nums.sort()

def getMolecule(self, index):
"""Return the molecule at the given index.
Expand Down Expand Up @@ -624,7 +612,7 @@ def getMolecules(self, group="all"):
raise ValueError("No molecules in group '%s'" % group)

# Return a molecules container.
return _Molecules(molgrp.molecules())
return _Molecules(molgrp)

def getWaterMolecules(self):
"""Return a list containing all of the water molecules in the system.
Expand All @@ -636,7 +624,7 @@ def getWaterMolecules(self):
A container of water molecule objects.
"""

return _Molecules(self._sire_object.search("water").toMolecules())
return _Molecules(self._sire_object.search("water").toGroup())

def nWaterMolecules(self):
"""Return the number of water molecules in the system.
Expand Down Expand Up @@ -1360,8 +1348,12 @@ def _updateCoordinates(self, system, property_map0={}, property_map1={},
# Try to update the coordinates property.
try:
mol0 = mol0.edit().setProperty(prop, mol1.property(prop1)).molecule().commit()
except:
raise _IncompatibleError("Unable to update 'coordinates' for molecule index '%d'" % idx)
except Exception as e:
msg = "Unable to update 'coordinates' for molecule index '%d'" % idx
if _isVerbose():
raise _IncompatibleError(msg) from e
else:
raise _IncompatibleError(msg) from None

# Update the molecule in the original system.
self._sire_object.update(mol0)
Expand Down Expand Up @@ -1444,7 +1436,7 @@ def _getRelativeIndices(self, abs_index):
mol_num_last = mol_num
mol_num = num
else:
mol_idx = self._mol_num_to_idx[mol_num_last]
mol_idx = self._molecule_index[mol_num_last]
return mol_idx, abs_index - tally_last

raise ValueError("'abs_index' exceeded system atom tally!")
Expand Down Expand Up @@ -1533,26 +1525,40 @@ def _set_water_topology(self, format):
else:
new_waters = _SireIO.setGromacsWater(self._sire_object.search("water"), water_model)

# Loop over all of the renamed water molecules, delete the old one
# from the system, then add the renamed one back in.
# TODO: This is a hack since the "update" method of Sire.System
# doesn't work properly at present.
self.removeWaterMolecules()
for idx, wat in enumerate(new_waters):
# Renumber the new molecule to match the original.
edit_mol = wat.edit()
edit_mol.renumber(waters[idx]._sire_object.number())
wat = edit_mol.commit()
# Create a new system and molecule group.
system = _SireSystem.System("BioSimSpace System")
molgrp = _SireMol.MoleculeGroup("all")

# Create a dictionary of water molecule numbers to list index.
wat_num_to_idx = {}
for idx, wat in enumerate(waters):
wat_num_to_idx[wat.number()] = idx

# Add all moleclules from the original system into the new one,
# preserving the order.
for mol in self.getMolecules():
try:
# Water molecule.
molgrp.add(new_waters[wat_num_to_idx[mol.number()]])
except:
# Non-water molecule.
molgrp.add(mol._sire_object)

# Add the group to the system.
system.add(molgrp)

# Copy across system properties.
for prop in self._sire_object.propertyKeys():
system.setProperty(prop, self._sire_object.property(prop))

# Re-add the water to the system.
self._sire_object.add(wat, _SireMol.MGName("all"))
# Set the system.
self._sire_object = system

# Reset the index mappings.
self._reset_mappings()

# Update the molecule numbers.
self._mol_nums = self._sire_object.molNums()
self._mol_nums.sort()

# Import at bottom of module to avoid circular dependency.
from ._atom import Atom as _Atom
Expand Down

0 comments on commit da1e68d

Please sign in to comment.