From f50a9b555a4f7f5ce70db2c1a99ba434ebad405a Mon Sep 17 00:00:00 2001 From: chrisjonesBSU Date: Wed, 18 Feb 2026 16:23:11 +0000 Subject: [PATCH 1/6] Add slicing methods, set number as option in iter_sites methods --- gmso/core/topology.py | 36 +++++++++++++++---- gmso/tests/test_topology.py | 11 ++++++ gmso/tests/test_utils.py | 38 ++++++++++++++++++++ gmso/utils/slicing.py | 71 +++++++++++++++++++++++++++++++++++++ 4 files changed, 150 insertions(+), 6 deletions(-) create mode 100644 gmso/utils/slicing.py diff --git a/gmso/core/topology.py b/gmso/core/topology.py index ddfec08ff..02b91e060 100644 --- a/gmso/core/topology.py +++ b/gmso/core/topology.py @@ -1447,7 +1447,7 @@ def iter_sites(self, key, value): if getattr(site, key) == value: yield site - def iter_sites_by_residue(self, residue_tag): + def iter_sites_by_residue(self, residue_tag, residue_number=None): """Iterate through this topology's sites which contain this specific residue name. See Also @@ -1457,12 +1457,24 @@ def iter_sites_by_residue(self, residue_tag): """ if isinstance(residue_tag, str): for site in self._sites: - if site.residue and getattr(site, "residue").name == residue_tag: - yield site + if residue_tag and residue_number is None: + if site.molecule and getattr(site, "residue").name == residue_tag: + yield site + elif residue_number and residue_tag is None: + if site.molecule and getattr(site, "residue").number == residue_number: + yield site + else: + if all( + [ + site.molecule and getattr(site, "residue").name == residue_tag, + site.molecule and getattr(site, "residue").number == residue_number + ] + ): + yield site else: return self.iter_sites("residue", residue_tag) - def iter_sites_by_molecule(self, molecule_tag): + def iter_sites_by_molecule(self, molecule_tag, molecule_number=None): """Iterate through this topology's sites which contain this specific molecule name. See Also @@ -1472,8 +1484,20 @@ def iter_sites_by_molecule(self, molecule_tag): """ if isinstance(molecule_tag, str): for site in self._sites: - if site.molecule and getattr(site, "molecule").name == molecule_tag: - yield site + if molecule_tag and molecule_number is None: + if site.molecule and getattr(site, "molecule").name == molecule_tag: + yield site + elif molecule_number and molecule_tag is None: + if site.molecule and getattr(site, "molecule").number == molecule_number: + yield site + else: + if all( + [ + site.molecule and getattr(site, "molecule").name == molecule_tag, + site.molecule and getattr(site, "molecule").number == molecule_number + ] + ): + yield site else: return self.iter_sites("molecule", molecule_tag) diff --git a/gmso/tests/test_topology.py b/gmso/tests/test_topology.py index bde0de76f..fbd0a1b62 100644 --- a/gmso/tests/test_topology.py +++ b/gmso/tests/test_topology.py @@ -886,6 +886,17 @@ def test_iter_sites_by_molecule(self, labeled_top): for site in labeled_top.iter_sites_by_molecule(molecule_name): assert site.molecule.name == molecule_name + def test_iter_sites_by_molecule_tag_and_number(self, labeled_top): + molecules = labeled_top.unique_site_labels("molecule", name_only=False) + for molecule in molecules: + for site in labeled_top.iter_sites_by_molecule(molecule): + assert site.residue == molecule + + molecule_names = labeled_top.unique_site_labels("molecule", name_only=True) + for molecule_name in molecule_names: + for site in labeled_top.iter_sites_by_molecule(molecule_name): + assert site.molecule.name == molecule_name + @pytest.mark.parametrize( "connections", ["bonds", "angles", "dihedrals", "impropers"], diff --git a/gmso/tests/test_utils.py b/gmso/tests/test_utils.py index c23b7ef0a..c72f8907e 100644 --- a/gmso/tests/test_utils.py +++ b/gmso/tests/test_utils.py @@ -1,3 +1,4 @@ +import mbuild as mb import numpy as np import pytest import unyt as u @@ -8,6 +9,7 @@ from gmso.utils.io import run_from_ipython from gmso.utils.misc import unyt_to_hashable from gmso.utils.sorting import sort_connection_members, sort_connection_strings +from gmso.utils.slicing import slice_topology_by_molecule, slice_topology_by_residue def test_unyt_to_hashable(): @@ -81,3 +83,39 @@ def test_moment_of_inertia(): masses=np.array([1.0 for i in xyz]), ) assert np.array_equal(tensor, np.array([1, 1, 2])) + + +def test_slice_by_molecule(): + benzene = mb.load("c1ccccc1", smiles=True) + benzene.name = "Benzene" + ethane = mb.load("CC", smiles=True) + ethane.name = "Ethane" + + system = mb.fill_box(compound=[benzene, ethane], n_compounds=[2, 2], box=[2,2,2]) + topology = system.to_gmso() + topology.identify_connections() + + single_benzene_top = slice_topology_by_molecule(topology, "Benzene", 0) + assert single_benzene_top.n_sites == 12 + + all_benzene_top = slice_topology_by_molecule(topology, "Benzene") + assert all_benzene_top.n_sites == 24 + + +def test_slice_by_residue(): + benzene = mb.load("c1ccccc1", smiles=True) + benzene.name = "Benzene" + ethane = mb.load("CC", smiles=True) + ethane.name = "Ethane" + + system = mb.fill_box(compound=[benzene, ethane], n_compounds=[2, 2], box=[2,2,2]) + topology = system.to_gmso() + topology.identify_connections() + + single_ethane_top = slice_topology_by_residue(topology, "Ethane", 0) + assert single_ethane_top.n_sites == 8 + + all_ethane_top = slice_topology_by_residue(topology, "Ethane") + assert all_ethane_top.n_sites == 16 + + diff --git a/gmso/utils/slicing.py b/gmso/utils/slicing.py new file mode 100644 index 000000000..7e8181774 --- /dev/null +++ b/gmso/utils/slicing.py @@ -0,0 +1,71 @@ +from gmso.core.topology import Topology + + +def slice_topology_by_molecule(topology, molecule_tag, molecule_number=None): + """Create a Topology that contains a subset of molecules from another Topology. + + Parameters + ---------- + topology : gmso.core.topology.Topology + The gmso Topology to perform the slice on + molecule_tag : str + The name of the gmso.abstract_site.Molecule object to include in the slice + molecule_number : int, default None + If given, only include a single molecule's sites + If None, then all sites in every molecule matching `molecule_tag` are included + in the sliced topology. + + Returns + ------- + gmso.core.topology.Topology + A new Topology instance containing only sites and connections from matching molecules. + """ + sites = [s for s in topology.iter_sites_by_molecule(molecule_tag=molecule_tag, molecule_number=molecule_number)] + return slice_by_sites(topology=topology, sites=sites) + + +def slice_topology_by_residue(topology, residue_tag, residue_number=None): + """Create a Topology that contains a subset of residues from another Topology. + + Parameters + ---------- + topology : gmso.core.topology.Topology + The gmso Topology to perform the slice on. + residue_tag : str + The name of the gmso.abstract_site.Residue object to include in the slice. + residue_number : int, default None + If given, only include a single residue's sites + If None, then all sites in every residue matching `residue_tag` are included + in the sliced topology. + + Returns + ------- + gmso.core.topology.Topology + A new Topology instance containing only sites and connections from matching molecules. + """ + sites = [s for s in topology.iter_sites_by_residue(residue_tag=residue_tag, residue_number=residue_number)] + return slice_by_sites(topology=topology, sites=sites) + + +def slice_by_sites(topology, sites): + """Used by slice_topology_by_molecule() and slice_topology_by_residue() + + Parameters + ---------- + sites : list of gmso.core.atom.Atom + List of sites to include in the sub-topology. + topology : gmso.core.topology.Topology + The topology being sliced. + """ + new_topology = Topology() + + connections = set() + for site in sites: + new_topology.add_site(site) + for connection in topology.iter_connections_by_site(site): + connections.add(connection) + + for connection in connections: + new_topology.add_connection(connection) + + return new_topology From 284c88b474839be0aa1de8dba1aeda5f27d31ceb Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Wed, 18 Feb 2026 16:38:34 +0000 Subject: [PATCH 2/6] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- gmso/core/topology.py | 22 ++++++++++++++++------ gmso/tests/test_utils.py | 17 +++++++++-------- gmso/utils/slicing.py | 20 +++++++++++++++----- 3 files changed, 40 insertions(+), 19 deletions(-) diff --git a/gmso/core/topology.py b/gmso/core/topology.py index 02b91e060..1a6101774 100644 --- a/gmso/core/topology.py +++ b/gmso/core/topology.py @@ -1461,13 +1461,18 @@ def iter_sites_by_residue(self, residue_tag, residue_number=None): if site.molecule and getattr(site, "residue").name == residue_tag: yield site elif residue_number and residue_tag is None: - if site.molecule and getattr(site, "residue").number == residue_number: + if ( + site.molecule + and getattr(site, "residue").number == residue_number + ): yield site else: if all( [ - site.molecule and getattr(site, "residue").name == residue_tag, - site.molecule and getattr(site, "residue").number == residue_number + site.molecule + and getattr(site, "residue").name == residue_tag, + site.molecule + and getattr(site, "residue").number == residue_number, ] ): yield site @@ -1488,13 +1493,18 @@ def iter_sites_by_molecule(self, molecule_tag, molecule_number=None): if site.molecule and getattr(site, "molecule").name == molecule_tag: yield site elif molecule_number and molecule_tag is None: - if site.molecule and getattr(site, "molecule").number == molecule_number: + if ( + site.molecule + and getattr(site, "molecule").number == molecule_number + ): yield site else: if all( [ - site.molecule and getattr(site, "molecule").name == molecule_tag, - site.molecule and getattr(site, "molecule").number == molecule_number + site.molecule + and getattr(site, "molecule").name == molecule_tag, + site.molecule + and getattr(site, "molecule").number == molecule_number, ] ): yield site diff --git a/gmso/tests/test_utils.py b/gmso/tests/test_utils.py index c72f8907e..e2e5c5e9a 100644 --- a/gmso/tests/test_utils.py +++ b/gmso/tests/test_utils.py @@ -8,8 +8,11 @@ from gmso.utils.geometry import moment_of_inertia from gmso.utils.io import run_from_ipython from gmso.utils.misc import unyt_to_hashable +from gmso.utils.slicing import ( + slice_topology_by_molecule, + slice_topology_by_residue, +) from gmso.utils.sorting import sort_connection_members, sort_connection_strings -from gmso.utils.slicing import slice_topology_by_molecule, slice_topology_by_residue def test_unyt_to_hashable(): @@ -91,7 +94,7 @@ def test_slice_by_molecule(): ethane = mb.load("CC", smiles=True) ethane.name = "Ethane" - system = mb.fill_box(compound=[benzene, ethane], n_compounds=[2, 2], box=[2,2,2]) + system = mb.fill_box(compound=[benzene, ethane], n_compounds=[2, 2], box=[2, 2, 2]) topology = system.to_gmso() topology.identify_connections() @@ -99,7 +102,7 @@ def test_slice_by_molecule(): assert single_benzene_top.n_sites == 12 all_benzene_top = slice_topology_by_molecule(topology, "Benzene") - assert all_benzene_top.n_sites == 24 + assert all_benzene_top.n_sites == 24 def test_slice_by_residue(): @@ -108,14 +111,12 @@ def test_slice_by_residue(): ethane = mb.load("CC", smiles=True) ethane.name = "Ethane" - system = mb.fill_box(compound=[benzene, ethane], n_compounds=[2, 2], box=[2,2,2]) + system = mb.fill_box(compound=[benzene, ethane], n_compounds=[2, 2], box=[2, 2, 2]) topology = system.to_gmso() topology.identify_connections() single_ethane_top = slice_topology_by_residue(topology, "Ethane", 0) - assert single_ethane_top.n_sites == 8 + assert single_ethane_top.n_sites == 8 all_ethane_top = slice_topology_by_residue(topology, "Ethane") - assert all_ethane_top.n_sites == 16 - - + assert all_ethane_top.n_sites == 16 diff --git a/gmso/utils/slicing.py b/gmso/utils/slicing.py index 7e8181774..166b948d2 100644 --- a/gmso/utils/slicing.py +++ b/gmso/utils/slicing.py @@ -11,7 +11,7 @@ def slice_topology_by_molecule(topology, molecule_tag, molecule_number=None): molecule_tag : str The name of the gmso.abstract_site.Molecule object to include in the slice molecule_number : int, default None - If given, only include a single molecule's sites + If given, only include a single molecule's sites If None, then all sites in every molecule matching `molecule_tag` are included in the sliced topology. @@ -20,7 +20,12 @@ def slice_topology_by_molecule(topology, molecule_tag, molecule_number=None): gmso.core.topology.Topology A new Topology instance containing only sites and connections from matching molecules. """ - sites = [s for s in topology.iter_sites_by_molecule(molecule_tag=molecule_tag, molecule_number=molecule_number)] + sites = [ + s + for s in topology.iter_sites_by_molecule( + molecule_tag=molecule_tag, molecule_number=molecule_number + ) + ] return slice_by_sites(topology=topology, sites=sites) @@ -34,8 +39,8 @@ def slice_topology_by_residue(topology, residue_tag, residue_number=None): residue_tag : str The name of the gmso.abstract_site.Residue object to include in the slice. residue_number : int, default None - If given, only include a single residue's sites - If None, then all sites in every residue matching `residue_tag` are included + If given, only include a single residue's sites + If None, then all sites in every residue matching `residue_tag` are included in the sliced topology. Returns @@ -43,7 +48,12 @@ def slice_topology_by_residue(topology, residue_tag, residue_number=None): gmso.core.topology.Topology A new Topology instance containing only sites and connections from matching molecules. """ - sites = [s for s in topology.iter_sites_by_residue(residue_tag=residue_tag, residue_number=residue_number)] + sites = [ + s + for s in topology.iter_sites_by_residue( + residue_tag=residue_tag, residue_number=residue_number + ) + ] return slice_by_sites(topology=topology, sites=sites) From 0d02c4aadf08b027ca2869e89984b8cdd6c8d671 Mon Sep 17 00:00:00 2001 From: chrisjonesBSU Date: Mon, 9 Mar 2026 10:15:37 +0000 Subject: [PATCH 3/6] Add method to convert topology to OpenFF Molecule object --- gmso/external/convert_openmm.py | 42 +++++++++++++++++++++++++++++++++ 1 file changed, 42 insertions(+) diff --git a/gmso/external/convert_openmm.py b/gmso/external/convert_openmm.py index b3ea139a6..83b275220 100644 --- a/gmso/external/convert_openmm.py +++ b/gmso/external/convert_openmm.py @@ -7,12 +7,54 @@ from gmso.utils.io import has_openmm, has_openmm_unit, import_ +from openff.toolkit import Molecule as openFF_Molecule + if has_openmm & has_openmm_unit: openmm_unit = import_("openmm.unit") from openmm import * from openmm.app import * + +def to_openff_molecule(topology): + """Converts a GMSO topology into an OpenFF Molecule.""" + mol = openFF_Molecule() + # Bond order not available at the site level + # But, it can be set during openFF's add_atom' + # Build up a list of sites first to track their bond orders later + sites = [] + for site in topology.sites: + sites.append(site) + + # Index-to-index mapping with sites list + is_aromatic = [False for i in sites] + bond_pairs = [] + bond_orders = [] + + for bond in topology.bonds: + bond_pair = [] + bond_orders.append(bond.bond_order) + for site in bond.connection_members: + site_index = sites.index(site) + if bond.bond_order == 1.5: + is_aromatic[site_index] = True + bond_pair.append(site_index) + bond_pairs.append(bond_pair) + # Now that we have parsed all sites for aromaticity, add each to Molecule + for site, aromatic in zip(sites, is_aromatic): + mol.add_atom(atomic_number=site.element.atomic_number, is_aromatic=aromatic, formal_charge=0) + + for bond, border in zip(bond_pairs, bond_orders): + mol.add_bond( + atom1=bond[0], + atom2=bond[1], + bond_order=int(border), + fractional_bond_order=border, + is_aromatic=True if border==1.5 else False + ) + return mol + + def to_openmm(topology, openmm_object="topology"): """Convert an untyped topology object to an untyped OpenMM modeller or topology. From 176bb784f1ba0a5d728cf99e862f446916bee280 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Mon, 9 Mar 2026 10:15:47 +0000 Subject: [PATCH 4/6] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- gmso/external/convert_openmm.py | 20 +++++++++++--------- 1 file changed, 11 insertions(+), 9 deletions(-) diff --git a/gmso/external/convert_openmm.py b/gmso/external/convert_openmm.py index 83b275220..b1a5c65c9 100644 --- a/gmso/external/convert_openmm.py +++ b/gmso/external/convert_openmm.py @@ -4,18 +4,16 @@ """Convert to and from an OpenMM Topology or System object.""" import unyt as u +from openff.toolkit import Molecule as openFF_Molecule from gmso.utils.io import has_openmm, has_openmm_unit, import_ -from openff.toolkit import Molecule as openFF_Molecule - if has_openmm & has_openmm_unit: openmm_unit = import_("openmm.unit") from openmm import * from openmm.app import * - def to_openff_molecule(topology): """Converts a GMSO topology into an OpenFF Molecule.""" mol = openFF_Molecule() @@ -25,12 +23,12 @@ def to_openff_molecule(topology): sites = [] for site in topology.sites: sites.append(site) - + # Index-to-index mapping with sites list is_aromatic = [False for i in sites] bond_pairs = [] bond_orders = [] - + for bond in topology.bonds: bond_pair = [] bond_orders.append(bond.bond_order) @@ -40,17 +38,21 @@ def to_openff_molecule(topology): is_aromatic[site_index] = True bond_pair.append(site_index) bond_pairs.append(bond_pair) - # Now that we have parsed all sites for aromaticity, add each to Molecule + # Now that we have parsed all sites for aromaticity, add each to Molecule for site, aromatic in zip(sites, is_aromatic): - mol.add_atom(atomic_number=site.element.atomic_number, is_aromatic=aromatic, formal_charge=0) - + mol.add_atom( + atomic_number=site.element.atomic_number, + is_aromatic=aromatic, + formal_charge=0, + ) + for bond, border in zip(bond_pairs, bond_orders): mol.add_bond( atom1=bond[0], atom2=bond[1], bond_order=int(border), fractional_bond_order=border, - is_aromatic=True if border==1.5 else False + is_aromatic=True if border == 1.5 else False, ) return mol From 7f2f9e74a61b9ed29fee3d8d41a3e6affe51c625 Mon Sep 17 00:00:00 2001 From: chrisjonesBSU Date: Thu, 30 Apr 2026 11:26:44 +0100 Subject: [PATCH 5/6] start charges.py file to calc partial charges --- gmso/core/topology.py | 17 +++++++++++++++++ gmso/utils/charges.py | 19 +++++++++++++++++++ gmso/utils/slicing.py | 5 ++--- 3 files changed, 38 insertions(+), 3 deletions(-) create mode 100644 gmso/utils/charges.py diff --git a/gmso/core/topology.py b/gmso/core/topology.py index 1a6101774..dda386e47 100644 --- a/gmso/core/topology.py +++ b/gmso/core/topology.py @@ -1251,6 +1251,23 @@ def get_index(self, member): return index + def calculate_charges(self, method="am1bcc", slice_by="Molecule"): + print("Not implemented yet. We are working on it.") + # Get all of the moleucles types in the topology + molecules = set() + for site in self.sites: + molecules.add(site.molecule.name) + + # Figure out how to handle the tags + for mol in molecules: + mol_slice = slice_topology_by_molecule(topology=self, molecule_tag=0) + openff_mol = to_openff_molecule(mol_slice) + openff_mol.assign_partial_charges(method=method) + for atom in openff_mol.atoms: + print(atom.partial_charge) + + + def write_forcefield(self, filename, overwrite=False): """Save an xml file for all parameters found in the topology. diff --git a/gmso/utils/charges.py b/gmso/utils/charges.py new file mode 100644 index 000000000..893af987a --- /dev/null +++ b/gmso/utils/charges.py @@ -0,0 +1,19 @@ +"""Calculate charges for a gmso Topology.""" +from gmso.external.convert_openmm import to_openff_molecule +from gmso.utils.slicing import slice_topology_by_molecule + +def calculate_charges(topology, method="am1bcc", slice_by="Molecule"): + # Get all of the moleucles types in the topology + charges_dict = dict() + + molecules = set() + for site in topology.sites: + molecules.add(site.molecule.name) + + # Figure out how to handle the tags + for mol in molecules: + mol_slice = slice_topology_by_molecule(topology=topology, molecule_tag=0) + openff_mol = to_openff_molecule(mol_slice) + openff_mol.assign_partial_charges(partial_charge_method=method) + for atom in openff_mol.atoms: + print(atom.partial_charge) diff --git a/gmso/utils/slicing.py b/gmso/utils/slicing.py index 166b948d2..f7b726da4 100644 --- a/gmso/utils/slicing.py +++ b/gmso/utils/slicing.py @@ -1,5 +1,4 @@ -from gmso.core.topology import Topology - +import gmso def slice_topology_by_molecule(topology, molecule_tag, molecule_number=None): """Create a Topology that contains a subset of molecules from another Topology. @@ -67,7 +66,7 @@ def slice_by_sites(topology, sites): topology : gmso.core.topology.Topology The topology being sliced. """ - new_topology = Topology() + new_topology = gmso.core.topology.Topology() connections = set() for site in sites: From d7b075549ab2948bb05e685d0f86d05fd3d8a5c8 Mon Sep 17 00:00:00 2001 From: "pre-commit-ci[bot]" <66853113+pre-commit-ci[bot]@users.noreply.github.com> Date: Thu, 30 Apr 2026 11:13:51 +0000 Subject: [PATCH 6/6] [pre-commit.ci] auto fixes from pre-commit.com hooks for more information, see https://pre-commit.ci --- gmso/core/topology.py | 4 +--- gmso/utils/charges.py | 2 ++ gmso/utils/slicing.py | 1 + 3 files changed, 4 insertions(+), 3 deletions(-) diff --git a/gmso/core/topology.py b/gmso/core/topology.py index 45784417d..1306b3b52 100644 --- a/gmso/core/topology.py +++ b/gmso/core/topology.py @@ -1262,7 +1262,7 @@ def calculate_charges(self, method="am1bcc", slice_by="Molecule"): molecules = set() for site in self.sites: molecules.add(site.molecule.name) - + # Figure out how to handle the tags for mol in molecules: mol_slice = slice_topology_by_molecule(topology=self, molecule_tag=0) @@ -1271,8 +1271,6 @@ def calculate_charges(self, method="am1bcc", slice_by="Molecule"): for atom in openff_mol.atoms: print(atom.partial_charge) - - def write_forcefield(self, filename, overwrite=False): """Save an xml file for all parameters found in the topology. diff --git a/gmso/utils/charges.py b/gmso/utils/charges.py index 893af987a..e336f72a1 100644 --- a/gmso/utils/charges.py +++ b/gmso/utils/charges.py @@ -1,7 +1,9 @@ """Calculate charges for a gmso Topology.""" + from gmso.external.convert_openmm import to_openff_molecule from gmso.utils.slicing import slice_topology_by_molecule + def calculate_charges(topology, method="am1bcc", slice_by="Molecule"): # Get all of the moleucles types in the topology charges_dict = dict() diff --git a/gmso/utils/slicing.py b/gmso/utils/slicing.py index f7b726da4..d94b73eb2 100644 --- a/gmso/utils/slicing.py +++ b/gmso/utils/slicing.py @@ -1,5 +1,6 @@ import gmso + def slice_topology_by_molecule(topology, molecule_tag, molecule_number=None): """Create a Topology that contains a subset of molecules from another Topology.