Update topology/system particle bookkeeping in OpenMM #1051
Conversation
|
Check out this pull request on See visual diffs & provide feedback on Jupyter Notebooks. Powered by ReviewNB |
|
|
||
| # make a new "residue" for each molecule which has virtual sites | ||
| if len(virtual_sites_in_this_molecule) > 0: | ||
| this_virtual_site_residue = openmm_topology.addResidue("VS", virtual_site_chain) |
There was a problem hiding this comment.
Would it make sense to name these residues by their parent residue name? The use case I'm thinking of here is let's say you want to select all your waters, you'd go "resname WAT HOH SOL" and know you're picking up all the waters and their virtual sites. Similarly if you knew you had a ligand named LIG, then your user experience (especially when you output to a PDB from an openmm topology via something like mdtraj) is a lot cleaner.
|
|
||
| # make a new "residue" for each molecule which has virtual sites | ||
| if len(virtual_sites_in_this_molecule) > 0: | ||
| this_virtual_site_residue = openmm_topology.addResidue("VS", virtual_site_chain) |
There was a problem hiding this comment.
If we're not picking up the resid from anywhere, I believe that doesn't need to be unique - so it might make sense to set the resid to of the virtual site residue to the same value of its parent residue?
There was a problem hiding this comment.
Sure - with a21f2ce,
In [1]: from openff.toolkit import *
In [2]: mol = Molecule.from_smiles("O")
In [3]: for atom in mol.atoms:
...: atom.metadata["residue_name"] = "SOL"
In [4]: In [30]: x = (
...: ...: ForceField("opc.offxml")
...: ...: .create_interchange(
...: ...: Topology.from_molecules(
...: ...: [
...: ...: mol,#Molecule.from_smiles("O"),
...: ...: mol,#Molecule.from_smiles("O"),
...: mol,
...: ...: ],
...: ...: ),
...: ...: )
...: ...: .to_openmm_topology()
...: ...: )
/Users/mattthompson/micromamba/envs/new-models/lib/python3.11/site-packages/smirnoff99frosst/smirnoff99frosst.py:11: UserWarning: Module openff was already imported from None, but /Users/mattthompson/software/openff-interchange is being added to sys.path
from pkg_resources import resource_filename
/Users/mattthompson/micromamba/envs/new-models/lib/python3.11/site-packages/mdtraj/formats/__init__.py:13: DeprecationWarning: 'xdrlib' is deprecated and slated for removal in Python 3.13
from mdtraj.formats.trr import TRRTrajectoryFile
In [5]: for residue in x.residues():
...: print(residue.index, [*residue.atoms()])
0 [<Atom 0 (O1x) of chain 0 residue 0 (SOL)>, <Atom 1 (H1x) of chain 0 residue 0 (SOL)>, <Atom 2 (H2x) of chain 0 residue 0 (SOL)>]
1 [<Atom 3 (O1x) of chain 1 residue 1 (SOL)>, <Atom 4 (H1x) of chain 1 residue 1 (SOL)>, <Atom 5 (H2x) of chain 1 residue 1 (SOL)>]
2 [<Atom 6 (O1x) of chain 2 residue 2 (SOL)>, <Atom 7 (H1x) of chain 2 residue 2 (SOL)>, <Atom 8 (H2x) of chain 2 residue 2 (SOL)>]
3 [<Atom 9 (EP) of chain 3 residue 3 (SOL)>]
4 [<Atom 10 (EP) of chain 3 residue 4 (SOL)>]
5 [<Atom 11 (EP) of chain 3 residue 5 (SOL)>]There was a problem hiding this comment.
Ah, yes, IDs and names are different nouns, after all ...
There was a problem hiding this comment.
Nice.
- Having
to_openmm_topologymatchto_openmm_systemseems pretty essential - Maintaining the guarantee that we'll keep atom indices constant no matter what you do is extremely valuable
- This is a great solution
I have one blocking bug fix in to_pdb, a suggestion of how to accomplish Irfan's suggestion of preserving parent atom metadata (which I think is very valuable for analysis), and a few nits.
| @@ -680,6 +697,7 @@ def to_pdb(self, file_path: Path | str, include_virtual_sites: bool = False): | |||
| ) | |||
|
|
|||
| openmm_topology = self.to_openmm_topology( | |||
| collate=True, | |||
| ensure_unique_atom_names=False, | |||
| ) | |||
| positions = get_positions_with_virtual_sites(self) | |||
There was a problem hiding this comment.
| positions = get_positions_with_virtual_sites(self) | |
| positions = get_positions_with_virtual_sites(self, collate=True) |
This was producing some very confusing visualizations for me! Noting that if you decide to go with the new behavior as I suggest above, then this line is correct.
Co-authored-by: Josh A. Mitchell <[email protected]>
Description
Resolves #1049
Checklist