T12

Author: alongd

arc/job/adapters/ts/linear.py

@@ -632,6 +632,15 @@ def interpolate_addition(rxn: 'ARCReaction',
                 # product compositions, partially displace them now.
                 ts_xyz = _migrate_h_between_fragments(
                     ts_xyz, uni_mol, cut, multi_species, weight)
+                # Revalidate: the migration may have introduced collisions
+                # or other geometry defects not present in the pre-migration guess.
+                is_valid, reason = _validate_ts_guess(
+                    ts_xyz, set(), cut, uni_mol,
+                    label=f'rxn={rxn.label}, frag-fallback-post-migrate')
+                if not is_valid:
+                    logger.debug(f'Linear (rxn={rxn.label}, frag-fallback): '
+                                 f'post-migration guess rejected — {reason}.')
+                    continue
                 ts_xyzs.append(ts_xyz)
 
     # Deduplicate near-identical guesses.
@@ -801,8 +810,10 @@ def _migrate_h_between_fragments(ts_xyz: dict,
     2. Computes element compositions for each fragment.
     3. Compares with product compositions to find H surplus/deficit.
     4. For each surplus fragment, finds the H atom closest to the deficit
-       fragment and displaces it partially (scaled by ``weight``) toward
-       the nearest heavy atom in the deficit fragment.
+       fragment and places it on the donor→acceptor axis at a TS-like
+       distance (interpolated by ``weight``).  Using the donor→acceptor
+       axis instead of a direct H→acceptor line avoids near-collisions
+       with other atoms in the source fragment (e.g. the C in a CO₂ group).
 
     Args:
         ts_xyz: TS guess XYZ from ``_stretch_bond`` (already stretched).
@@ -930,19 +941,72 @@ def heavy_formula(f: Dict[str, int]) -> Dict[str, int]:
         h_dists.sort(key=lambda x: x[1])
 
         for h_idx, _, nearest_heavy in h_dists[:n_to_move]:
-            # Displace H partially toward the nearest heavy atom in the deficit fragment.
-            # At weight=0.5 (default TS), move halfway. The H is between the two fragments.
-            direction = ts_coords[nearest_heavy] - ts_coords[h_idx]
-            dist = float(np.linalg.norm(direction))
-            if dist < 1e-6:
+            # Find the donor heavy atom: the heavy atom bonded to this H in the
+            # source fragment.
+            donor_heavy = None
+            for nbr in uni_mol.atoms[h_idx].bonds.keys():
+                nbr_idx = atom_to_idx[nbr]
+                if symbols[nbr_idx] != 'H' and nbr_idx in fragments[s_fi]:
+                    donor_heavy = nbr_idx
+                    break
+
+            if donor_heavy is None:
                 continue
-            direction /= dist
-            # Target: the H should be at a TS-like distance from the heavy atom.
-            target_dist = get_single_bond_length(symbols[nearest_heavy], 'H') + _PAULING_DELTA
-            # Move the H so it's at target_dist from the heavy atom, scaled by weight.
-            move_dist = (dist - target_dist) * 2.0 * (1.0 - weight)
-            if move_dist > 0:
-                ts_coords[h_idx] += direction * move_dist
+
+            # Triangulate: place H at the intersection of two spheres centred
+            # on donor and acceptor with TS-like radii, choosing the point
+            # closest to the current H position.  This produces a non-collinear
+            # D-H-A geometry that avoids passing through atoms between donor
+            # and acceptor (e.g. the C in a CO₂ group).
+            d_pos = ts_coords[donor_heavy]
+            a_pos = ts_coords[nearest_heavy]
+            h_pos = ts_coords[h_idx]
+            da_vec = a_pos - d_pos
+            da_dist = float(np.linalg.norm(da_vec))
+            if da_dist < 1e-6:
+                continue
+            da_hat = da_vec / da_dist
+
+            d_DH = get_single_bond_length(symbols[donor_heavy], 'H') + _PAULING_DELTA
+            d_AH = get_single_bond_length(symbols[nearest_heavy], 'H') + _PAULING_DELTA
+
+            if da_dist <= d_DH + d_AH:
+                # Spheres overlap → triangulate.
+                x = (da_dist ** 2 + d_DH ** 2 - d_AH ** 2) / (2.0 * da_dist)
+                h_sq = d_DH ** 2 - x ** 2
+                h_perp = np.sqrt(max(h_sq, 0.0))
+
+                proj = d_pos + np.dot(h_pos - d_pos, da_hat) * da_hat
+                perp = h_pos - proj
+                perp_norm = float(np.linalg.norm(perp))
+                if perp_norm > 1e-8:
+                    n_perp = perp / perp_norm
+                else:
+                    ref = np.array([1.0, 0.0, 0.0]) if abs(da_hat[0]) < 0.9 \
+                        else np.array([0.0, 1.0, 0.0])
+                    n_perp = np.cross(da_hat, ref)
+                    n_perp /= np.linalg.norm(n_perp)
+
+                cand_plus = d_pos + x * da_hat + h_perp * n_perp
+                cand_minus = d_pos + x * da_hat - h_perp * n_perp
+
+                # Pick the candidate with greater clearance from source-fragment
+                # heavy atoms (avoids crowding the TS ring interior).
+                def _min_frag_dist(pos):
+                    md = float('inf')
+                    for idx in fragments[s_fi]:
+                        if idx == h_idx or idx == donor_heavy or symbols[idx] == 'H':
+                            continue
+                        md = min(md, float(np.linalg.norm(pos - ts_coords[idx])))
+                    return md
+
+                new_h = cand_plus if _min_frag_dist(cand_plus) >= _min_frag_dist(cand_minus) \
+                    else cand_minus
+            else:
+                # Spheres don't overlap → collinear placement at d_DH from donor.
+                new_h = d_pos + d_DH * da_hat
+
+            ts_coords[h_idx] = new_h
 
     result = {
         'symbols': ts_xyz['symbols'],
@@ -2231,21 +2295,55 @@ def _postprocess_h_migration(xyz: dict,
     return xyz, migrating_hs
 
 
+def _postprocess_group_shift(xyz: dict,
+                              r_mol: 'Molecule',
+                              forming_bonds: List[Tuple[int, int]],
+                              breaking_bonds: List[Tuple[int, int]],
+                              r_label_map: Optional[Dict[str, int]] = None,
+                              ) -> Tuple[dict, Set[int]]:
+    """
+    Post-processing pipeline for group-migration and sigmatropic families.
+
+    These families involve a non-H group migrating between heavy atoms
+    (e.g. 1,2-shiftC, 1,2-shiftS, 1,3-sigmatropic rearrangements).
+    The pipeline applies forming-bond distance correction, umbrella
+    orientation, and universal H fixes, but omits H-transfer-specific
+    steps (donor terminal H staggering) that are not applicable when
+    the migrating atom is not hydrogen.
+
+    Args:
+        xyz: Raw TS guess XYZ dictionary.
+        r_mol: Reactant RMG Molecule (bond topology reference).
+        forming_bonds: List of (i, j) index pairs for bonds that form.
+        breaking_bonds: List of (i, j) index pairs for bonds that break.
+        r_label_map: Family-labelled atom map (unused in this handler).
+
+    Returns:
+        Tuple of (corrected XYZ dict, empty set — no migrating H atoms).
+    """
+    xyz = fix_forming_bond_distances(xyz, r_mol, forming_bonds)
+    xyz = fix_nonreactive_h_distances(xyz, r_mol, migrating_h_indices=set())
+    xyz = fix_crowded_h_atoms(xyz, r_mol, skip_h_indices=set())
+    xyz = fix_migrating_group_umbrella(xyz, r_mol, breaking_bonds, forming_bonds)
+    return xyz, set()
+
+
 # Register H-migration families.
 for _fam in ('intra_H_migration', 'Ketoenol', 'Intra_Disproportionation',
              'Intra_ene_reaction', 'intra_OH_migration',
              'Intra_RH_Add_Endocyclic', 'Intra_RH_Add_Exocyclic',
              'Concerted_Intra_Diels_alder_monocyclic_1,2_shiftH'):
     _FAMILY_POSTPROCESSORS[_fam] = _postprocess_h_migration
 
-# Register group-migration and sigmatropic families — the umbrella
-# flip and forming-bond fix steps from the H-migration pipeline also
-# help these related transfer/rearrangement reactions.
+# Register group-migration and sigmatropic families — these involve
+# non-H group transfer (e.g. CH3, S, ring shifts).  They benefit from
+# forming-bond distance correction and umbrella orientation but do not
+# need H-transfer-specific steps like donor terminal H staggering.
 for _fam in ('1,2_shiftC', '1,2_shiftS',
              '1,3_sigmatropic_rearrangement',
              '6_membered_central_C-C_shift',
              'intra_substitutionCS_isomerization'):
-    _FAMILY_POSTPROCESSORS[_fam] = _postprocess_h_migration
+    _FAMILY_POSTPROCESSORS[_fam] = _postprocess_group_shift
 
 
 # ---------------------------------------------------------------------------
@@ -2318,9 +2416,19 @@ def _postprocess_ts_guess(xyz: dict,
     """
     Dispatch to the appropriate family-specific post-processing handler.
 
-    Looks up ``family`` in :data:`_FAMILY_POSTPROCESSORS`.  If a handler is
-    registered, it is called; otherwise :func:`_postprocess_generic` is used
-    as the safe default (applies only universally safe H fixes).
+    Three handler tiers are available:
+
+    * :func:`_postprocess_h_migration` — for H-transfer families (intra_H_migration,
+      Ketoenol, etc.).  Full pipeline: forming-bond triangulation, donor terminal H
+      staggering, non-reactive H distance fix, crowded-H redistribution, umbrella flip.
+    * :func:`_postprocess_group_shift` — for non-H group migrations (1,2_shiftC/S,
+      sigmatropic, etc.).  Applies forming-bond fix, umbrella, and universal H fixes
+      but omits H-transfer-specific donor staggering.
+    * :func:`_postprocess_generic` — safe default for unknown families.  Only universal
+      H fixes (distance and crowding).
+
+    Looks up ``family`` in :data:`_FAMILY_POSTPROCESSORS`.  If no handler is
+    registered, :func:`_postprocess_generic` is used.
 
     Args:
         xyz: Raw TS guess XYZ dictionary.

arc/job/adapters/ts/linear_test.py

@@ -1540,6 +1540,34 @@ def test_interpolate_1_3_insertion_co2(self):
         for ts_xyz in ts_xyzs:
             print('\n\n***********')
             print(xyz_to_str(ts_xyz))
+        expected_ts = """C      -2.49526563   -1.71744655   -0.05070502
+C      -1.02439736   -1.43442148   -0.09903098
+C      -0.03369873   -2.05724101    0.58071768
+C      -0.33706853   -3.21011130    1.53261903
+C      -0.85997546   -2.75301776    2.88970062
+C       1.47028993   -1.69663277    0.40874963
+C       2.07012721   -1.30878861    1.77489924
+C       1.74154080   -0.52957843   -0.56586530
+C       2.40316640   -3.23990858   -0.35913140
+O       1.93702700   -4.09412667   -1.09797433
+O       3.71891372   -3.28501831   -0.05846324
+H      -3.10560600   -0.84943191    0.20522867
+H      -2.69351927   -2.49442175    0.69448241
+H      -2.82934183   -2.08783564   -1.02499122
+H      -0.78561006   -0.62690380   -0.78689260
+H       0.55596394   -3.82539962    1.69076067
+H      -1.06313426   -3.89101221    1.07201521
+H      -1.79380282   -2.19600714    2.76086637
+H      -0.14816985   -2.09644367    3.39930225
+H      -1.06676022   -3.60034589    3.54974021
+H       2.00787775   -2.13013981    2.49728135
+H       3.13273869   -1.05110182    1.68998787
+H       1.55881292   -0.43721490    2.19959957
+H       1.26204457    0.39674910   -0.22955829
+H       2.81639214   -0.32714424   -0.65140849
+H       1.38896771   -0.75762301   -1.57886997
+H       2.85557208   -2.29609554    0.36706348"""
+        self.assertTrue(almost_equal_coords(ts_xyzs[0], str_to_xyz(expected_ts)))
 
     def test_interpolate_1_3_insertion_ror(self):
         """Test the interpolate_isomerization() function for 1,3_Insertion_ROR: CCOCOCC <=> C=C + CCOCO"""