Merge pull request #29 from openforcefield/clear-conformers-before-co…

…mparison

Clear conformers before adding new conformers for analysis

yammbs/_tests/conftest.py

@@ -97,3 +97,21 @@ def tiny_cache() -> CachedResultCollection:
 def diphenylvinylbenzene():
     """Return 1,2-diphenylvinylbenzene"""
     return Molecule.from_smiles("c1ccc(cc1)C=C(c2ccccc2)c3ccccc3")
+
+
+@pytest.fixture()
+def allicin():
+    """Return allicin, inspired by PQR"""
+    return Molecule.from_smiles(
+        "C=CCSS(=O)CC=C",
+        allow_undefined_stereo=True,
+    )
+
+
+@pytest.fixture()
+def conformers(allicin):
+    other_allicin = Molecule(allicin)
+
+    other_allicin.generate_conformers(n_conformers=10)
+
+    return other_allicin.conformers

yammbs/_tests/unit_tests/test_analysis.py

@@ -2,7 +2,77 @@
 import pytest
 from openff.toolkit import Molecule
 
-from yammbs.analysis import get_internal_coordinate_rmsds
+from yammbs.analysis import get_internal_coordinate_rmsds, get_rmsd, get_tfd
+
+
+class TestAnalysis:
+    def test_rmsd(self, allicin, conformers):
+        # Passing the same conformers should return 0.0
+        last_last = get_rmsd(
+            molecule=allicin,
+            reference=conformers[-1],
+            target=conformers[-1],
+        )
+
+        assert last_last == 0.0
+
+        first_last = get_rmsd(
+            molecule=allicin,
+            reference=conformers[0],
+            target=conformers[-1],
+        )
+
+        assert isinstance(first_last, float)
+
+        first_second = get_rmsd(
+            molecule=allicin,
+            reference=conformers[0],
+            target=conformers[1],
+        )
+
+        assert first_second != first_last
+
+        last_first = get_rmsd(
+            molecule=allicin,
+            reference=conformers[-1],
+            target=conformers[0],
+        )
+
+        assert last_first == first_last
+
+    def test_tfd(self, allicin, conformers):
+        # Passing the same conformers should return 0.0
+        last_last = get_tfd(
+            molecule=allicin,
+            reference=conformers[-1],
+            target=conformers[-1],
+        )
+
+        assert last_last == 0.0
+
+        first_last = get_tfd(
+            molecule=allicin,
+            reference=conformers[0],
+            target=conformers[-1],
+        )
+
+        assert isinstance(first_last, float)
+
+        first_second = get_tfd(
+            molecule=allicin,
+            reference=conformers[0],
+            target=conformers[1],
+        )
+
+        assert first_second != first_last
+
+        last_first = get_tfd(
+            molecule=allicin,
+            reference=conformers[-1],
+            target=conformers[0],
+        )
+
+        assert last_first == first_last
 
 
 class TestInternalCoordinateRMSD:

yammbs/analysis.py

@@ -97,9 +97,14 @@ def get_rmsd(
     from openff.units import Quantity, unit
 
     molecule1 = Molecule(molecule)
+    molecule2 = Molecule(molecule)
+
+    for molecule in (molecule1, molecule2):
+        if molecule.conformers is not None:
+            molecule.conformers.clear()
+
     molecule1.add_conformer(Quantity(reference, unit.angstrom))
 
-    molecule2 = Molecule(molecule)
     molecule2.add_conformer(Quantity(target, unit.angstrom))
 
     # oechem appears to not support named arguments, but it's hard to tell
@@ -210,8 +215,6 @@ def _rdmol(
         molecule: Molecule,
         conformer: Array,
     ):
-        from copy import deepcopy
-
         from openff.units import Quantity, unit
 
         # TODO: Do we need to remap indices?
@@ -236,11 +239,11 @@ def _rdmol(
 
             molecule.remap(mapping_dict=atom_map)
 
-        molecule = deepcopy(molecule)
+        molecule = Molecule(molecule)
+        if molecule.conformers is not None:
+            molecule.conformers.clear()
 
-        molecule.add_conformer(
-            Quantity(conformer, unit.angstrom),
-        )
+        molecule.add_conformer(Quantity(conformer, unit.angstrom))
 
         return molecule.to_rdkit()