diff --git a/arc/scheduler.py b/arc/scheduler.py index c56fae7d72..adb04bb8ed 100644 --- a/arc/scheduler.py +++ b/arc/scheduler.py @@ -3100,9 +3100,12 @@ def check_directed_scan(self, label, pivots, scan, energies): # a lower conformation was found deg_increment = actions[1] self.species_dict[label].set_dihedral(scan=scan, index=1, deg_increment=deg_increment) - is_isomorphic = self.species_dict[label].check_xyz_isomorphism( - allow_nonisomorphic_2d=self.allow_nonisomorphic_2d, - xyz=self.species_dict[label].initial_xyz) + if self.species_dict[label].is_ts: + is_isomorphic = True + else: + is_isomorphic = self.species_dict[label].check_xyz_isomorphism( + allow_nonisomorphic_2d=self.allow_nonisomorphic_2d, + xyz=self.species_dict[label].initial_xyz) if is_isomorphic: self.delete_all_species_jobs(label) # Remove all completed rotor calculation information @@ -3162,7 +3165,10 @@ def check_directed_scan_job(self, label: str, job: 'JobAdapter'): """ if job.job_status[1]['status'] == 'done': xyz = parser.parse_geometry(log_file_path=job.local_path_to_output_file) - is_isomorphic = self.species_dict[label].check_xyz_isomorphism(xyz=xyz, verbose=False) + if self.species_dict[label].is_ts: + is_isomorphic = True + else: + is_isomorphic = self.species_dict[label].check_xyz_isomorphism(xyz=xyz, verbose=False) for rotor_dict in self.species_dict[label].rotors_dict.values(): if rotor_dict['pivots'] == job.pivots: key = tuple(f'{dihedral:.2f}' for dihedral in job.dihedrals) @@ -3395,9 +3401,12 @@ def troubleshoot_scan_job(self, break else: # If 'change conformer' is not used, check for isomorphism. - is_isomorphic = self.species_dict[label].check_xyz_isomorphism( - allow_nonisomorphic_2d=self.allow_nonisomorphic_2d, - xyz=new_xyz) + if self.species_dict[label].is_ts: + is_isomorphic = True + else: + is_isomorphic = self.species_dict[label].check_xyz_isomorphism( + allow_nonisomorphic_2d=self.allow_nonisomorphic_2d, + xyz=new_xyz) if is_isomorphic: self.species_dict[label].final_xyz = new_xyz # Remove all completed rotor calculation information. diff --git a/arc/scheduler_test.py b/arc/scheduler_test.py index cdc4b17f07..a1a9f9c710 100644 --- a/arc/scheduler_test.py +++ b/arc/scheduler_test.py @@ -6,7 +6,7 @@ """ import unittest -from unittest.mock import patch +from unittest.mock import MagicMock, patch import os import shutil @@ -1005,6 +1005,91 @@ def test_switch_ts_rotors_reset(self, mock_run_opt): # rotors_dict=None must be preserved — do not re-enable rotor scans. self.assertIsNone(sched2.species_dict[ts_label2].rotors_dict) + def test_check_directed_scan_job_skips_isomorphism_for_ts(self): + """check_directed_scan_job must not call check_xyz_isomorphism for a TS; is_isomorphic is recorded as True.""" + ts_xyz = str_to_xyz("""N 0.91779059 0.51946178 0.00000000 + H 1.81402049 1.03819414 0.00000000 + H 0.00000000 0.00000000 0.00000000 + H 0.91779059 1.22790192 0.72426890""") + ts_spc = ARCSpecies(label='TS_dirscan', is_ts=True, xyz=ts_xyz, multiplicity=1, charge=0, + compute_thermo=False) + ts_spc.rotors_dict = {0: {'pivots': [1, 2], 'directed_scan': {}}} + + project_directory = os.path.join(ARC_PATH, 'Projects', 'arc_project_ts_iso_dirscan') + self.addCleanup(shutil.rmtree, project_directory, ignore_errors=True) + sched = Scheduler(project='test_ts_iso_dirscan', ess_settings=self.ess_settings, + species_list=[ts_spc], + opt_level=Level(repr=default_levels_of_theory['opt']), + freq_level=Level(repr=default_levels_of_theory['freq']), + sp_level=Level(repr=default_levels_of_theory['sp']), + ts_guess_level=Level(repr=default_levels_of_theory['ts_guesses']), + project_directory=project_directory, + testing=True, + job_types=self.job_types1, + ) + + job_mock = MagicMock() + job_mock.job_status = [None, {'status': 'done'}] + job_mock.local_path_to_output_file = '/fake/path.log' + job_mock.pivots = [1, 2] + job_mock.dihedrals = [45.0] + job_mock.ess_trsh_methods = [] + + with patch('arc.species.species.ARCSpecies.check_xyz_isomorphism') as mock_iso, \ + patch('arc.scheduler.parser.parse_geometry', return_value=ts_xyz), \ + patch('arc.scheduler.parser.parse_e_elect', return_value=-123.45): + sched.check_directed_scan_job(label='TS_dirscan', job=job_mock) + + mock_iso.assert_not_called() + recorded = sched.species_dict['TS_dirscan'].rotors_dict[0]['directed_scan'][('45.00',)] + self.assertTrue(recorded['is_isomorphic']) + + @patch('arc.scheduler.Scheduler.run_opt_job') + def test_troubleshoot_scan_job_skips_isomorphism_for_ts(self, mock_run_opt): + """troubleshoot_scan_job must not call check_xyz_isomorphism for a TS when applying 'change conformer'.""" + ts_xyz = str_to_xyz("""N 0.91779059 0.51946178 0.00000000 + H 1.81402049 1.03819414 0.00000000 + H 0.00000000 0.00000000 0.00000000 + H 0.91779059 1.22790192 0.72426890""") + new_xyz = str_to_xyz("""N 0.91000000 0.52000000 0.00000000 + H 1.81000000 1.04000000 0.00000000 + H 0.00000000 0.00000000 0.00000000 + H 0.91000000 1.23000000 0.72000000""") + ts_spc = ARCSpecies(label='TS_trsh', is_ts=True, xyz=ts_xyz, multiplicity=1, charge=0, + compute_thermo=False) + ts_spc.rotors_dict = {0: {'pivots': [1, 2], 'scan': [3, 1, 2, 4], 'scan_path': '', + 'invalidation_reason': '', 'success': None, 'symmetry': None, + 'times_dihedral_set': 0, 'trsh_methods': [], 'trsh_counter': 0}} + + project_directory = os.path.join(ARC_PATH, 'Projects', 'arc_project_ts_iso_trsh') + self.addCleanup(shutil.rmtree, project_directory, ignore_errors=True) + sched = Scheduler(project='test_ts_iso_trsh', ess_settings=self.ess_settings, + species_list=[ts_spc], + opt_level=Level(repr=default_levels_of_theory['opt']), + freq_level=Level(repr=default_levels_of_theory['freq']), + sp_level=Level(repr=default_levels_of_theory['sp']), + ts_guess_level=Level(repr=default_levels_of_theory['ts_guesses']), + project_directory=project_directory, + testing=True, + job_types=self.job_types1, + ) + sched.trsh_ess_jobs = True + sched.trsh_rotors = True + + job_mock = MagicMock() + job_mock.species_label = 'TS_trsh' + job_mock.rotor_index = 0 + job_mock.torsions = [[3, 1, 2, 4]] + job_mock.job_name = 'scan_a200' + + with patch('arc.species.species.ARCSpecies.check_xyz_isomorphism') as mock_iso, \ + patch('arc.scheduler.Scheduler.delete_all_species_jobs'): + sched.troubleshoot_scan_job(job=job_mock, methods={'change conformer': new_xyz}) + + mock_iso.assert_not_called() + self.assertEqual(sched.species_dict['TS_trsh'].final_xyz, new_xyz) + mock_run_opt.assert_called_once() + @classmethod def tearDownClass(cls): """ diff --git a/arc/species/species.py b/arc/species/species.py index 76a9105dfd..880f2e03f2 100644 --- a/arc/species/species.py +++ b/arc/species/species.py @@ -1617,6 +1617,8 @@ def mol_from_xyz(self, Important for TS searches and for identifying rotor indices. This works by generating a molecule from xyz and using the 2D structure to confirm that the perceived molecule is correct. + For TSs, the perceived molecule is accepted without enforcing + 2D-graph isomorphism, since TS connectivity is not strictly defined. If ``xyz`` is not given, the species xyz attribute will be used. Args: @@ -1639,28 +1641,32 @@ def mol_from_xyz(self, n_fragments=self.get_n_fragments(), ) if perceived_mol is not None: - allow_nonisomorphic_2d = (self.charge is not None and self.charge) \ - or self.mol.has_charge() or perceived_mol.has_charge() \ - or (self.multiplicity is not None and self.multiplicity >= 3) \ - or self.mol.multiplicity >= 3 or perceived_mol.multiplicity >= 3 - isomorphic = self.check_xyz_isomorphism(mol=perceived_mol, - xyz=xyz, - allow_nonisomorphic_2d=allow_nonisomorphic_2d) - if not isomorphic: - logger.warning(f'XYZ and the 2D graph representation for {self.label} are not isomorphic.\nGot ' - f'xyz:\n{xyz}\n\nwhich corresponds to {self.mol.copy(deep=True).to_smiles()}\n' - f'{self.mol.copy(deep=True).to_adjacency_list()}\n\nand: ' - f'{self.mol.copy(deep=True).to_smiles()}\n' - f'{self.mol.copy(deep=True).to_adjacency_list()}') - raise SpeciesError(f'XYZ and the 2D graph representation for {self.label} are not compliant.') - if not self.keep_mol: - if is_mol_valid(perceived_mol, charge=self.charge, multiplicity=self.multiplicity, n_radicals=self.number_of_radicals): + if self.is_ts: + if not self.keep_mol: self.mol = perceived_mol - else: - try: - order_atoms(ref_mol=perceived_mol, mol=self.mol) - except SanitizationError: + else: + allow_nonisomorphic_2d = (self.charge is not None and self.charge) \ + or self.mol.has_charge() or perceived_mol.has_charge() \ + or (self.multiplicity is not None and self.multiplicity >= 3) \ + or self.mol.multiplicity >= 3 or perceived_mol.multiplicity >= 3 + isomorphic = self.check_xyz_isomorphism(mol=perceived_mol, + xyz=xyz, + allow_nonisomorphic_2d=allow_nonisomorphic_2d) + if not isomorphic: + logger.warning(f'XYZ and the 2D graph representation for {self.label} are not isomorphic.\nGot ' + f'xyz:\n{xyz}\n\nwhich corresponds to {self.mol.copy(deep=True).to_smiles()}\n' + f'{self.mol.copy(deep=True).to_adjacency_list()}\n\nand: ' + f'{perceived_mol.copy(deep=True).to_smiles()}\n' + f'{perceived_mol.copy(deep=True).to_adjacency_list()}') + raise SpeciesError(f'XYZ and the 2D graph representation for {self.label} are not compliant.') + if not self.keep_mol: + if is_mol_valid(perceived_mol, charge=self.charge, multiplicity=self.multiplicity, n_radicals=self.number_of_radicals): self.mol = perceived_mol + else: + try: + order_atoms(ref_mol=perceived_mol, mol=self.mol) + except SanitizationError: + self.mol = perceived_mol else: perceived_mol = perceive_molecule_from_xyz(xyz, charge=self.charge, diff --git a/arc/species/species_test.py b/arc/species/species_test.py index 466217ff36..3e98186121 100644 --- a/arc/species/species_test.py +++ b/arc/species/species_test.py @@ -1684,6 +1684,43 @@ def test_mol_from_xyz(self): radical_count += atom.radical_electrons self.assertEqual(radical_count, 2) + def test_ts_mol_from_xyz_skips_isomorphism_enforcement(self): + """Test that a TS accepts the perceived xyz-derived molecule even if a stored 2D graph disagrees.""" + xyz = {'symbols': ('O', 'O', 'H', 'C', 'C', 'C', 'C', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H', 'H'), + 'isotopes': (16, 16, 1, 12, 12, 12, 12, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1), + 'coords': ((1.570004, 0.919237, -0.063628), (2.346384, -0.215837, 0.116826), + (2.578713, -0.467096, -0.784756), (-0.673058, -1.117816, -1.045216), + (-0.76037, -0.036132, 0.001231), (-0.850886, -0.54288, 1.418092), + (-1.694638, 1.099253, -0.333714), (-1.612555, -1.681809, -1.088997), + (-0.491563, -0.700452, -2.03735), (0.122945, -1.827698, -0.811666), + (0.427835, 0.508898, -0.034066), (-0.034676, -1.231859, 1.641118), + (-1.795224, -1.080235, 1.568607), (-0.814118, 0.27711, 2.136337), + (-2.733958, 0.749029, -0.320335), (-1.610541, 1.910407, 0.391085), + (-1.494252, 1.501954, -1.327915))} + adj = """multiplicity 2 +1 O u0 p2 c0 {2,S} {17,S} +2 O u1 p2 c0 {1,S} +3 C u0 p0 c0 {4,S} {5,S} {6,S} {7,S} +4 C u0 p0 c0 {3,S} {8,S} {9,S} {10,S} +5 C u0 p0 c0 {3,S} {11,S} {12,S} {13,S} +6 C u0 p0 c0 {3,S} {14,S} {15,S} {16,S} +7 H u0 p0 c0 {3,S} +8 H u0 p0 c0 {4,S} +9 H u0 p0 c0 {4,S} +10 H u0 p0 c0 {4,S} +11 H u0 p0 c0 {5,S} +12 H u0 p0 c0 {5,S} {17,vdW} +13 H u0 p0 c0 {5,S} +14 H u0 p0 c0 {6,S} +15 H u0 p0 c0 {6,S} +16 H u0 p0 c0 {6,S} +17 H u0 p0 c0 {1,S} {12,vdW}""" + + spc = ARCSpecies(label='TS0', adjlist=adj, xyz=xyz, is_ts=True, multiplicity=2, charge=0) + self.assertIn('3 H u0 p0 c0 {2,S}', spc.mol.to_adjacency_list()) + self.assertIn('11 H u0 p0 c0 {1,S} {5,S}', spc.mol.to_adjacency_list()) + self.assertNotIn('{17,vdW}', spc.mol.to_adjacency_list()) + def test_consistent_atom_order(self): """Test that the atom order is preserved whether starting from SMILES or from xyz""" xyz9 = """O -1.17310019 -0.30822930 0.16269772 diff --git a/docs/source/TS_search.rst b/docs/source/TS_search.rst index 73513c9afb..87199e4b8c 100644 --- a/docs/source/TS_search.rst +++ b/docs/source/TS_search.rst @@ -47,6 +47,9 @@ Outputs and validation """""""""""""""""""""" Validated TS results are reported in the project output (log files and generated artifacts), together with the supporting calculations (optimization, frequency, and IRC). +ARC does not require TS geometries to be isomorphic with a stored 2D adjacency list, since a TS does not have a +single strict graph representation. Instead, TS validation relies on TS-specific checks such as the imaginary +frequency, normal mode displacement analysis, IRC results, and energetic consistency. Reference """"""""" diff --git a/docs/source/advanced.rst b/docs/source/advanced.rst index 091fa90b41..8393a1452d 100644 --- a/docs/source/advanced.rst +++ b/docs/source/advanced.rst @@ -481,14 +481,18 @@ __ Truhlar1_ Isomorphism checks ^^^^^^^^^^^^^^^^^^ -When a species is defined using a 2D graph (i.e., SMILES, AdjList, or InChI), an isomorphism check +For non-TS species defined using a 2D graph (i.e., SMILES, AdjList, or InChI), an isomorphism check is performed on the optimized geometry (all conformers and final optimization). -If the molecule perceived from the 3D coordinate is not isomorphic +If the molecule perceived from the 3D coordinates is not isomorphic with the input 2D graph, ARC will not spawn any additional jobs for the species, and will not use it further -(for thermo and/or rates calculations). However, sometimes the perception algorithm doesn't work as expected (e.g., -issues with charged species and triplets are known). To continue spawning jobs for all species in an ARC +(for thermo and/or rates calculations). However, sometimes the perception algorithm does not work as expected (e.g., +issues with charged species and triplets are known). To continue spawning jobs for all non-TS species in an ARC project, pass ``True`` to the ``allow_nonisomorphic_2d`` argument (it is ``False`` by default). +Transition states are handled differently. ARC does not enforce 2D-graph isomorphism for TS species, since TS +connectivity and bond orders are not uniquely defined. TS validation is instead based on TS-specific criteria such as +the imaginary frequency, normal mode displacement checks, IRC calculations, and energetic consistency. + .. _directory: