Fixes stereo correction and generalizes to TSs

TODO

@@ -10,6 +10,7 @@
 !!! [ ] Stereo-correct join for bond stereo
     LEFT OFF AT automol.geom.ts.join_at_forming_bond()
     -> SEE test_ts_geometry.ipynb
+    [ ] Figure out and TEST reagent reordering
     [ ] _4geom.py
     [ ] Include zmatrix data in Reaction string
     [ ] _5zmat.py

automol/geom/_conv.py

@@ -675,6 +675,9 @@ def change_zmatrix_row_values(
     tol = 2.0 * phycon.DEG2RAD
     lin = numpy.pi
 
+    print(automol.graph.string(gra))
+    print("idx1", idx1)
+    print("idx", idx)
     idxs = automol.graph.branch_atom_keys(gra, idx1, idx)
 
     # Note that the coordinates will change throughout, but the change

automol/geom/base/_core.py

@@ -110,7 +110,7 @@ def coordinates(geo, idxs=None, angstrom=False):
     else:
         xyzs = ()
     xyzs = xyzs if not angstrom else numpy.multiply(xyzs, phycon.BOHR2ANG)
-    xyzs = tuple(tuple(xyz) for idx, xyz in enumerate(xyzs) if idx in idxs)
+    xyzs = tuple(map(tuple, map(xyzs.__getitem__, idxs)))
 
     return xyzs
 

automol/geom/ts.py

@@ -48,10 +48,8 @@ def geometry_from_reactants(
     :return: TS geometry
     :rtype: automol geom data structure
     """
-    # 1. If there are multiple reactants, combine the geometries
-    if len(geos) == 1:
-        (geo,) = geos
-    else:
+    # 1. Join geometries for bimolecular reactions, yielding a single starting structure
+    if len(geos) > 1:
         assert (
             len(geos) == 2 and len(automol.graph.base.ts.forming_bond_keys(tsg)) == 1
         ), "Generating a TS geometry for this case is not implemented:\n{tsg}"
@@ -62,12 +60,12 @@ def geometry_from_reactants(
             fdist_factor=fdist_factor,
             debug_visualize=debug_visualize,
         )
+    else:
+        (geo,) = geos
 
     # 2. Correct the stereochemistry against the TS graph, so it is consistent with both
     # reactants and products
-    ts_geo = geo
-    # ts_geo = automol.graph.stereo_corrected_geometry(tsg, geo,
-    # geo_idx_dct=geo_idx_dct)
+    ts_geo = automol.graph.stereo_corrected_geometry(tsg, geo, geo_idx_dct=geo_idx_dct)
 
     print(bdist_factor)
 
@@ -106,6 +104,14 @@ def join_at_forming_bond(
     :returns: The joined geometry
     :rtype: automol geom data structure
     """
+    rcts_gra = automol.graph.ts.reactants_graph(tsg, stereo=False)
+    geos_gra = automol.graph.base.union_from_sequence(
+        list(map(automol.geom.graph, geos)), shift_keys=True
+    )
+    assert geos_gra == rcts_gra, (
+        f"The geometries don't match the TS graph. Have they been reordered?"
+        f"\ngeos:\n{geos}\ntsg:\n{tsg}"
+    )
     akeys = sorted(automol.graph.base.atom_keys(tsg))
     geo_idx_dct = (
         {k: i for i, k in enumerate(akeys)} if geo_idx_dct is None else geo_idx_dct
@@ -121,6 +127,8 @@ def join_at_forming_bond(
     len1, len2 = map(count, geos)
     idxs1 = tuple(range(len1))
     idxs2 = tuple(range(len1, len1 + len2))
+    print("idxs1", idxs1)
+    print("idxs2", idxs2)
     (fidx1,) = frm_key & set(idxs1)
     (fidx2,) = frm_key & set(idxs2)
     geo = sum(geos, ())
@@ -146,7 +154,8 @@ def join_at_forming_bond(
     if debug_visualize:
         view = py3dmol_.create_view()
         # Visualize the geomtry
-        view = py3dmol_view(geo, view=view)
+        rcts_gra = automol.graph.ts.reactants_graph(tsg)
+        view = py3dmol_view(geo, gra=rcts_gra, view=view)
         # Visualize the first direction vector
         view = py3dmol_.view_vector(rvec1, orig_xyz=fxyz1, view=view)
         # Visualize the rotational axis
@@ -208,6 +217,7 @@ def _direction_vector(key):
         # If this is a vinyl radical, the reacting electron is pointed in-plane at a 120
         # degree angle to the double bond
         if key in hyd_dct:
+            print("AA HERE?")
             (nkey,) = hyd_dct[key]
             nxyz = xyzs[nkey]
             rvec = vec.unit_norm(numpy.subtract(xyz, nxyz))
@@ -252,9 +262,11 @@ def _direction_vector(key):
         # Otherwise, assume this is a reacting pi-electron, which is pointed
         # perpendicular to the plane of the neighboring atoms
         else:
+            print("BB OR HERE?")
             nbh = automol.graph.base.atom_neighborhood(no_rxn_bnd_gra, key)
             # Use this to find the
             nkeys = automol.graph.base.atom_keys(nbh)
+            print("nkeys", nkeys)
             nxyzs = coordinates(geo, idxs=nkeys)
             rvec = vec.best_unit_perpendicular(xyzs=nxyzs)
 

automol/graph/__init__.py

@@ -172,6 +172,8 @@
 from automol.graph.base._algo import ring_systems
 from automol.graph.base._algo import ring_systems_atom_keys
 from automol.graph.base._algo import ring_systems_bond_keys
+from automol.graph.base._algo import spiro_atoms_grouped_neighbor_keys
+from automol.graph.base._algo import spiro_atom_keys
 from automol.graph.base._algo import is_ring_system
 from automol.graph.base._algo import ring_system_decomposed_atom_keys
 from automol.graph.base._algo import ring_systems_decomposed_atom_keys
@@ -222,6 +224,7 @@
 from automol.graph.base._geom import geometries_parity_mismatches
 # # corrections
 from automol.graph.base._geom import linear_vinyl_corrected_geometry
+from automol.graph.base._geom import geometry_pseudorotate_atom
 from automol.graph.base._geom import geometry_rotate_bond
 from automol.graph.base._geom import geometry_dihedrals_near_value
 from automol.graph.base._geom import perturb_geometry_planar_dihedrals
@@ -258,7 +261,6 @@
 # stereo functions:
 # # core functions
 from automol.graph.base._stereo import expand_stereo
-from automol.graph.base._stereo import expand_stereo_with_priorities_and_amchis
 # # stereo correction
 from automol.graph.base._stereo import stereo_corrected_geometry
 # functional groups code:
@@ -473,6 +475,8 @@
     'ring_systems',
     'ring_systems_atom_keys',
     'ring_systems_bond_keys',
+    'spiro_atoms_grouped_neighbor_keys',
+    'spiro_atom_keys',
     'is_ring_system',
     'ring_system_decomposed_atom_keys',
     'ring_systems_decomposed_atom_keys',
@@ -523,6 +527,7 @@
     'geometries_parity_mismatches',
     # # corrections
     'linear_vinyl_corrected_geometry',
+    'geometry_pseudorotate_atom',
     'geometry_rotate_bond',
     'geometry_dihedrals_near_value',
     'perturb_geometry_planar_dihedrals',
@@ -559,7 +564,6 @@
     # stereo functions:
     # # core functions
     'expand_stereo',
-    'expand_stereo_with_priorities_and_amchis',
     # # stereo correction
     'stereo_corrected_geometry',
     # functional groups code:

automol/graph/base/__init__.py

@@ -179,6 +179,8 @@
 from automol.graph.base._algo import ring_systems
 from automol.graph.base._algo import ring_systems_atom_keys
 from automol.graph.base._algo import ring_systems_bond_keys
+from automol.graph.base._algo import spiro_atoms_grouped_neighbor_keys
+from automol.graph.base._algo import spiro_atom_keys
 from automol.graph.base._algo import is_ring_system
 from automol.graph.base._algo import ring_system_decomposed_atom_keys
 from automol.graph.base._algo import ring_systems_decomposed_atom_keys
@@ -229,6 +231,7 @@
 from automol.graph.base._geom import geometries_parity_mismatches
 # # corrections
 from automol.graph.base._geom import linear_vinyl_corrected_geometry
+from automol.graph.base._geom import geometry_pseudorotate_atom
 from automol.graph.base._geom import geometry_rotate_bond
 from automol.graph.base._geom import geometry_dihedrals_near_value
 from automol.graph.base._geom import perturb_geometry_planar_dihedrals
@@ -265,7 +268,6 @@
 # stereo functions:
 # # core functions
 from automol.graph.base._stereo import expand_stereo
-from automol.graph.base._stereo import expand_stereo_with_priorities_and_amchis
 # # stereo correction
 from automol.graph.base._stereo import stereo_corrected_geometry
 # functional groups code:
@@ -463,6 +465,8 @@
     'ring_systems',
     'ring_systems_atom_keys',
     'ring_systems_bond_keys',
+    'spiro_atoms_grouped_neighbor_keys',
+    'spiro_atom_keys',
     'is_ring_system',
     'ring_system_decomposed_atom_keys',
     'ring_systems_decomposed_atom_keys',
@@ -513,6 +517,7 @@
     'geometries_parity_mismatches',
     # # corrections
     'linear_vinyl_corrected_geometry',
+    'geometry_pseudorotate_atom',
     'geometry_rotate_bond',
     'geometry_dihedrals_near_value',
     'perturb_geometry_planar_dihedrals',
@@ -549,7 +554,6 @@
     # stereo functions:
     # # core functions
     'expand_stereo',
-    'expand_stereo_with_priorities_and_amchis',
     # # stereo correction
     'stereo_corrected_geometry',
     # functional groups code: