Use the converter.check_isomorphic

instead of isIsomorphic()

arc/scheduler.py

@@ -26,11 +26,11 @@
 import arc.rmgdb as rmgdb
 from arc import plotter
 from arc import parser
-from arc.species.converter import get_xyz_string, molecules_from_xyz
 from arc.job.job import Job
 from arc.arc_exceptions import SpeciesError, SchedulerError
 from arc.job.ssh import SSH_Client
 from arc.species.species import ARCSpecies, TSGuess, determine_rotor_symmetry
+from arc.species.converter import get_xyz_string, molecules_from_xyz, check_isomorphism
 from arc.ts.atst import autotst
 from arc.settings import rotor_scan_resolution, inconsistency_ab, inconsistency_az, maximum_barrier
 
@@ -717,10 +717,7 @@ def determine_most_stable_conformer(self, label):
                 for i, xyz in enumerate(xyzs):
                     _, b_mol = molecules_from_xyz(xyz)
                     if b_mol is not None:
-                        # make copies of the molecules, since isIsomorphic() changes atom orders
-                        b_mol_copy, mol_copy = b_mol.copy(deep=True), self.species_dict[label].mol.copy(deep=True)
-                        match = mol_copy.isIsomorphic(b_mol_copy)
-                        if match:
+                        if check_isomorphism(self.species_dict[label].mol, b_mol):
                             if i == 0:
                                 logging.info('Most stable conformer for species {0} was found to be isomorphic '
                                              'with the 2D graph representation {1}\n'.format(label, b_mol.toSMILES()))
@@ -732,7 +729,8 @@ def determine_most_stable_conformer(self, label):
                                              'with the 2D graph representation {1}. This conformer is {2} kJ/mol '
                                              'above the most stable one (which is not isomorphic). Using the '
                                              'isomorphic conformer for further geometry optimization.'.format(
-                                              label, mol_copy.toSMILES(), (energies[i] - energies[0]) * 2625.50))
+                                              label, self.species_dict[label].mol.toSMILES(),
+                                              (energies[i] - energies[0]) * 2625.50))
                                 conformer_xyz = xyz
                                 self.output[label]['status'] += 'passed isomorphism check but not for the most stable conformer; '
                             break

arc/species/species.py

@@ -27,7 +27,7 @@
 from arc.settings import arc_path, default_ts_methods, valid_chars, minimum_barrier
 from arc.parser import parse_xyz_from_file
 from arc.species.converter import get_xyz_string, get_xyz_matrix, rdkit_conf_from_mol, check_xyz,\
-    molecules_from_xyz, rmg_mol_from_inchi, order_atoms_in_mol_list
+    molecules_from_xyz, rmg_mol_from_inchi, order_atoms_in_mol_list, check_isomorphism
 from arc.ts import atst
 
 ##################################################################
@@ -811,19 +811,14 @@ def mol_from_xyz(self, xyz=None):
             xyz = self.initial_xyz
         original_mol = self.mol.copy(deep=True)
         _, self.mol = molecules_from_xyz(xyz)
-        if not original_mol.isIsomorphic(self.mol):
+        if not check_isomorphism(original_mol, self.mol):
             raise InputError('XYZ and the 2D graph representation of the Molecule are not isomorphic.\n'
                              'Got xyz:\n{0}\n\nwhich corresponds to {1}\n{2}\n\nand: {3}\n{4}'.format(
                                xyz, self.mol.toSMILES(), self.mol.toAdjacencyList(),
                                original_mol.toSMILES(), original_mol.toAdjacencyList()))
         if self.mol_list is None:
             self.mol.assignAtomIDs()
             self.mol_list = self.mol.generate_resonance_structures(keep_isomorphic=False, filter_structures=True)
-        # The isIsomorphic and generate_resonance_structures methods changes atom order, reorder atoms in self.mol:
-        id_mol = self.mol.copy(deep=True)
-        _, self.mol = molecules_from_xyz(xyz)
-        for i, atom in enumerate(self.mol.atoms):
-            atom.id = id_mol.atoms[i].id
         order_atoms_in_mol_list(ref_mol=self.mol, mol_list=self.mol_list)