Repository: openforcefield/openff-toolkit

Issue:
Make `offmol.smirnoff_impropers` and `offmol.single_impropers`
**Describe the bug**
Not really a bug, possibly naive question about intended behavior in enumerating indices for improper torsions.

I was surprised at the large number of torsions included in `offmol.impropers` for some molecules (e.g. `len(methane.impropers) == 24`).

[The way `offmol.impropers` are enumerated](https://github.com/openforcefield/openforcefield/blob/ce77ebafc29eccd53c3fba4ea011a750d2a2ef81/openforcefield/topology/molecule.py#L4798-L4838) appears maximally permissive, in that it doesn't place a constraint on the connectivity or identity of the central atom (except that it must have connectivity of >= 3). I think it's equivalent to returning all chemical environment matches for the SMARTS pattern `"[*:1]~[*:2](~[*:3])~[*:4]"` ("a central atom and 3 tagged neighbors").

Are there other constraints we wish to impose on the indices returned by `offmol.impropers`? For example, would it be desirable to impose a constraint that the connectivity of the central atom is _exactly_ 3 (`"[*:1]~[X3:2](~[*:3])~[*:4]"` )?

**To Reproduce**
```python
from openforcefield.topology import Molecule

improper_smarts = "[*:1]~[X3:2](~[*:3])~[*:4]"
#permissive_improper_smarts = "[*:1]~[*:2](~[*:3])~[*:4]"

smiles_list = ["C", "N", "CCCCC", 'C'*1000]

def compare_improper_enumeration_schemes(smiles):
    offmol = Molecule.from_smiles(smiles)

    smarts_inds = offmol.chemical_environment_matches(improper_smarts)
    offmol_inds = [tuple([atom.molecule_atom_index for atom in improper]) for improper in offmol.impropers]

    print(f'smiles:\n\t{smiles}')
    print(f'matches for substructure search "{improper_smarts}":\n\t{len(smarts_inds)}')
    print(f'offmol.impropers:\n\t{len(offmol_inds)}')
    print(f'set(smarts_matches) is subset of set(offmol.impropers):\n\t{set(smarts_inds).issubset(offmol_inds)}\n\n')

for smiles in smiles_list:
    compare_improper_enumeration_schemes(smiles)
```

**Output**
```
smiles:
	C
matches for substructure search "[*:1]~[X3:2](~[*:3])~[*:4]":
	0
offmol.impropers:
	24
set(smarts_matches) is subset of set(offmol.impropers):
	True


smiles:
	N
matches for substructure search "[*:1]~[X3:2](~[*:3])~[*:4]":
	6
offmol.impropers:
	6
set(smarts_matches) is subset of set(offmol.impropers):
	True


smiles:
	CCCCC
matches for substructure search "[*:1]~[X3:2](~[*:3])~[*:4]":
	0
offmol.impropers:
	120
set(smarts_matches) is subset of set(offmol.impropers):
	True


smiles:
	CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC
matches for substructure search "[*:1]~[X3:2](~[*:3])~[*:4]":
	0
offmol.impropers:
	24000
set(smarts_matches) is subset of set(offmol.impropers):
	True
```

**Computing environment:**
`openforcefield==0.7.1`

Discussion:

Patch:
diff --git a/docs/releasehistory.rst b/docs/releasehistory.rst
index 72b6c9d5..5d77ab43 100644
--- a/docs/releasehistory.rst
+++ b/docs/releasehistory.rst
@@ -16,6 +16,14 @@ New features
   OpenEye wrapper.
 - `PR #793 <https://github.com/openforcefield/openforcefield/pull/793>`_: Add an initial ELF conformer selection
   implementation which uses RDKit.
+- `PR #799 <https://github.com/openforcefield/openforcefield/pull/799>`_: Closes
+  `Issue #746 <https://github.com/openforcefield/openforcefield/issues/746>`_ by adding
+  :py:meth:`Molecule.smirnoff_impropers <openff.toolkit.topology.FrozenMolecule.smirnoff_impropers>`,
+  :py:meth:`Molecule.amber_impropers <openff.toolkit.topology.FrozenMolecule.amber_impropers>`,
+  :py:meth:`TopologyMolecule.smirnoff_impropers <openff.toolkit.topology.TopologyMolecule.smirnoff_impropers>`,
+  :py:meth:`TopologyMolecule.amber_impropers <openff.toolkit.topology.TopologyMolecule.amber_impropers>`,
+  :py:meth:`Topology.smirnoff_impropers <openff.toolkit.topology.Topology.smirnoff_impropers>`, and
+  :py:meth:`Topology.amber_impropers <openff.toolkit.topology.Topology.amber_impropers>`.
 
 Behavior changed
 """"""""
diff --git a/openff/toolkit/topology/molecule.py b/openff/toolkit/topology/molecule.py
index 77d591c0..d20bc502 100644
--- a/openff/toolkit/topology/molecule.py
+++ b/openff/toolkit/topology/molecule.py
@@ -3518,7 +3518,7 @@ class FrozenMolecule(Serializable):
 
     @property
     def n_impropers(self):
-        """int: number of improper torsions in the Molecule."""
+        """int: number of possible improper torsions in the Molecule."""
         self._construct_torsions()
         return len(self._impropers)
 
@@ -3624,15 +3624,105 @@ class FrozenMolecule(Serializable):
     @property
     def impropers(self):
         """
-        Iterate over all proper torsions in the molecule
-
-        .. todo::
+        Iterate over all improper torsions in the molecule.
 
+        .. todo ::
            * Do we need to return a ``Torsion`` object that collects information about fractional bond orders?
+
+        Returns
+        -------
+        impropers : set of tuple
+            An iterator of tuples, each containing the indices of atoms making
+            up a possible improper torsion.
+
+        See Also
+        --------
+        smirnoff_impropers, amber_impropers
         """
         self._construct_torsions()
         return self._impropers
 
+    @property
+    def smirnoff_impropers(self):
+        """
+        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.
+
+        The exact three orderings that will be applied during parameterization can not be
+        determined in this method, since it requires sorting the particle indices, and
+        those indices may change when this molecule is added to a Topology.
+
+        For more details on the use of three-fold ('trefoil') impropers, see
+        https://open-forcefield-toolkit.readthedocs.io/en/latest/smirnoff.html#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
+        --------
+        impropers, amber_impropers
+
+        """
+        # TODO: Replace with non-cheminformatics-toolkit method
+        #       (ie. just looping over all atoms and finding ones that have 3 bonds?)
+
+        smirnoff_improper_smarts = "[*:1]~[X3:2](~[*:3])~[*:4]"
+        improper_idxs = self.chemical_environment_matches(smirnoff_improper_smarts)
+        smirnoff_impropers = {
+            tuple(self.atoms[idx] for idx in imp) for imp in improper_idxs
+        }
+        return smirnoff_impropers
+
+    @property
+    def amber_impropers(self):
+        """
+        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
+        --------
+        impropers, smirnoff_impropers
+
+        """
+        # TODO: Replace with non-cheminformatics-toolkit method
+        #       (ie. just looping over all atoms and finding ones that have 3 bonds?)
+        amber_improper_smarts = "[X3:1](~[*:2])(~[*:3])~[*:4]"
+        improper_idxs = self.chemical_environment_matches(amber_improper_smarts)
+        amber_impropers = {
+            tuple(self.atoms[idx] for idx in imp) for imp in improper_idxs
+        }
+        return amber_impropers
+
     @property
     def total_charge(self):
         """
@@ -4901,7 +4991,6 @@ class FrozenMolecule(Serializable):
         if not hasattr(self, "_torsions"):
             self._construct_bonded_atoms_list()
 
-            # self._torsions = set()
             self._propers = set()
             self._impropers = set()
             for atom1 in self._atoms:
@@ -4933,7 +5022,6 @@ class FrozenMolecule(Serializable):
                             self._impropers.add(improper)
 
             self._torsions = self._propers | self._impropers
-        # return iter(self._torsions)
 
     def _construct_bonded_atoms_list(self):
         """
diff --git a/openff/toolkit/topology/topology.py b/openff/toolkit/topology/topology.py
index 20cc9830..183e49bb 100644
--- a/openff/toolkit/topology/topology.py
+++ b/openff/toolkit/topology/topology.py
@@ -1119,34 +1119,94 @@ class TopologyMolecule:
 
     @property
     def angles(self):
-        """Iterable of Tuple[TopologyAtom]: iterator over the angles in this Topology."""
+        """Iterable of Tuple[TopologyAtom]: iterator over the angles in this TopologyMolecule."""
         return self._convert_to_topology_atom_tuples(self._reference_molecule.angles)
 
     @property
     def n_angles(self):
-        """int: number of angles in this Topology."""
+        """int: number of angles in this TopologyMolecule."""
         return self._reference_molecule.n_angles
 
     @property
     def propers(self):
-        """Iterable of Tuple[TopologyAtom]: iterator over the proper torsions in this Topology."""
+        """Iterable of Tuple[TopologyAtom]: iterator over the proper torsions in this TopologyMolecule."""
         return self._convert_to_topology_atom_tuples(self._reference_molecule.propers)
 
     @property
     def n_propers(self):
-        """int: number of proper torsions in this Topology."""
+        """int: number of proper torsions in this TopologyMolecule."""
         return self._reference_molecule.n_propers
 
     @property
     def impropers(self):
-        """Iterable of Tuple[TopologyAtom]: iterator over the improper torsions in this Topology."""
+        """Iterable of Tuple[TopologyAtom]: iterator over the possible improper torsions in this TopologyMolecule."""
         return self._convert_to_topology_atom_tuples(self._reference_molecule.impropers)
 
     @property
     def n_impropers(self):
-        """int: number of proper torsions in this Topology."""
+        """int: number of possible improper torsions in this TopologyMolecule."""
         return self._reference_molecule.n_impropers
 
+    @property
+    def smirnoff_impropers(self):
+        """
+        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://open-forcefield-toolkit.readthedocs.io/en/latest/smirnoff.html#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.
+
+        .. todo::
+           * Offer a way to do the keytransform and get the final 3 orderings in this
+             method? How can we keep this logic synced up with the parameterization
+             machinery?
+        """
+        return self._convert_to_topology_atom_tuples(
+            self._reference_molecule.smirnoff_impropers
+        )
+
+    @property
+    def amber_impropers(self):
+        """
+        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.
+        """
+        return self._convert_to_topology_atom_tuples(
+            self._reference_molecule.amber_impropers
+        )
+
     @property
     def virtual_site_start_topology_index(self):
         """Get the topology index of the first virtual site in this TopologyMolecule"""
@@ -1706,16 +1766,85 @@ class Topology(Serializable):
 
     @property
     def n_impropers(self):
-        """int: number of improper torsions in this Topology."""
+        """int: number of possible improper torsions in this Topology."""
         return sum(mol.n_impropers for mol in self._topology_molecules)
 
     @property
     def impropers(self):
-        """Iterable of Tuple[TopologyAtom]: iterator over the improper torsions in this Topology."""
+        """Iterable of Tuple[TopologyAtom]: iterator over the possible improper torsions in this Topology."""
         for topology_molecule in self._topology_molecules:
             for improper in topology_molecule.impropers:
                 yield improper
 
+    @property
+    def smirnoff_impropers(self):
+        """
+        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://open-forcefield-toolkit.readthedocs.io/en/latest/smirnoff.html#impropertorsions
+
+        .. todo:: Offer a way to do the keytransform and get the final 3 orderings in this
+                  method? How can we keep this logic synced up with the parameterization
+                  machinery?
+
+        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
+        --------
+        impropers, amber_impropers
+
+        """
+        for topology_molecule in self._topology_molecules:
+            for smirnoff_improper in topology_molecule.smirnoff_impropers:
+                yield smirnoff_improper
+
+    @property
+    def amber_impropers(self):
+        """
+        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
+        --------
+        impropers, smirnoff_impropers
+        """
+        for topology_molecule in self._topology_molecules:
+            for amber_improper in topology_molecule.amber_impropers:
+                yield amber_improper
+
     class _ChemicalEnvironmentMatch:
         """Represents the match of a given chemical environment query, storing
         both the matched topology atom indices and the indices of the corresponding
diff --git a/openff/toolkit/utils/toolkits.py b/openff/toolkit/utils/toolkits.py
index b9a51281..c84dde31 100644
--- a/openff/toolkit/utils/toolkits.py
+++ b/openff/toolkit/utils/toolkits.py
@@ -2135,7 +2135,7 @@ class OpenEyeToolkitWrapper(ToolkitWrapper):
 
         See Also
         --------
-        ``RDKitToolkitWrapper.apply_elf_conformer_selection``
+        RDKitToolkitWrapper.apply_elf_conformer_selection
 
         Parameters
         ----------
@@ -3840,7 +3840,7 @@ class RDKitToolkitWrapper(ToolkitWrapper):
 
         See Also
         --------
-        ``RDKitToolkitWrapper._elf_select_diverse_conformers``
+        RDKitToolkitWrapper._elf_select_diverse_conformers
 
         Parameters
         ----------


Test patch:
diff --git a/openff/toolkit/tests/test_molecule.py b/openff/toolkit/tests/test_molecule.py
index 34c802a9..f2799792 100644
--- a/openff/toolkit/tests/test_molecule.py
+++ b/openff/toolkit/tests/test_molecule.py
@@ -2016,6 +2016,32 @@ class TestMolecule:
             )
             assert is_not_cyclic or is_three_memebered_ring_torsion(improper)
 
+    @pytest.mark.parametrize(
+        ("molecule", "n_impropers", "n_pruned"),
+        [
+            ("C", 24, 0),
+            ("CC", 48, 0),
+            ("N", 6, 6),
+        ],
+    )
+    def test_pruned_impropers(self, molecule, n_impropers, n_pruned):
+        """Test the amber_impropers and smirnoff_impropers properties"""
+        mol = Molecule.from_smiles(molecule)
+        assert mol.n_impropers == n_impropers
+        assert len(mol.smirnoff_impropers) == n_pruned
+        assert len(mol.amber_impropers) == n_pruned
+
+        # Order not guaranteed, so cannot zip and compare directly
+        for smirnoff_imp in mol.smirnoff_impropers:
+            # Convert SMIRNOFF-style improper into AMBER-style
+            mod_imp = (
+                smirnoff_imp[1],
+                smirnoff_imp[0],
+                smirnoff_imp[2],
+                smirnoff_imp[3],
+            )
+            assert mod_imp in mol.amber_impropers
+
     @pytest.mark.parametrize("molecule", mini_drug_bank())
     def test_torsions(self, molecule):
         """Test torsions property"""
diff --git a/openff/toolkit/tests/test_topology.py b/openff/toolkit/tests/test_topology.py
index fc375a8e..ec5b37df 100644
--- a/openff/toolkit/tests/test_topology.py
+++ b/openff/toolkit/tests/test_topology.py
@@ -587,6 +587,26 @@ class TestTopology(TestCase):
             top_improper_atoms3, mol_improper_atoms2, transformed_dict_cls=ImproperDict
         )
 
+    def test_pruned_impropers(self):
+        """Test {smirnoff|amber}_impropers from the Topology API"""
+        top = Topology.from_molecules(
+            [Molecule.from_smiles(smi) for smi in ["N", "C=C"]]
+        )
+
+        assert len([*top.smirnoff_impropers]) == 18
+        assert len([*top.amber_impropers]) == 18
+
+        # Order not guaranteed, so cannot zip and compare directly
+        for smirnoff_imp in top.smirnoff_impropers:
+            # Convert SMIRNOFF-style improper into AMBER-style
+            mod_imp = (
+                smirnoff_imp[1],
+                smirnoff_imp[0],
+                smirnoff_imp[2],
+                smirnoff_imp[3],
+            )
+            assert mod_imp in top.amber_impropers
+
     # test_get_fractional_bond_order
     # test_two_of_same_molecule
     # test_two_different_molecules