Merge pull request #545 from avcopan/dev

Fix: Strained bridgehead filtering

automol/graph/base/_05stereo.py

@@ -239,9 +239,9 @@ def stereoatom_bridgehead_pairs(
     bhp_dct = {}
     for key1, key2 in itertools.combinations(cand_atm_dct.keys(), r=2):
         paths = simple_paths_between_atoms(gra, key1, key2)
-        conn_nkeys1 = tuple(sorted(set(p[1] for p in paths)))
-        conn_nkeys2 = tuple(sorted(set(p[-2] for p in paths)))
-        if len(set(conn_nkeys1)) == len(set(conn_nkeys2)) == 3:
+        conn_dct = {p[1]: p[-2] for p in paths}
+        if len(set(conn_dct.keys())) == len(set(conn_dct.values())) == 3:
+            conn_nkeys1, conn_nkeys2 = zip(*sorted(conn_dct.items()), strict=True)
             bhp_dct[(key1, key2)] = (conn_nkeys1, conn_nkeys2)
 
     return bhp_dct

automol/graph/base/_11stereo.py

@@ -156,11 +156,12 @@ def _select_ts_canonical_direction_priorities(gprs):
 
 def _remove_strained_stereoisomers_from_expansion(gps, cand_dct):
     """Remove strained stereoisomers from an expansion"""
-    gps = list(gps)
-    gra = without_stereo(gps[0][0])
+    gps0 = list(gps)
+    gra = without_stereo(gps0[0][0])
     bhp_dct = stereoatom_bridgehead_pairs(gra, cand_dct)
 
-    for gra, pri_dct in gps:
+    gps = gps0.copy()
+    for gra, pri_dct in gps0:
         par_dct = stereo_parities(gra)
         can_nkeys_dct, _ = stereocenter_candidates_grouped(cand_dct, pri_dct=pri_dct)
         for (key1, key2), (conn_nkeys1, conn_nkeys2) in bhp_dct.items():
@@ -176,11 +177,12 @@ def _remove_strained_stereoisomers_from_expansion(gps, cand_dct):
 
             # If the parities relative to the above ordering are not opposite, then the
             # configuration of the bridgehead pair is strained
+            # Exception: Adjacent bridgehead pairs will have groups on the same side
             if par1 is not None and par2 is not None:
                 sgn1 = util.is_odd_permutation(srt_nkeys1, can_nkeys1)
                 sgn2 = util.is_odd_permutation(srt_nkeys2, can_nkeys2)
                 is_strained = not par1 ^ par2 ^ sgn1 ^ sgn2
-                if is_strained and (gra, pri_dct) in gps:
+                if is_strained:
                     gps.remove((gra, pri_dct))
 
     return gps
@@ -548,9 +550,9 @@ def geometry_pseudorotate_atom(
 
     rxn_keys = ts_reacting_atom_keys(gra)
     rsy_keys_lst = []
-    for rgra in (gra_reac,gra_prod):
+    for rgra in (gra_reac, gra_prod):
         rsy_keys_lst.extend(ring_systems_atom_keys(rgra, lump_spiro=False))
-    rsy_keys_lst=list(set(rsy_keys_lst))
+    rsy_keys_lst = list(set(rsy_keys_lst))
     nkeys = atom_neighbor_atom_keys(gra, key)
     # Gather neighbors connected in a ring system
     ring_nkey_sets = [nkeys & ks for ks in rsy_keys_lst if nkeys & ks]
@@ -584,9 +586,7 @@ def geometry_pseudorotate_atom(
 
     # Identify the remaining keys to be rotated
     rot_nkeys = nkeys - {nkey1, nkey2}
-    rot_keys = set(
-        itertools.chain(*(branch_atom_keys(gra, key, k) for k in rot_nkeys))
-    )
+    rot_keys = set(itertools.chain(*(branch_atom_keys(gra, key, k) for k in rot_nkeys)))
 
     geo = geom_base.rotate(geo, rot_axis, ang, orig_xyz=xyz, idxs=rot_keys)
 

automol/tests/test_graph.py

@@ -1868,6 +1868,7 @@ def test__stereogenic_keys():
 
 def test__expand_stereo():
     """test graph.expand_stereo"""
+    print(graph.smiles(C2H2CL2F2_CGR))
     assert graph.expand_stereo(C2H2CL2F2_CGR) == C2H2CL2F2_SGRS
     assert graph.expand_stereo(C3H3CL2F3_CGR) == C3H3CL2F3_SGRS
     # When symmetry equivalents are filtered out, we can't guarantee that the
@@ -1953,42 +1954,32 @@ def test__expand_stereo():
     assert len(graph.expand_stereo(gra, enant=False, symeq=True)) == 2
 
 
-def test__expand_stereo__strained():
+@pytest.mark.parametrize(
+    "smi,npars1,npars2,par_dct",
+    [
+        ("C1C2OC2CC1", 3, 1, {1: False, 3: True}),
+        ("C1(O2)CC2CC1", 3, 1, {0: False, 3: True}),
+        ("C=1C2CC(C=1)O2", 3, 1, {1: False, 3: True}),
+    ]
+)
+def test__expand_stereo__strained(smi, npars1, npars2, par_dct):
     """test removal of strained stereoisomers from stereoexpansion"""
-    # C=1C2CC(C=1)O2
-    gra = (
-        {
-            0: ("C", 1, None),
-            1: ("C", 1, None),
-            2: ("C", 2, None),
-            3: ("C", 1, None),
-            4: ("C", 1, None),
-            5: ("O", 0, None),
-        },
-        {
-            frozenset({1, 4}): (1, None),
-            frozenset({2, 3}): (1, None),
-            frozenset({0, 3}): (1, None),
-            frozenset({4, 5}): (1, None),
-            frozenset({0, 1}): (1, None),
-            frozenset({2, 4}): (1, None),
-            frozenset({3, 5}): (1, None),
-        },
-    )
+    print(f"smi = {smi}")
+    gra = automol.smiles.graph(smi)
 
     # Without removing strained stereoisomers, there are three distinct possibilities
     sgras = graph.expand_stereo(gra, strained=True)
-    assert len(sgras) == 3, sgras
+    assert len(sgras) == npars1, f"{len(sgras)} != {npars1}"
 
     # With removing strained stereoisomers, there is only one distinct possibility
     sgras = graph.expand_stereo(gra, strained=False)
-    assert len(sgras) == 1, sgras
+    assert len(sgras) == npars2, f"{len(sgras)} != {npars2}"
 
     (sgra,) = sgras
-    par_dct = automol.util.dict_.filter_by_value(
+    par_dct_ = automol.util.dict_.filter_by_value(
         automol.graph.stereo_parities(sgra), lambda x: x is not None
     )
-    assert par_dct == {3: False, 4: True}, par_dct
+    assert par_dct_ == par_dct, f"{par_dct_} != {par_dct}"
 
 
 def test__ring_systems():
@@ -2606,4 +2597,6 @@ def test__atom_hypervalencies():
     # test__stereo_corrected_geometry()
     # test__atom_hypervalencies()
     # test__align_with_geometry()
-    test__embed__clean_geometry()
+    # test__embed__clean_geometry()
+    # test__expand_stereo__strained("C=1C2CC(C=1)O2", 3, 1, {1: False, 3: True})
+    test__expand_stereo()