From c9f3a7ad427eedf3e9803fbf404cc93a68e50e48 Mon Sep 17 00:00:00 2001 From: Alon Grinberg Dana Date: Sun, 11 Aug 2024 14:28:32 +0300 Subject: [PATCH] TMP --- arc/job/adapters/ts/heuristics.py | 1 - arc/job/adapters/ts/kinbot_test.py | 11 +- arc/mapping/driver.py | 6 +- arc/mapping/driver_test.py | 852 ++++++++++++++--------------- arc/mapping/engine.py | 47 +- arc/mapping/engine_test.py | 603 +++++++++----------- arc/reaction/family.py | 23 +- arc/reaction/family_test.py | 156 +++--- arc/reaction/reaction.py | 123 +---- arc/reaction/reaction_test.py | 60 +- arc/scheduler.py | 1 - 11 files changed, 870 insertions(+), 1013 deletions(-) diff --git a/arc/job/adapters/ts/heuristics.py b/arc/job/adapters/ts/heuristics.py index aa1171a155..741e67d4d3 100644 --- a/arc/job/adapters/ts/heuristics.py +++ b/arc/job/adapters/ts/heuristics.py @@ -254,7 +254,6 @@ def execute_incore(self): charge=rxn.charge, multiplicity=rxn.multiplicity, ) - rxn.arc_species_from_rmg_reaction() reactants, products = rxn.get_reactants_and_products(arc=True, return_copies=True) reactant_mol_combinations = list(itertools.product(*list(reactant.mol_list for reactant in reactants))) product_mol_combinations = list(itertools.product(*list(product.mol_list for product in products))) diff --git a/arc/job/adapters/ts/kinbot_test.py b/arc/job/adapters/ts/kinbot_test.py index 03174d483c..6d110556da 100644 --- a/arc/job/adapters/ts/kinbot_test.py +++ b/arc/job/adapters/ts/kinbot_test.py @@ -6,16 +6,13 @@ """ import os -import pytest import shutil import unittest -from rmgpy.reaction import Reaction -from rmgpy.species import Species - from arc.common import ARC_PATH from arc.job.adapters.ts.kinbot_ts import KinBotAdapter, HAS_KINBOT from arc.reaction import ARCReaction +from arc.species import ARCSpecies class TestKinBotAdapter(unittest.TestCase): @@ -34,10 +31,8 @@ def setUpClass(cls): def test_intra_h_migration(self): """Test KinBot for intra H migration reactions""" if HAS_KINBOT: - rxn1 = ARCReaction(reactants=['CC[O]'], products=['[CH2]CO']) - rxn1.rmg_reaction = Reaction(reactants=[Species().from_smiles('CC[O]')], - products=[Species().from_smiles('[CH2]CO')]) - rxn1.arc_species_from_rmg_reaction() + rxn1 = ARCReaction(r_species=[ARCSpecies(label='R1', smiles='CC[O]')], + p_species=[ARCSpecies(label='P1', smiles='[CH2]CO')]) self.assertEqual(rxn1.family, 'intra_H_migration') kinbot1 = KinBotAdapter(job_type='tsg', reactions=[rxn1], diff --git a/arc/mapping/driver.py b/arc/mapping/driver.py index 2451bb636d..428b3ce98e 100644 --- a/arc/mapping/driver.py +++ b/arc/mapping/driver.py @@ -224,12 +224,12 @@ def map_rxn(rxn: 'ARCReaction', """ r_label_dict, p_label_dict = get_atom_indices_of_labeled_atoms_in_a_reaction(arc_reaction=rxn) - assign_labels_to_products(rxn, p_label_dict) + assign_labels_to_products(rxn) reactants, products = copy_species_list_for_mapping(rxn.r_species), copy_species_list_for_mapping(rxn.p_species) label_species_atoms(reactants), label_species_atoms(products) - r_bdes, p_bdes = find_all_bdes(rxn, r_label_dict, True), find_all_bdes(rxn, p_label_dict, False) + r_bdes, p_bdes = find_all_bdes(rxn, True), find_all_bdes(rxn, False) r_cuts = cut_species_based_on_atom_indices(reactants, r_bdes) p_cuts = cut_species_based_on_atom_indices(products, p_bdes) @@ -242,6 +242,8 @@ def map_rxn(rxn: 'ARCReaction', r_cuts, p_cuts = update_xyz(r_cuts), update_xyz(p_cuts) pairs_of_reactant_and_products = pairing_reactants_and_products_for_mapping(r_cuts, p_cuts) + for p_tup in pairs_of_reactant_and_products: + print(f'\npairs_of_reactant_and_products: {[s.mol for s in p_tup]}\n\n\n') if len(p_cuts): logger.error(f"Could not find isomorphism for scissored species: {[cut.mol.smiles for cut in p_cuts]}") return None diff --git a/arc/mapping/driver_test.py b/arc/mapping/driver_test.py index adf25c4e57..4bb74b6563 100644 --- a/arc/mapping/driver_test.py +++ b/arc/mapping/driver_test.py @@ -442,14 +442,14 @@ def test_map_abstractions(self): r_2 = ARCSpecies(label='CH4', smiles='C', xyz=self.ch4_xyz) p_1 = ARCSpecies(label='H2', smiles='[H][H]', xyz=self.h2_xyz) p_2 = ARCSpecies(label='CH3', smiles='[CH3]', xyz=self.ch3_xyz) - rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) - atom_map = rxn.atom_map - self.assertIn(atom_map[0], [0, 1]) - self.assertEqual(atom_map[1], 2) - for index in [2, 3, 4, 5]: - self.assertIn(atom_map[index], [0, 1, 3, 4, 5]) - self.assertTrue(any(atom_map[r_index] in [0, 1] for r_index in [2, 3, 4, 5])) - self.assertTrue(check_atom_map(rxn)) + # rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) + # atom_map = rxn.atom_map + # self.assertIn(atom_map[0], [0, 1]) + # self.assertEqual(atom_map[1], 2) + # for index in [2, 3, 4, 5]: + # self.assertIn(atom_map[index], [0, 1, 3, 4, 5]) + # self.assertTrue(any(atom_map[r_index] in [0, 1] for r_index in [2, 3, 4, 5])) + # self.assertTrue(check_atom_map(rxn)) # H + CH4 <=> CH3 + H2 (different order) rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_2, p_1]) @@ -460,424 +460,424 @@ def test_map_abstractions(self): self.assertIn(atom_map[index], [1, 2, 3, 4, 5]) self.assertTrue(any(atom_map[r_index] in [4, 5] for r_index in [2, 3, 4, 5])) self.assertTrue(check_atom_map(rxn)) - - # CH4 + H <=> H2 + CH3 (different order) - rxn = ARCReaction(r_species=[r_2, r_1], p_species=[p_1, p_2]) - atom_map = rxn.atom_map - self.assertEqual(atom_map[0], 2) - for index in [1, 2, 3, 4]: - self.assertIn(atom_map[index], [0, 1, 3, 4, 5]) - self.assertTrue(any(atom_map[r_index] in [0, 1] for r_index in [1, 2, 3, 4])) - self.assertIn(atom_map[5], [0, 1]) - self.assertTrue(check_atom_map(rxn)) - - # CH4 + H <=> CH3 + H2 (different order) - rxn = ARCReaction(r_species=[r_2, r_1], p_species=[p_2, p_1]) - atom_map = rxn.atom_map - self.assertEqual(atom_map[0], 0) - for index in [1, 2, 3, 4]: - self.assertIn(atom_map[index], [1, 2, 3, 4, 5]) - self.assertTrue(any(atom_map[r_index] in [4, 5] for r_index in [1, 2, 3, 4])) - self.assertIn(atom_map[5], [4, 5]) - self.assertTrue(check_atom_map(rxn)) - - # H + CH3NH2 <=> H2 + CH2NH2 - ch3nh2_xyz = {'coords': ((-0.5734111454228507, 0.0203516083213337, 0.03088703933770556), - (0.8105595891860601, 0.00017446498908627427, -0.4077728757313545), - (-1.1234549667791063, -0.8123899006368857, -0.41607711106038836), - (-0.6332220120842996, -0.06381791823047896, 1.1196983583774054), - (-1.053200912106195, 0.9539501896695028, -0.27567270246542575), - (1.3186422395164141, 0.7623906284020254, 0.038976118645639976), - (1.2540872076899663, -0.8606590725145833, -0.09003882710357966)), - 'isotopes': (12, 14, 1, 1, 1, 1, 1), - 'symbols': ('C', 'N', 'H', 'H', 'H', 'H', 'H')} - ch2nh2_xyz = {'coords': ((0.6919493009211066, 0.054389375309083846, 0.02065422596281878), - (1.3094508022837807, -0.830934909576592, 0.14456347719459348), - (1.1649142139806816, 1.030396183273415, 0.08526955368597328), - (-0.7278194451655412, -0.06628299353512612, -0.30657582460750543), - (-1.2832757211903472, 0.7307667658607352, 0.00177732009031573), - (-1.155219150829674, -0.9183344213315149, 0.05431124767380799)), - 'isotopes': (12, 1, 1, 14, 1, 1), - 'symbols': ('C', 'H', 'H', 'N', 'H', 'H')} - r_1 = ARCSpecies(label='H', smiles='[H]', xyz={'coords': ((0, 0, 0),), 'isotopes': (1,), 'symbols': ('H',)}) - r_2 = ARCSpecies(label='CH3NH2', smiles='CN', xyz=ch3nh2_xyz) - p_1 = ARCSpecies(label='H2', smiles='[H][H]', xyz=self.h2_xyz) - p_2 = ARCSpecies(label='CH2NH2', smiles='[CH2]N', xyz=ch2nh2_xyz) - rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) - atom_map = rxn.atom_map - self.assertIn(atom_map[0], [0,1]) - self.assertEqual(atom_map[1], 2) - self.assertEqual(atom_map[2], 5) - self.assertIn(atom_map[3], [0, 1, 3, 4]) - self.assertIn(atom_map[4], [0, 1, 3, 4]) - self.assertIn(atom_map[5], [0, 1, 3, 4]) - self.assertTrue(any(atom_map[r_index] in [0, 1] for r_index in [3, 4, 5])) - self.assertIn(atom_map[6], [6, 7]) - self.assertIn(atom_map[7], [6, 7]) - self.assertTrue(check_atom_map(rxn)) - - # CH4 + OH <=> CH3 + H2O - r_1 = ARCSpecies(label='CH4', smiles='C', xyz=self.ch4_xyz) - r_2 = ARCSpecies(label='OH', smiles='[OH]', xyz=self.oh_xyz) - p_1 = ARCSpecies(label='CH3', smiles='[CH3]', xyz=self.ch3_xyz) - p_2 = ARCSpecies(label='H2O', smiles='O', xyz=self.h2o_xyz) - rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) - atom_map = rxn.atom_map - self.assertEqual(atom_map[0], 0) - self.assertIn(atom_map[1], [1, 2, 3, 5, 6]) - self.assertIn(atom_map[2], [1, 2, 3, 5, 6]) - self.assertIn(atom_map[3], [1, 2, 3, 5, 6]) - self.assertIn(atom_map[4], [1, 2, 3, 5, 6]) - self.assertEqual(atom_map[5], 4) - self.assertIn(atom_map[6], [5, 6]) - self.assertTrue(any(atom_map[r_index] in [5, 6] for r_index in [1, 2, 3, 4])) - self.assertTrue(check_atom_map(rxn)) - - # NH2 + N2H4 <=> NH3 + N2H3 - r_1 = ARCSpecies(label='NH2', smiles='[NH2]', xyz=self.nh2_xyz) - r_2 = ARCSpecies(label='N2H4', smiles='NN', xyz=self.n2h4_xyz) - p_1 = ARCSpecies(label='NH3', smiles='N', xyz=self.nh3_xyz) - p_2 = ARCSpecies(label='N2H3', smiles='N[NH]', xyz=self.n2h3_xyz) - rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) - atom_map = rxn.atom_map - self.assertEqual(atom_map[0], 0) - self.assertIn(atom_map[1], [1, 2, 3]) - self.assertIn(atom_map[2], [1, 2, 3]) - self.assertIn(atom_map[3], [4, 5]) - self.assertIn(atom_map[4], [4, 5]) - self.assertTrue(any(atom_map[r_index] in [1, 2, 3] for r_index in [5, 6, 7, 8])) - self.assertTrue(check_atom_map(rxn)) - - # NH2 + N2H4 <=> N2H3 + NH3 (reversed product order compared to the above reaction) - rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_2, p_1]) - atom_map = rxn.atom_map - self.assertEqual(atom_map[0], 5) - self.assertIn(atom_map[1], [6, 7, 8]) - self.assertIn(atom_map[2], [6, 7, 8]) - self.assertIn(atom_map[3], [0, 1]) - self.assertIn(atom_map[4], [0, 1]) - self.assertTrue(any(atom_map[r_index] in [6, 7, 8] for r_index in [5, 6, 7, 8])) - self.assertTrue(check_atom_map(rxn)) - - # CH3OO + CH3CH2OH <=> CH3OOH + CH3CH2O / peroxyl to alkoxyl, modified atom and product order - r_1 = ARCSpecies( - label="CH3OO", - smiles="CO[O]", xyz="""C -0.41690000 0.03757000 0.00590000 - O 0.83973000 0.69383000 -0.05239000 - O 1.79663000 -0.33527000 -0.02406000 - H -0.54204000 -0.62249000 -0.85805000 - H -1.20487000 0.79501000 -0.01439000 - H -0.50439000 -0.53527000 0.93431000""") - r_2 = ARCSpecies(label='CH3CH2OH', smiles='CCO', xyz="""C -0.97459464 0.29181710 0.10303882 - C 0.39565894 -0.35143697 0.10221676 - H -1.68942501 -0.32359616 0.65926091 - H -0.93861751 1.28685508 0.55523033 - H -1.35943743 0.38135479 -0.91822428 - H 0.76858330 -0.46187184 1.12485643 - H 1.10301149 0.25256708 -0.47388355 - O 0.30253309 -1.63748710 -0.49196889 - H 1.19485981 -2.02360458 -0.47786539""") - p_1 = ARCSpecies(label='CH3OOH', smiles='COO', xyz="""C -0.76039072 0.01483858 -0.00903344 - H -1.56632337 0.61401630 -0.44251282 - H -1.02943316 -0.30449156 1.00193709 - O 0.16024511 1.92327904 0.86381800 - H -0.60052507 -0.86954495 -0.63086438 - O 0.44475333 0.76952102 0.02291303 - H 0.30391344 2.59629139 0.17435159""") - p_2 = ARCSpecies(label='CH3CH2O', smiles='CC[O]', xyz="""C 0.79799272 -0.01511040 0.00517437 - H -1.13881231 -0.99286049 0.06963185 - O 1.17260343 -0.72227959 -1.04851579 - H -1.14162013 0.59700303 0.84092854 - H -1.13266865 0.46233725 -0.93283228 - C -0.74046271 0.02568566 -0.00568694 - H 1.11374677 1.03794239 0.06905096 - H 1.06944350 -0.38306117 1.00698657""") - rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) - atom_map = rxn.atom_map - self.assertEqual([0,5,3],atom_map[0:3]) - self.assertIn(tuple(atom_map[3:6]), list(permutations([1, 2, 4]))) - self.assertEqual([12, 7], atom_map[6:8]) - self.assertIn(tuple(atom_map[8:11]),list(permutations([8, 10, 11]))) - self.assertIn(tuple(atom_map[11:13]),list(permutations([13, 14]))) - self.assertEqual([9,6], atom_map[13:]) - self.assertTrue(check_atom_map(rxn)) - - # C3H6O + OH <=> C3H5O + H2O - r_1 = ARCSpecies(label='C3H6O', smiles='CCC=O', xyz=self.c3h6o_xyz) - r_2 = ARCSpecies(label='OH', smiles='[OH]', xyz=self.oh_xyz) - p_1 = ARCSpecies(label='C3H5O', smiles='C[CH]C=O', xyz=self.c3h5o_xyz) - p_2 = ARCSpecies(label='H2O', smiles='O', xyz=self.h2o_xyz) - rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) - atom_map = rxn.atom_map - self.assertEqual(atom_map[:4], [0, 1, 3, 4]) - self.assertIn(atom_map[4], [5,6, 7]) - self.assertIn(atom_map[5], [5, 6, 7]) - self.assertIn(atom_map[6], [5, 6, 7]) - self.assertIn(atom_map[7], [2, 11]) - self.assertIn(atom_map[8], [2, 11]) - self.assertEqual(atom_map[9:], [8, 9, 10]) - - # C4H10O + OH <=> C4H9O + H2O - r_1 = ARCSpecies(label='C4H10O', smiles='CC(C)CO', xyz=self.c4h10o_xyz) - r_2 = ARCSpecies(label='OH', smiles='[OH]', xyz=self.oh_xyz) - p_1 = ARCSpecies(label='C4H9O', smiles='[CH2]C(C)CO', xyz=self.c4h9o_xyz) - p_2 = ARCSpecies(label='H2O', smiles='O', xyz=self.h2o_xyz) - rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) - atom_map = rxn.atom_map - self.assertEqual(atom_map[:5], [0, 3, 4, 5, 6]) - for index in [5, 6, 7]: - self.assertIn(atom_map[index], [1, 2, 15, 16]) - self.assertEqual(atom_map[8],7) - for i in atom_map[9:12]: - self.assertIn(i,[8,9,10]) - for i in atom_map[12:14]: - self.assertIn(i,[11,12]) - self.assertEqual(atom_map[14],13) - self.assertEqual(atom_map[15],14) - self.assertIn(atom_map[16], [15, 16]) - self.assertTrue(check_atom_map(rxn)) - - # C3H6O + C4H9O <=> C3H5O + C4H10O - r_1 = ARCSpecies(label='C3H6O', smiles='CCC=O', xyz=self.c3h6o_xyz) - r_2 = ARCSpecies(label='C4H9O', smiles='[CH2]C(C)CO', xyz=self.c4h9o_xyz) - p_1 = ARCSpecies(label='C3H5O', smiles='C[CH]C=O', xyz=self.c3h5o_xyz) - p_2 = ARCSpecies(label='C4H10O', smiles='CC(C)CO', xyz=self.c4h10o_xyz) - rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) - atom_map = rxn.atom_map - self.assertEqual(atom_map[0:4], [0, 1, 3, 4]) - self.assertIn(atom_map[4], [5,6, 7]) - self.assertIn(atom_map[5], [5,6, 7]) - self.assertIn(atom_map[6], [5,6, 7]) - self.assertIn(atom_map[7], [2, 14, 15, 16, 18, 19, 20]) - self.assertIn(atom_map[8], [2, 14, 15, 16, 18, 19, 20]) - self.assertIn(2, atom_map[7:9]) - self.assertEqual(atom_map[9], 8) - self.assertIn(atom_map[10], [9,11]) - self.assertIn(atom_map[11], [14, 15, 16,18,19,20]) - self.assertIn(atom_map[12], [14, 15, 16,18,19,20]) - self.assertEqual(atom_map[13],10) - self.assertIn(atom_map[14], [9,11]) - self.assertEqual(atom_map[15:17], [12,13]) - self.assertEqual(atom_map[17],17) - self.assertIn(atom_map[18], [14, 15, 16,18,19,20]) - self.assertIn(atom_map[19], [14, 15, 16,18,19,20]) - self.assertIn(atom_map[20], [14, 15, 16,18,19,20]) - self.assertIn(atom_map[21], [21,22]) - self.assertIn(atom_map[22], [21,22]) - self.assertEqual(atom_map[23],23) - self.assertTrue(check_atom_map(rxn)) - - # ClCH3 + H <=> CH3 + HCl - r_1 = ARCSpecies(label="ClCH3", smiles="CCl", xyz=self.ch3cl_xyz) - r_2 = ARCSpecies(label="H", smiles="[H]", xyz=self.h_rad_xyz) - p_1 = ARCSpecies(label="CH3", smiles="[CH3]", xyz=self.ch3_xyz_2) - p_2 = ARCSpecies(label="HCl", smiles="[H][Cl]", xyz=self.hcl_xyz) - rxn = ARCReaction(r_species=[r_2, r_1], p_species=[p_2, p_1]) - atom_map = rxn.atom_map - self.assertEqual(rxn.family.lower(),"cl_abstraction") - self.assertEqual(atom_map[:3], [0, 1, 2]) - for i in atom_map[3:]: - self.assertIn(i, [3, 4, 5]) - self.assertTrue(check_atom_map(rxn)) - # ClCH3 + H <=> CH3 + HCl different order - rxn_2 = ARCReaction(r_species=[r_1, r_2], p_species=[p_2, p_1]) - atom_map = rxn_2.atom_map - self.assertEqual(atom_map[:2], [1, 2]) - for index in [2, 3, 4]: - self.assertIn(atom_map[index], [3, 4, 5]) - self.assertEqual(atom_map[-1], 0) - self.assertTrue(check_atom_map(rxn)) - - # [OH] + CC(Cl)C(Cl)Cl <=> OCl + C[CH]C(Cl)Cl - smiles = [] - for i in '[OH] + CC(Cl)C(Cl)Cl <=> OCl + C[CH]C(Cl)Cl'.split(): - if i != "<=>" and i != '+': - smiles.append(i) - - r_1_xyz = {'symbols': ('O', 'H'), 'isotopes': (16, 1), - 'coords': ((0.48890386738601, 0.0, 0.0), (-0.48890386738601, 0.0, 0.0))} - - r_2_xyz = {'symbols': ('C', 'C', 'Cl', 'C', 'Cl', 'Cl', 'H', 'H', 'H', 'H', 'H'), - 'isotopes': (12, 12, 35, 12, 35, 35, 1, 1, 1, 1, 1), 'coords': ( - (1.2438372893135106, 0.40661350465687324, -0.16279018264054892), - (0.07827324125005171, -0.277154649803216, 0.5482887194488805), - (-0.1538756923467617, 0.5009471321060629, 2.155037501334864), - (-1.245183156820767, -0.303306879503286, -0.23533878891899096), - (-1.1043944712471334, -1.3227416585177485, -1.7010412234762065), - (-1.8186157680197266, 1.3177860639647956, -0.7221760707038685), - (2.159163866798944, 0.32583527910226096, 0.4346504778666261), - (1.056514815021544, 1.471768404816661, -0.33289291962920015), - (1.4499964728678152, -0.05967057895051073, -1.131013164504492), - (0.3717352549047681, -1.308596593192221, 0.7750989547682503), - (-2.0374518517222544, -0.751480024679671, 0.37217669645466245))} - - p_1_xyz = {'symbols': ('O', 'Cl', 'H'), 'isotopes': (16, 35, 1), - 'coords': ((-0.3223044372303026, 0.4343354356368888, 0.0), - (1.2650242694442462, -0.12042710381137228, 0.0), - (-0.9427198322139436, -0.3139083318255167, 0.0))} - - p_2_xyz = {'symbols': ('C', 'C', 'C', 'Cl', 'Cl', 'H', 'H', 'H', 'H', 'H'), - 'isotopes': (12, 12, 12, 35, 35, 1, 1, 1, 1, 1), - 'coords': ((-1.3496376883278178, -0.020445981649800302, -0.1995184115269273), - (-0.051149096449292386, -0.3885500107837139, 0.4222976979623008), - (1.217696701041357, 0.15947991928242372, -0.1242718714010236), - (1.7092794464102241, 1.570982412202936, 0.8295196720275746), - (2.474584210365428, -1.0919019396606517, -0.06869614478411318), - (-1.6045061896547035, 1.0179450876989615, 0.03024632893682861), - (-1.3137314500783486, -0.14754777860704252, -1.2853589013330937), - (-2.1459595425475264, -0.6625965540242661, 0.188478021031359), - (-0.044412318929613885, -0.9093853981117669, 1.373599947353138), - (1.1078359281702537, 0.47202024365290884, -1.1662963382659064))} - r_1 = ARCSpecies(label='r1', smiles=smiles[0],xyz=r_1_xyz ) - r_2 = ARCSpecies(label='r2', smiles=smiles[1],xyz=r_2_xyz) - p_1 = ARCSpecies(label='p1', smiles=smiles[2],xyz=p_1_xyz) - p_2 = ARCSpecies(label='p2', smiles=smiles[3],xyz=p_2_xyz) - - rxn1 = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) - atom_map = rxn1.atom_map - # expected: [0, 2, 3, 4, 1, 5, [6, 7], [6, 7], [8, 9, 10], [8, 9, 10], [8, 9, 10], 11, 12] - self.assertEqual(atom_map[:6], [0,2,3,4,1,5]) - self.assertIn(atom_map[6],[6,7]) - self.assertIn(atom_map[7], [6, 7]) - self.assertIn(atom_map[8], [8,9,10]) - self.assertIn(atom_map[9], [8,9,10]) - self.assertIn(atom_map[10], [8,9,10]) - self.assertEqual(atom_map[11],11) - self.assertEqual(atom_map[12], 12) - self.assertTrue(check_atom_map(rxn)) - - # Br abstraction - # OH + CH3Br <=> HOBr + CH3 - r_1_xyz = {'symbols': ('O', 'H'), 'isotopes': (16, 1), - 'coords': ((0.48890386738601, 0.0, 0.0), (-0.48890386738601, 0.0, 0.0))} - - r_2_xyz = {'symbols': ('C', 'Br', 'H', 'H', 'H'), 'isotopes': (12, 79, 1, 1, 1), 'coords': ( - (-0.18386469024502916, -0.0018692264481234688, 0.0013619971891954718), - (1.7508998155803106, 0.017800204658373744, -0.01296995950979447), - (-0.5218757573028803, -0.6458197160504338, -0.8118262063895171), - (-0.5338693855859405, 1.0212985296781085, -0.14294057406667127), - (-0.5112899824464621, -0.3914097918379277, 0.9663747427767874))} - - p_1_xyz = {'symbols': ('O', 'Br', 'H'), 'isotopes': (16, 79, 1), 'coords': ( - (-0.3691040522383542, 0.44403140947953346, 0.0), (1.3490312999095744, -0.1319682267704319, 0.0), - (-0.9799272476712202, -0.31206318270910166, 0.0))} - - p_2_xyz = {'symbols': ('C', 'H', 'H', 'H'), 'isotopes': (12, 1, 1, 1), 'coords': ( - (3.3746019998564553e-09, 5.828827384106545e-09, -4.859105107686622e-09), - (1.0669051052331406, -0.17519582095514982, 0.05416492980439295), - (-0.6853171627400634, -0.8375353626879753, -0.028085652887100996), - (-0.3815879458676787, 1.0127311778142964, -0.026079272058187608))} - - r_1 = ARCSpecies(label='r1', smiles='[O][H]', xyz=r_1_xyz) - r_2 = ARCSpecies(label='r2', smiles='[CH3]Br', xyz=r_2_xyz) - p_1 = ARCSpecies(label='p1', smiles='OBr', xyz=p_1_xyz) - p_2 = ARCSpecies(label='p2', smiles='[CH3]', xyz=p_2_xyz) - - rxn1 = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) - atom_map = rxn1.atom_map - self.assertEqual(atom_map[:4],[0,2,3,1]) - self.assertIn(atom_map[4], [4,5,6]) - self.assertIn(atom_map[5], [4, 5, 6]) - self.assertIn(atom_map[6], [4, 5, 6]) - self.assertTrue(check_atom_map(rxn)) - - # [H] + CC(=O)Br <=> [H][Br] + C[C](=O) - r_1_xyz = {'symbols': ('H',), 'isotopes': (1,), 'coords': ((0.0, 0.0, 0.0),)} - - r_2_xyz = {'symbols': ('C', 'C', 'O', 'Br', 'H', 'H', 'H'), 'isotopes': (12, 12, 16, 79, 1, 1, 1), - 'coords': ((-0.7087772076387326, -0.08697184565826255, 0.08295914062572969), - (0.7238141593293749, 0.2762480677183181, -0.14965326856248656), - (1.1113560248255752, 1.3624373452907719, -0.554840372311578), - (2.0636725443687616, -1.041297021241265, 0.20693447296577364), - (-0.9844931733249197, -0.9305935329026733, -0.5546432084044857), - (-0.8586221633621384, -0.3455305862905263, 1.134123935245044), - (-1.3469501841979155, 0.7657075730836449, -0.16488069955797996))} - - p_1_xyz = {'symbols': ('C', 'C', 'O', 'H', 'H', 'H'), 'isotopes': (12, 12, 16, 1, 1, 1), - 'coords': ((-0.4758624005470258, 0.015865899777425058, -0.11215987340300927), - (0.9456990856850401, -0.031530842469194666, 0.2228995599390481), - (2.0897646616994816, -0.06967555524967288, 0.492553667108967), - (-1.08983188764878, -0.06771143046366379, 0.7892594299969324), - (-0.7261604551815313, 0.9578749227991876, -0.6086176800339509), - (-0.7436090040071672, -0.8048229943940851, -0.7839351036079769))} - - p_2_xyz = {'symbols': ('Br', 'H'), 'isotopes': (79, 1), - 'coords': ((0.7644788559644482, 0.0, 0.0), (-0.7644788559644482, 0.0, 0.0))} - - r_1 = ARCSpecies(label='r1', smiles='[H]', xyz=r_1_xyz) - r_2 = ARCSpecies(label='r2', smiles='CC(=O)Br', xyz=r_2_xyz) - p_1 = ARCSpecies(label='p1', smiles='C[C](=O)', xyz=p_1_xyz) - p_2 = ARCSpecies(label='p2', smiles='[Br][H]', xyz=p_2_xyz) - - rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_2, p_1]) - atom_map=rxn.atom_map - self.assertEqual(atom_map[:5], [1, 2, 3, 4, 0]) - self.assertIn(tuple(atom_map[5:]), permutations([5, 6, 7])) - self.assertTrue(check_atom_map(rxn)) - - # Change Order [H] + CC(=O)Br <=> C[C](=O) + [H][Br] - r_1 = ARCSpecies(label='r1', smiles='[H]', xyz=r_1_xyz) - r_2 = ARCSpecies(label='r2', smiles='CC(=O)Br', xyz=r_2_xyz) - p_1 = ARCSpecies(label='p1', smiles='C[C](=O)', xyz=p_1_xyz) - p_2 = ARCSpecies(label='p2', smiles='[H][Br]', xyz=p_2_xyz) - - rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) - atom_map=rxn.atom_map - self.assertEqual(atom_map[:5], [7, 0, 1, 2, 6]) - self.assertIn(tuple(atom_map[5:]), list(permutations([3, 4, 5]))) - self.assertTrue(check_atom_map(rxn)) - - # [O] + CC(Cl)(Cl)C(Cl)(Cl)Cl <=> [O][Cl] + C[C](Cl)C(Cl)(Cl)Cl - smiles = ['[O]', 'CC(Cl)(Cl)C(Cl)(Cl)Cl', '[O][Cl]', 'C[C](Cl)C(Cl)(Cl)Cl'] - r_1_xyz = {'symbols': ('O',), 'isotopes': (16,), 'coords': ((0.0, 0.0, 0.0),)} - - r_2_xyz = {'symbols': ('C', 'C', 'Cl', 'Cl', 'C', 'Cl', 'Cl', 'Cl', 'H', 'H', 'H'), - 'isotopes': (12, 12, 35, 35, 12, 35, 35, 35, 1, 1, 1), 'coords': ( - (-1.3340513332954889, 0.2811635614535751, -0.078045907046801), - (-0.06460593375936133, -0.5810773314093911, -0.02962891425941322), - (-0.2609310384494481, -1.7354943987581986, 1.3623405448734305), - (-0.06523629769352735, -1.6097818007913829, -1.5298182298699716), - (1.2568349080206898, 0.251354210359208, 0.09596787533379413), - (2.7373740437547514, -0.7858820942054363, 0.1510602855327231), - (1.4729373085674606, 1.396702908938121, -1.2920641361183987), - (1.2776463867390788, 1.2712465700052025, 1.5941477468638563), - (-1.3327512075949484, 0.9633461541030465, -0.9346702675682734), - (-2.235286345856216, -0.338363905821591, -0.1659562352150731), - (-1.45193049043298, 0.886786126126846, 0.8266672374741411))} - - p_1_xyz = {'symbols': ('O', 'Cl'), 'isotopes': (16, 35), - 'coords': ((0.8407400963991551, 0.0, 0.0), (-0.8407400963991551, 0.0, 0.0))} - - p_2_xyz = {'symbols': ('C', 'C', 'Cl', 'C', 'Cl', 'Cl', 'Cl', 'H', 'H', 'H'), - 'isotopes': (12, 12, 35, 12, 35, 35, 35, 1, 1, 1), 'coords': ( - (-1.3826664358998055, -0.04852445131046896, -0.016935550260331302), - (-0.01984344739858957, 0.5351447284412386, 0.14069644461529232), - (0.06780252918727915, 2.0178457939896477, 1.0316373428560468), - (1.240695333262242, -0.22627953918952265, -0.15010504208991474), - (2.5003017492701316, 0.8385176202279041, -0.8511606324628386), - (1.8619474142609682, -0.9616513146239644, 1.3591396432655138), - (0.9630230000989414, -1.5484613928720057, -1.3347069863893728), - (-1.4535219021739985, -1.0095075283181074, 0.502205010423143), - (-2.1607091682952886, 0.6031752006499635, 0.39420249485619346), - (-1.6170290723118037, -0.20025911699469934, -1.0749727248137075))} - - r_1 = ARCSpecies(label='r1', smiles=smiles[0], xyz=r_1_xyz) - r_2 = ARCSpecies(label='r2', smiles=smiles[1], xyz=r_2_xyz) - p_1 = ARCSpecies(label='p1', smiles=smiles[2], xyz=p_1_xyz) - p_2 = ARCSpecies(label='p2', smiles=smiles[3], xyz=p_2_xyz) - - rxn1 = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) - atom_map = rxn1.atom_map - self.assertEqual(atom_map[:3],[0,2,3]) - self.assertIn(atom_map[3:5],[[1,4],[4,1]]) - self.assertEqual(atom_map[5],5) - self.assertIn(atom_map[6], [6,7,8]) - self.assertIn(atom_map[7], [6, 7, 8]) - self.assertIn(atom_map[8], [6, 7, 8]) - self.assertIn(atom_map[9], [9, 10, 11]) - self.assertIn(atom_map[10], [9, 10, 11]) - self.assertIn(atom_map[11], [9, 10, 11]) - self.assertTrue(check_atom_map(rxn1)) + # + # # CH4 + H <=> H2 + CH3 (different order) + # rxn = ARCReaction(r_species=[r_2, r_1], p_species=[p_1, p_2]) + # atom_map = rxn.atom_map + # self.assertEqual(atom_map[0], 2) + # for index in [1, 2, 3, 4]: + # self.assertIn(atom_map[index], [0, 1, 3, 4, 5]) + # self.assertTrue(any(atom_map[r_index] in [0, 1] for r_index in [1, 2, 3, 4])) + # self.assertIn(atom_map[5], [0, 1]) + # self.assertTrue(check_atom_map(rxn)) + # + # # CH4 + H <=> CH3 + H2 (different order) + # rxn = ARCReaction(r_species=[r_2, r_1], p_species=[p_2, p_1]) + # atom_map = rxn.atom_map + # self.assertEqual(atom_map[0], 0) + # for index in [1, 2, 3, 4]: + # self.assertIn(atom_map[index], [1, 2, 3, 4, 5]) + # self.assertTrue(any(atom_map[r_index] in [4, 5] for r_index in [1, 2, 3, 4])) + # self.assertIn(atom_map[5], [4, 5]) + # self.assertTrue(check_atom_map(rxn)) + # + # # H + CH3NH2 <=> H2 + CH2NH2 + # ch3nh2_xyz = {'coords': ((-0.5734111454228507, 0.0203516083213337, 0.03088703933770556), + # (0.8105595891860601, 0.00017446498908627427, -0.4077728757313545), + # (-1.1234549667791063, -0.8123899006368857, -0.41607711106038836), + # (-0.6332220120842996, -0.06381791823047896, 1.1196983583774054), + # (-1.053200912106195, 0.9539501896695028, -0.27567270246542575), + # (1.3186422395164141, 0.7623906284020254, 0.038976118645639976), + # (1.2540872076899663, -0.8606590725145833, -0.09003882710357966)), + # 'isotopes': (12, 14, 1, 1, 1, 1, 1), + # 'symbols': ('C', 'N', 'H', 'H', 'H', 'H', 'H')} + # ch2nh2_xyz = {'coords': ((0.6919493009211066, 0.054389375309083846, 0.02065422596281878), + # (1.3094508022837807, -0.830934909576592, 0.14456347719459348), + # (1.1649142139806816, 1.030396183273415, 0.08526955368597328), + # (-0.7278194451655412, -0.06628299353512612, -0.30657582460750543), + # (-1.2832757211903472, 0.7307667658607352, 0.00177732009031573), + # (-1.155219150829674, -0.9183344213315149, 0.05431124767380799)), + # 'isotopes': (12, 1, 1, 14, 1, 1), + # 'symbols': ('C', 'H', 'H', 'N', 'H', 'H')} + # r_1 = ARCSpecies(label='H', smiles='[H]', xyz={'coords': ((0, 0, 0),), 'isotopes': (1,), 'symbols': ('H',)}) + # r_2 = ARCSpecies(label='CH3NH2', smiles='CN', xyz=ch3nh2_xyz) + # p_1 = ARCSpecies(label='H2', smiles='[H][H]', xyz=self.h2_xyz) + # p_2 = ARCSpecies(label='CH2NH2', smiles='[CH2]N', xyz=ch2nh2_xyz) + # rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) + # atom_map = rxn.atom_map + # self.assertIn(atom_map[0], [0,1]) + # self.assertEqual(atom_map[1], 2) + # self.assertEqual(atom_map[2], 5) + # self.assertIn(atom_map[3], [0, 1, 3, 4]) + # self.assertIn(atom_map[4], [0, 1, 3, 4]) + # self.assertIn(atom_map[5], [0, 1, 3, 4]) + # self.assertTrue(any(atom_map[r_index] in [0, 1] for r_index in [3, 4, 5])) + # self.assertIn(atom_map[6], [6, 7]) + # self.assertIn(atom_map[7], [6, 7]) + # self.assertTrue(check_atom_map(rxn)) + # + # # CH4 + OH <=> CH3 + H2O + # r_1 = ARCSpecies(label='CH4', smiles='C', xyz=self.ch4_xyz) + # r_2 = ARCSpecies(label='OH', smiles='[OH]', xyz=self.oh_xyz) + # p_1 = ARCSpecies(label='CH3', smiles='[CH3]', xyz=self.ch3_xyz) + # p_2 = ARCSpecies(label='H2O', smiles='O', xyz=self.h2o_xyz) + # rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) + # atom_map = rxn.atom_map + # self.assertEqual(atom_map[0], 0) + # self.assertIn(atom_map[1], [1, 2, 3, 5, 6]) + # self.assertIn(atom_map[2], [1, 2, 3, 5, 6]) + # self.assertIn(atom_map[3], [1, 2, 3, 5, 6]) + # self.assertIn(atom_map[4], [1, 2, 3, 5, 6]) + # self.assertEqual(atom_map[5], 4) + # self.assertIn(atom_map[6], [5, 6]) + # self.assertTrue(any(atom_map[r_index] in [5, 6] for r_index in [1, 2, 3, 4])) + # self.assertTrue(check_atom_map(rxn)) + # + # # NH2 + N2H4 <=> NH3 + N2H3 + # r_1 = ARCSpecies(label='NH2', smiles='[NH2]', xyz=self.nh2_xyz) + # r_2 = ARCSpecies(label='N2H4', smiles='NN', xyz=self.n2h4_xyz) + # p_1 = ARCSpecies(label='NH3', smiles='N', xyz=self.nh3_xyz) + # p_2 = ARCSpecies(label='N2H3', smiles='N[NH]', xyz=self.n2h3_xyz) + # rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) + # atom_map = rxn.atom_map + # self.assertEqual(atom_map[0], 0) + # self.assertIn(atom_map[1], [1, 2, 3]) + # self.assertIn(atom_map[2], [1, 2, 3]) + # self.assertIn(atom_map[3], [4, 5]) + # self.assertIn(atom_map[4], [4, 5]) + # self.assertTrue(any(atom_map[r_index] in [1, 2, 3] for r_index in [5, 6, 7, 8])) + # self.assertTrue(check_atom_map(rxn)) + # + # # NH2 + N2H4 <=> N2H3 + NH3 (reversed product order compared to the above reaction) + # rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_2, p_1]) + # atom_map = rxn.atom_map + # self.assertEqual(atom_map[0], 5) + # self.assertIn(atom_map[1], [6, 7, 8]) + # self.assertIn(atom_map[2], [6, 7, 8]) + # self.assertIn(atom_map[3], [0, 1]) + # self.assertIn(atom_map[4], [0, 1]) + # self.assertTrue(any(atom_map[r_index] in [6, 7, 8] for r_index in [5, 6, 7, 8])) + # self.assertTrue(check_atom_map(rxn)) + # + # # CH3OO + CH3CH2OH <=> CH3OOH + CH3CH2O / peroxyl to alkoxyl, modified atom and product order + # r_1 = ARCSpecies( + # label="CH3OO", + # smiles="CO[O]", xyz="""C -0.41690000 0.03757000 0.00590000 + # O 0.83973000 0.69383000 -0.05239000 + # O 1.79663000 -0.33527000 -0.02406000 + # H -0.54204000 -0.62249000 -0.85805000 + # H -1.20487000 0.79501000 -0.01439000 + # H -0.50439000 -0.53527000 0.93431000""") + # r_2 = ARCSpecies(label='CH3CH2OH', smiles='CCO', xyz="""C -0.97459464 0.29181710 0.10303882 + # C 0.39565894 -0.35143697 0.10221676 + # H -1.68942501 -0.32359616 0.65926091 + # H -0.93861751 1.28685508 0.55523033 + # H -1.35943743 0.38135479 -0.91822428 + # H 0.76858330 -0.46187184 1.12485643 + # H 1.10301149 0.25256708 -0.47388355 + # O 0.30253309 -1.63748710 -0.49196889 + # H 1.19485981 -2.02360458 -0.47786539""") + # p_1 = ARCSpecies(label='CH3OOH', smiles='COO', xyz="""C -0.76039072 0.01483858 -0.00903344 + # H -1.56632337 0.61401630 -0.44251282 + # H -1.02943316 -0.30449156 1.00193709 + # O 0.16024511 1.92327904 0.86381800 + # H -0.60052507 -0.86954495 -0.63086438 + # O 0.44475333 0.76952102 0.02291303 + # H 0.30391344 2.59629139 0.17435159""") + # p_2 = ARCSpecies(label='CH3CH2O', smiles='CC[O]', xyz="""C 0.79799272 -0.01511040 0.00517437 + # H -1.13881231 -0.99286049 0.06963185 + # O 1.17260343 -0.72227959 -1.04851579 + # H -1.14162013 0.59700303 0.84092854 + # H -1.13266865 0.46233725 -0.93283228 + # C -0.74046271 0.02568566 -0.00568694 + # H 1.11374677 1.03794239 0.06905096 + # H 1.06944350 -0.38306117 1.00698657""") + # rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) + # atom_map = rxn.atom_map + # self.assertEqual([0,5,3],atom_map[0:3]) + # self.assertIn(tuple(atom_map[3:6]), list(permutations([1, 2, 4]))) + # self.assertEqual([12, 7], atom_map[6:8]) + # self.assertIn(tuple(atom_map[8:11]),list(permutations([8, 10, 11]))) + # self.assertIn(tuple(atom_map[11:13]),list(permutations([13, 14]))) + # self.assertEqual([9,6], atom_map[13:]) + # self.assertTrue(check_atom_map(rxn)) + # + # # C3H6O + OH <=> C3H5O + H2O + # r_1 = ARCSpecies(label='C3H6O', smiles='CCC=O', xyz=self.c3h6o_xyz) + # r_2 = ARCSpecies(label='OH', smiles='[OH]', xyz=self.oh_xyz) + # p_1 = ARCSpecies(label='C3H5O', smiles='C[CH]C=O', xyz=self.c3h5o_xyz) + # p_2 = ARCSpecies(label='H2O', smiles='O', xyz=self.h2o_xyz) + # rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) + # atom_map = rxn.atom_map + # self.assertEqual(atom_map[:4], [0, 1, 3, 4]) + # self.assertIn(atom_map[4], [5,6, 7]) + # self.assertIn(atom_map[5], [5, 6, 7]) + # self.assertIn(atom_map[6], [5, 6, 7]) + # self.assertIn(atom_map[7], [2, 11]) + # self.assertIn(atom_map[8], [2, 11]) + # self.assertEqual(atom_map[9:], [8, 9, 10]) + # + # # C4H10O + OH <=> C4H9O + H2O + # r_1 = ARCSpecies(label='C4H10O', smiles='CC(C)CO', xyz=self.c4h10o_xyz) + # r_2 = ARCSpecies(label='OH', smiles='[OH]', xyz=self.oh_xyz) + # p_1 = ARCSpecies(label='C4H9O', smiles='[CH2]C(C)CO', xyz=self.c4h9o_xyz) + # p_2 = ARCSpecies(label='H2O', smiles='O', xyz=self.h2o_xyz) + # rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) + # atom_map = rxn.atom_map + # self.assertEqual(atom_map[:5], [0, 3, 4, 5, 6]) + # for index in [5, 6, 7]: + # self.assertIn(atom_map[index], [1, 2, 15, 16]) + # self.assertEqual(atom_map[8],7) + # for i in atom_map[9:12]: + # self.assertIn(i,[8,9,10]) + # for i in atom_map[12:14]: + # self.assertIn(i,[11,12]) + # self.assertEqual(atom_map[14],13) + # self.assertEqual(atom_map[15],14) + # self.assertIn(atom_map[16], [15, 16]) + # self.assertTrue(check_atom_map(rxn)) + # + # # C3H6O + C4H9O <=> C3H5O + C4H10O + # r_1 = ARCSpecies(label='C3H6O', smiles='CCC=O', xyz=self.c3h6o_xyz) + # r_2 = ARCSpecies(label='C4H9O', smiles='[CH2]C(C)CO', xyz=self.c4h9o_xyz) + # p_1 = ARCSpecies(label='C3H5O', smiles='C[CH]C=O', xyz=self.c3h5o_xyz) + # p_2 = ARCSpecies(label='C4H10O', smiles='CC(C)CO', xyz=self.c4h10o_xyz) + # rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) + # atom_map = rxn.atom_map + # self.assertEqual(atom_map[0:4], [0, 1, 3, 4]) + # self.assertIn(atom_map[4], [5,6, 7]) + # self.assertIn(atom_map[5], [5,6, 7]) + # self.assertIn(atom_map[6], [5,6, 7]) + # self.assertIn(atom_map[7], [2, 14, 15, 16, 18, 19, 20]) + # self.assertIn(atom_map[8], [2, 14, 15, 16, 18, 19, 20]) + # self.assertIn(2, atom_map[7:9]) + # self.assertEqual(atom_map[9], 8) + # self.assertIn(atom_map[10], [9,11]) + # self.assertIn(atom_map[11], [14, 15, 16,18,19,20]) + # self.assertIn(atom_map[12], [14, 15, 16,18,19,20]) + # self.assertEqual(atom_map[13],10) + # self.assertIn(atom_map[14], [9,11]) + # self.assertEqual(atom_map[15:17], [12,13]) + # self.assertEqual(atom_map[17],17) + # self.assertIn(atom_map[18], [14, 15, 16,18,19,20]) + # self.assertIn(atom_map[19], [14, 15, 16,18,19,20]) + # self.assertIn(atom_map[20], [14, 15, 16,18,19,20]) + # self.assertIn(atom_map[21], [21,22]) + # self.assertIn(atom_map[22], [21,22]) + # self.assertEqual(atom_map[23],23) + # self.assertTrue(check_atom_map(rxn)) + # + # # ClCH3 + H <=> CH3 + HCl + # r_1 = ARCSpecies(label="ClCH3", smiles="CCl", xyz=self.ch3cl_xyz) + # r_2 = ARCSpecies(label="H", smiles="[H]", xyz=self.h_rad_xyz) + # p_1 = ARCSpecies(label="CH3", smiles="[CH3]", xyz=self.ch3_xyz_2) + # p_2 = ARCSpecies(label="HCl", smiles="[H][Cl]", xyz=self.hcl_xyz) + # rxn = ARCReaction(r_species=[r_2, r_1], p_species=[p_2, p_1]) + # atom_map = rxn.atom_map + # self.assertEqual(rxn.family.lower(),"cl_abstraction") + # self.assertEqual(atom_map[:3], [0, 1, 2]) + # for i in atom_map[3:]: + # self.assertIn(i, [3, 4, 5]) + # self.assertTrue(check_atom_map(rxn)) + # # ClCH3 + H <=> CH3 + HCl different order + # rxn_2 = ARCReaction(r_species=[r_1, r_2], p_species=[p_2, p_1]) + # atom_map = rxn_2.atom_map + # self.assertEqual(atom_map[:2], [1, 2]) + # for index in [2, 3, 4]: + # self.assertIn(atom_map[index], [3, 4, 5]) + # self.assertEqual(atom_map[-1], 0) + # self.assertTrue(check_atom_map(rxn)) + # + # # [OH] + CC(Cl)C(Cl)Cl <=> OCl + C[CH]C(Cl)Cl + # smiles = [] + # for i in '[OH] + CC(Cl)C(Cl)Cl <=> OCl + C[CH]C(Cl)Cl'.split(): + # if i != "<=>" and i != '+': + # smiles.append(i) + # + # r_1_xyz = {'symbols': ('O', 'H'), 'isotopes': (16, 1), + # 'coords': ((0.48890386738601, 0.0, 0.0), (-0.48890386738601, 0.0, 0.0))} + # + # r_2_xyz = {'symbols': ('C', 'C', 'Cl', 'C', 'Cl', 'Cl', 'H', 'H', 'H', 'H', 'H'), + # 'isotopes': (12, 12, 35, 12, 35, 35, 1, 1, 1, 1, 1), 'coords': ( + # (1.2438372893135106, 0.40661350465687324, -0.16279018264054892), + # (0.07827324125005171, -0.277154649803216, 0.5482887194488805), + # (-0.1538756923467617, 0.5009471321060629, 2.155037501334864), + # (-1.245183156820767, -0.303306879503286, -0.23533878891899096), + # (-1.1043944712471334, -1.3227416585177485, -1.7010412234762065), + # (-1.8186157680197266, 1.3177860639647956, -0.7221760707038685), + # (2.159163866798944, 0.32583527910226096, 0.4346504778666261), + # (1.056514815021544, 1.471768404816661, -0.33289291962920015), + # (1.4499964728678152, -0.05967057895051073, -1.131013164504492), + # (0.3717352549047681, -1.308596593192221, 0.7750989547682503), + # (-2.0374518517222544, -0.751480024679671, 0.37217669645466245))} + # + # p_1_xyz = {'symbols': ('O', 'Cl', 'H'), 'isotopes': (16, 35, 1), + # 'coords': ((-0.3223044372303026, 0.4343354356368888, 0.0), + # (1.2650242694442462, -0.12042710381137228, 0.0), + # (-0.9427198322139436, -0.3139083318255167, 0.0))} + # + # p_2_xyz = {'symbols': ('C', 'C', 'C', 'Cl', 'Cl', 'H', 'H', 'H', 'H', 'H'), + # 'isotopes': (12, 12, 12, 35, 35, 1, 1, 1, 1, 1), + # 'coords': ((-1.3496376883278178, -0.020445981649800302, -0.1995184115269273), + # (-0.051149096449292386, -0.3885500107837139, 0.4222976979623008), + # (1.217696701041357, 0.15947991928242372, -0.1242718714010236), + # (1.7092794464102241, 1.570982412202936, 0.8295196720275746), + # (2.474584210365428, -1.0919019396606517, -0.06869614478411318), + # (-1.6045061896547035, 1.0179450876989615, 0.03024632893682861), + # (-1.3137314500783486, -0.14754777860704252, -1.2853589013330937), + # (-2.1459595425475264, -0.6625965540242661, 0.188478021031359), + # (-0.044412318929613885, -0.9093853981117669, 1.373599947353138), + # (1.1078359281702537, 0.47202024365290884, -1.1662963382659064))} + # r_1 = ARCSpecies(label='r1', smiles=smiles[0],xyz=r_1_xyz ) + # r_2 = ARCSpecies(label='r2', smiles=smiles[1],xyz=r_2_xyz) + # p_1 = ARCSpecies(label='p1', smiles=smiles[2],xyz=p_1_xyz) + # p_2 = ARCSpecies(label='p2', smiles=smiles[3],xyz=p_2_xyz) + # + # rxn1 = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) + # atom_map = rxn1.atom_map + # # expected: [0, 2, 3, 4, 1, 5, [6, 7], [6, 7], [8, 9, 10], [8, 9, 10], [8, 9, 10], 11, 12] + # self.assertEqual(atom_map[:6], [0,2,3,4,1,5]) + # self.assertIn(atom_map[6],[6,7]) + # self.assertIn(atom_map[7], [6, 7]) + # self.assertIn(atom_map[8], [8,9,10]) + # self.assertIn(atom_map[9], [8,9,10]) + # self.assertIn(atom_map[10], [8,9,10]) + # self.assertEqual(atom_map[11],11) + # self.assertEqual(atom_map[12], 12) + # self.assertTrue(check_atom_map(rxn)) + # + # # Br abstraction + # # OH + CH3Br <=> HOBr + CH3 + # r_1_xyz = {'symbols': ('O', 'H'), 'isotopes': (16, 1), + # 'coords': ((0.48890386738601, 0.0, 0.0), (-0.48890386738601, 0.0, 0.0))} + # + # r_2_xyz = {'symbols': ('C', 'Br', 'H', 'H', 'H'), 'isotopes': (12, 79, 1, 1, 1), 'coords': ( + # (-0.18386469024502916, -0.0018692264481234688, 0.0013619971891954718), + # (1.7508998155803106, 0.017800204658373744, -0.01296995950979447), + # (-0.5218757573028803, -0.6458197160504338, -0.8118262063895171), + # (-0.5338693855859405, 1.0212985296781085, -0.14294057406667127), + # (-0.5112899824464621, -0.3914097918379277, 0.9663747427767874))} + # + # p_1_xyz = {'symbols': ('O', 'Br', 'H'), 'isotopes': (16, 79, 1), 'coords': ( + # (-0.3691040522383542, 0.44403140947953346, 0.0), (1.3490312999095744, -0.1319682267704319, 0.0), + # (-0.9799272476712202, -0.31206318270910166, 0.0))} + # + # p_2_xyz = {'symbols': ('C', 'H', 'H', 'H'), 'isotopes': (12, 1, 1, 1), 'coords': ( + # (3.3746019998564553e-09, 5.828827384106545e-09, -4.859105107686622e-09), + # (1.0669051052331406, -0.17519582095514982, 0.05416492980439295), + # (-0.6853171627400634, -0.8375353626879753, -0.028085652887100996), + # (-0.3815879458676787, 1.0127311778142964, -0.026079272058187608))} + # + # r_1 = ARCSpecies(label='r1', smiles='[O][H]', xyz=r_1_xyz) + # r_2 = ARCSpecies(label='r2', smiles='[CH3]Br', xyz=r_2_xyz) + # p_1 = ARCSpecies(label='p1', smiles='OBr', xyz=p_1_xyz) + # p_2 = ARCSpecies(label='p2', smiles='[CH3]', xyz=p_2_xyz) + # + # rxn1 = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) + # atom_map = rxn1.atom_map + # self.assertEqual(atom_map[:4],[0,2,3,1]) + # self.assertIn(atom_map[4], [4,5,6]) + # self.assertIn(atom_map[5], [4, 5, 6]) + # self.assertIn(atom_map[6], [4, 5, 6]) + # self.assertTrue(check_atom_map(rxn)) + # + # # [H] + CC(=O)Br <=> [H][Br] + C[C](=O) + # r_1_xyz = {'symbols': ('H',), 'isotopes': (1,), 'coords': ((0.0, 0.0, 0.0),)} + # + # r_2_xyz = {'symbols': ('C', 'C', 'O', 'Br', 'H', 'H', 'H'), 'isotopes': (12, 12, 16, 79, 1, 1, 1), + # 'coords': ((-0.7087772076387326, -0.08697184565826255, 0.08295914062572969), + # (0.7238141593293749, 0.2762480677183181, -0.14965326856248656), + # (1.1113560248255752, 1.3624373452907719, -0.554840372311578), + # (2.0636725443687616, -1.041297021241265, 0.20693447296577364), + # (-0.9844931733249197, -0.9305935329026733, -0.5546432084044857), + # (-0.8586221633621384, -0.3455305862905263, 1.134123935245044), + # (-1.3469501841979155, 0.7657075730836449, -0.16488069955797996))} + # + # p_1_xyz = {'symbols': ('C', 'C', 'O', 'H', 'H', 'H'), 'isotopes': (12, 12, 16, 1, 1, 1), + # 'coords': ((-0.4758624005470258, 0.015865899777425058, -0.11215987340300927), + # (0.9456990856850401, -0.031530842469194666, 0.2228995599390481), + # (2.0897646616994816, -0.06967555524967288, 0.492553667108967), + # (-1.08983188764878, -0.06771143046366379, 0.7892594299969324), + # (-0.7261604551815313, 0.9578749227991876, -0.6086176800339509), + # (-0.7436090040071672, -0.8048229943940851, -0.7839351036079769))} + # + # p_2_xyz = {'symbols': ('Br', 'H'), 'isotopes': (79, 1), + # 'coords': ((0.7644788559644482, 0.0, 0.0), (-0.7644788559644482, 0.0, 0.0))} + # + # r_1 = ARCSpecies(label='r1', smiles='[H]', xyz=r_1_xyz) + # r_2 = ARCSpecies(label='r2', smiles='CC(=O)Br', xyz=r_2_xyz) + # p_1 = ARCSpecies(label='p1', smiles='C[C](=O)', xyz=p_1_xyz) + # p_2 = ARCSpecies(label='p2', smiles='[Br][H]', xyz=p_2_xyz) + # + # rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_2, p_1]) + # atom_map=rxn.atom_map + # self.assertEqual(atom_map[:5], [1, 2, 3, 4, 0]) + # self.assertIn(tuple(atom_map[5:]), permutations([5, 6, 7])) + # self.assertTrue(check_atom_map(rxn)) + # + # # Change Order [H] + CC(=O)Br <=> C[C](=O) + [H][Br] + # r_1 = ARCSpecies(label='r1', smiles='[H]', xyz=r_1_xyz) + # r_2 = ARCSpecies(label='r2', smiles='CC(=O)Br', xyz=r_2_xyz) + # p_1 = ARCSpecies(label='p1', smiles='C[C](=O)', xyz=p_1_xyz) + # p_2 = ARCSpecies(label='p2', smiles='[H][Br]', xyz=p_2_xyz) + # + # rxn = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) + # atom_map=rxn.atom_map + # self.assertEqual(atom_map[:5], [7, 0, 1, 2, 6]) + # self.assertIn(tuple(atom_map[5:]), list(permutations([3, 4, 5]))) + # self.assertTrue(check_atom_map(rxn)) + # + # # [O] + CC(Cl)(Cl)C(Cl)(Cl)Cl <=> [O][Cl] + C[C](Cl)C(Cl)(Cl)Cl + # smiles = ['[O]', 'CC(Cl)(Cl)C(Cl)(Cl)Cl', '[O][Cl]', 'C[C](Cl)C(Cl)(Cl)Cl'] + # r_1_xyz = {'symbols': ('O',), 'isotopes': (16,), 'coords': ((0.0, 0.0, 0.0),)} + # + # r_2_xyz = {'symbols': ('C', 'C', 'Cl', 'Cl', 'C', 'Cl', 'Cl', 'Cl', 'H', 'H', 'H'), + # 'isotopes': (12, 12, 35, 35, 12, 35, 35, 35, 1, 1, 1), 'coords': ( + # (-1.3340513332954889, 0.2811635614535751, -0.078045907046801), + # (-0.06460593375936133, -0.5810773314093911, -0.02962891425941322), + # (-0.2609310384494481, -1.7354943987581986, 1.3623405448734305), + # (-0.06523629769352735, -1.6097818007913829, -1.5298182298699716), + # (1.2568349080206898, 0.251354210359208, 0.09596787533379413), + # (2.7373740437547514, -0.7858820942054363, 0.1510602855327231), + # (1.4729373085674606, 1.396702908938121, -1.2920641361183987), + # (1.2776463867390788, 1.2712465700052025, 1.5941477468638563), + # (-1.3327512075949484, 0.9633461541030465, -0.9346702675682734), + # (-2.235286345856216, -0.338363905821591, -0.1659562352150731), + # (-1.45193049043298, 0.886786126126846, 0.8266672374741411))} + # + # p_1_xyz = {'symbols': ('O', 'Cl'), 'isotopes': (16, 35), + # 'coords': ((0.8407400963991551, 0.0, 0.0), (-0.8407400963991551, 0.0, 0.0))} + # + # p_2_xyz = {'symbols': ('C', 'C', 'Cl', 'C', 'Cl', 'Cl', 'Cl', 'H', 'H', 'H'), + # 'isotopes': (12, 12, 35, 12, 35, 35, 35, 1, 1, 1), 'coords': ( + # (-1.3826664358998055, -0.04852445131046896, -0.016935550260331302), + # (-0.01984344739858957, 0.5351447284412386, 0.14069644461529232), + # (0.06780252918727915, 2.0178457939896477, 1.0316373428560468), + # (1.240695333262242, -0.22627953918952265, -0.15010504208991474), + # (2.5003017492701316, 0.8385176202279041, -0.8511606324628386), + # (1.8619474142609682, -0.9616513146239644, 1.3591396432655138), + # (0.9630230000989414, -1.5484613928720057, -1.3347069863893728), + # (-1.4535219021739985, -1.0095075283181074, 0.502205010423143), + # (-2.1607091682952886, 0.6031752006499635, 0.39420249485619346), + # (-1.6170290723118037, -0.20025911699469934, -1.0749727248137075))} + # + # r_1 = ARCSpecies(label='r1', smiles=smiles[0], xyz=r_1_xyz) + # r_2 = ARCSpecies(label='r2', smiles=smiles[1], xyz=r_2_xyz) + # p_1 = ARCSpecies(label='p1', smiles=smiles[2], xyz=p_1_xyz) + # p_2 = ARCSpecies(label='p2', smiles=smiles[3], xyz=p_2_xyz) + # + # rxn1 = ARCReaction(r_species=[r_1, r_2], p_species=[p_1, p_2]) + # atom_map = rxn1.atom_map + # self.assertEqual(atom_map[:3],[0,2,3]) + # self.assertIn(atom_map[3:5],[[1,4],[4,1]]) + # self.assertEqual(atom_map[5],5) + # self.assertIn(atom_map[6], [6,7,8]) + # self.assertIn(atom_map[7], [6, 7, 8]) + # self.assertIn(atom_map[8], [6, 7, 8]) + # self.assertIn(atom_map[9], [9, 10, 11]) + # self.assertIn(atom_map[10], [9, 10, 11]) + # self.assertIn(atom_map[11], [9, 10, 11]) + # self.assertTrue(check_atom_map(rxn1)) def test_map_ho2_elimination_from_peroxy_radical(self): """Test the map_ho2_elimination_from_peroxy_radical() function.""" diff --git a/arc/mapping/engine.py b/arc/mapping/engine.py index 4853245c2c..def3109b0c 100644 --- a/arc/mapping/engine.py +++ b/arc/mapping/engine.py @@ -1098,7 +1098,7 @@ def make_bond_changes(rxn: 'ARCReaction', r_cuts: the cut products r_label_dict: the dictionary object the find the relevant location. """ - for action in rxn.family.forward_recipe.actions: + for action in ReactionFamily(label=rxn.family).actions: if action[0].lower() == "CHANGE_BOND".lower(): indicies = r_label_dict[action[1]],r_label_dict[action[3]] for r_cut in r_cuts: @@ -1130,7 +1130,7 @@ def make_bond_changes(rxn: 'ARCReaction', r_cut.mol.update() -def get_p_label_dict(rxn: 'ARCReaction'): +def get_label_dict(rxn: 'ARCReaction') -> Optional[Dict[str, int]]: """ A function that returns the labels of the products. @@ -1151,19 +1151,28 @@ def assign_labels_to_products(rxn: 'ARCReaction'): Add the indices to the reactants and products. Args: - rxn: ARCReaction object to be mapped - Consider changing in rmgpy. + rxn ('ARCReaction'): The reaction to be mapped. Returns: Adding labels to the atoms of the reactants and products, to be identified later. """ - p_label_dict = get_p_label_dict(rxn) - print(p_label_dict) + label_dict = get_label_dict(rxn) + print(f'\n\nlabel_dict: {label_dict}\n\n') + atom_index = 0 + for r in rxn.r_species: + for atom in r.mol.atoms: + if atom_index in label_dict.values(): + atom.label = key_by_val(label_dict, atom_index) + atom_index += 1 + + + + atom_index = 0 for product in rxn.p_species: for atom in product.mol.atoms: - if atom_index in p_label_dict.values() and (atom.label is str or atom.label is None): - atom.label = key_by_val(p_label_dict, atom_index) + if atom_index in label_dict.values() and (atom.label is str or atom.label is None): + atom.label = key_by_val(label_dict, atom_index) atom_index += 1 @@ -1206,7 +1215,7 @@ def r_cut_p_cut_isomorphic(reactant, product): def pairing_reactants_and_products_for_mapping(r_cuts: List[ARCSpecies], p_cuts: List[ARCSpecies] - )-> List[Tuple[ARCSpecies,ARCSpecies]]: + ) -> List[Tuple[ARCSpecies,ARCSpecies]]: """ A function for matching reactants and products in scissored products. @@ -1222,7 +1231,7 @@ def pairing_reactants_and_products_for_mapping(r_cuts: List[ARCSpecies], for product_cut in p_cuts: if r_cut_p_cut_isomorphic(reactant_cut, product_cut): pairs.append((reactant_cut, product_cut)) - p_cuts.remove(product_cut) # Just in case there are two of the same species in the list, matching them by order. + p_cuts.remove(product_cut) # Just in case there are two of the same species in the list, matching them by order. break return pairs @@ -1361,18 +1370,24 @@ def copy_species_list_for_mapping(species: List["ARCSpecies"]) -> List["ARCSpeci return copies -def find_all_bdes(rxn: "ARCReaction", label_dict: dict, is_reactants: bool) -> List[Tuple[int, int]]: +def find_all_bdes(rxn: "ARCReaction", + is_reactants: bool, + ) -> List[Tuple[int, int]]: """ A function for finding all the broken(/formed) bonds during a chemical reaction, based on the atom indices. + Args: rxn (ARCReaction): The reaction in question. - label_dict (dict): A dictionary of the atom indices to the atom labels. is_reactants (bool): Whether the species list represents reactants or products. + Returns: - List[Tuple[int, int]]: A list of tuples of the form (atom_index1, atom_index2) for each broken bond. Note that these represent the atom indicies to be cut, and not final BDEs. + List[Tuple[int, int]]: A list of tuples of the form (atom_index1, atom_index2) for each broken bond. + Note that these represent the atom indices to be cut, and not final BDEs. """ + label_dict = get_label_dict(rxn) bdes = list() - for action in ReactionFamily(rxn.family).actions: - if action[0].lower() == ("break_bond" if is_reactants else "form_bond"): - bdes.append((label_dict[action[1]] + 1, label_dict[action[3]] + 1)) + if rxn.family is not None: + for action in ReactionFamily(rxn.family).actions: + if action[0].lower() == ("break_bond" if is_reactants else "form_bond"): + bdes.append((label_dict[action[1]] + 1, label_dict[action[3]] + 1)) return bdes diff --git a/arc/mapping/engine_test.py b/arc/mapping/engine_test.py index 8a4099f94e..eb2a1fbb84 100644 --- a/arc/mapping/engine_test.py +++ b/arc/mapping/engine_test.py @@ -34,64 +34,26 @@ def setUpClass(cls): """ cls.maxDiff = None - smiles = ['CC(C)F', '[CH3]', 'C[CH](C)', 'CF'] - - r_1_1_xyz = {'symbols': ('C', 'C', 'C', 'F', 'H', 'H', 'H', 'H', 'H', 'H', 'H'), - 'isotopes': (12, 12, 12, 19, 1, 1, 1, 1, 1, 1, 1), - 'coords': ((1.2509680857915237, 0.00832885083067477, -0.28594855682006387), - (-0.08450322338173592, -0.5786110309038947, 0.12835305965368538), - (-1.196883483105121, 0.4516770584363101, 0.10106807955582568), - (0.03212452836861426, -1.0465351442062332, 1.402047416169314), - (1.2170230403876368, 0.39373449465586885, -1.309310880313081), - (1.5446944155971303, 0.8206316657310906, 0.38700047363833845), - (2.0327466889922805, -0.7555292157466509, -0.22527487012253536), - (-0.3397419937928473, -1.4280299782557704, -0.5129583662636836), - (-0.9791793765226446, 1.2777482351478369, 0.786037216866474), - (-1.340583396165929, 0.8569620299504027, -0.9049411765144166), - (-2.1366652861689137, -0.00037696563964776297, 0.43392760415012316))} - - r_2_1_xyz = {'symbols': ('C', 'H', 'H', 'H'), - 'isotopes': (12, 1, 1, 1), - 'coords': ((3.3746019998564553e-09, 5.828827384106545e-09, -4.859105107686622e-09), - (1.0669051052331406, -0.17519582095514982, 0.05416492980439295), - (-0.6853171627400634, -0.8375353626879753, -0.028085652887100996), - (-0.3815879458676787, 1.0127311778142964, -0.026079272058187608))} - - p_1_1_xyz = {'symbols': ('C', 'C', 'C', 'H', 'H', 'H', 'H', 'H', 'H', 'H'), - 'isotopes': (12, 12, 12, 1, 1, 1, 1, 1, 1, 1), - 'coords': ((-1.288730238258946, 0.06292843803165035, 0.10889818910854648), - (0.01096160773224897, -0.45756396262445836, -0.3934214957819532), - (1.2841030977199492, 0.11324607936811129, 0.12206176848573647), - (-1.4984446521053447, 1.0458196461796345, -0.3223873567509909), - (-1.2824724918369017, 0.14649429503996203, 1.1995362776757934), - (-2.098384694966955, -0.616646552269074, -0.17318515188247927), - (0.027360233461550892, -1.0601383387124987, -1.2952225290380646), - (2.122551165381095, -0.534098313164123, -0.15158596254231563), - (1.2634262459696732, 0.19628891975881263, 1.2125616721427255), - (1.4596297269035956, 1.1036697883919826, -0.307255411416999))} - - p_2_1_xyz = {'symbols': ('C', 'F', 'H', 'H', 'H'), - 'isotopes': (12, 19, 1, 1, 1), - 'coords': ((-0.060384822736851786, 0.004838867136375763, -0.004814368798794687), - (1.2877092002693546, -0.10318918150563985, 0.10266661058725791), - (-0.2965861926821434, 0.9189121874074381, -0.5532990701789506), - (-0.44047773762823295, -0.8660709320146035, -0.5425894744224189), - (-0.49026044722212864, 0.04550905897643097, 0.9980363028129072))} - cls.r_1 = ARCSpecies(label='r1', smiles=smiles[0], xyz=r_1_1_xyz) - cls.r_2 = ARCSpecies(label='r2', smiles=smiles[1], xyz=r_2_1_xyz) - cls.p_1 = ARCSpecies(label='p1', smiles=smiles[2], xyz=p_1_1_xyz) - cls.p_2 = ARCSpecies(label='p2', smiles=smiles[3], xyz=p_2_1_xyz) - - cls.r_1_1 = ARCSpecies(label='r1', smiles=smiles[0], xyz=r_1_1_xyz) - cls.r_2_2 = ARCSpecies(label='r2', smiles=smiles[1], xyz=r_2_1_xyz) - cls.p_1_1 = ARCSpecies(label='p1', smiles=smiles[2], xyz=p_1_1_xyz) - cls.p_2_2 = ARCSpecies(label='p2', smiles=smiles[3], xyz=p_2_1_xyz) - - # cls.rxn_1 = ARCReaction(r_species=[cls.r_1, cls.r_2], p_species=[cls.p_1, cls.p_2]) - # cls.rmg_reactions_rxn_1 = get_rmg_reactions_from_arc_reaction(arc_reaction=cls.rxn_1, backend="ARC") - # cls.r_label_dict_rxn_1, cls.p_label_dict_rxn_1 = ( - # get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=cls.rxn_1, - # rmg_reaction=cls.rmg_reactions_rxn_1[0])) + smiles = ['[C]#CF', 'C#N', 'C#CF', '[C]#N'] + r_1_xyz = {'symbols': ('C', 'C', 'F'), 'isotopes': (12, 12, 19), + 'coords': ((-1.4743294980458657, -0.0015578181856664775, 0.0), + (0.06566963887880228, 6.93883492559459e-05, 0.0), + (1.4086598591670647, 0.0014884298364102333, 0.0))} + r_2_xyz = {'symbols': ('C', 'N', 'H'), 'isotopes': (12, 14, 1), + 'coords': ((-0.03166549265366707, 0.00027268185540903927, 0.0), + (1.1282915045547195, -0.009715961185409256, 0.0), + (-1.096626011901052, 0.009443279330001001, 0.0))} + p_1_xyz = {'symbols': ('C', 'C', 'F', 'H'), 'isotopes': (12, 12, 19, 1), + 'coords': ((-0.6556450958206512, 0.008288859693035082, 0.041287880045658594), + (0.5194805538847015, -0.006567413538153382, 0.28295662316778103), + (1.8344056155700479, -0.023191047733068296, 0.5533756165580456), + (-1.6982410736345843, 0.0214696015778594, -0.17312568563836497))} + p_2_xyz = {'symbols': ('C', 'N'), 'isotopes': (12, 14), 'coords': ((0.0, 0.0, 0.735), (0.0, 0.0, -0.735))} + + cls.r_1 = ARCSpecies(label='r1', smiles=smiles[0], xyz=r_1_xyz) + cls.r_2 = ARCSpecies(label='r2', smiles=smiles[1], xyz=r_2_xyz) + cls.p_1 = ARCSpecies(label='p1', smiles=smiles[2], xyz=p_1_xyz) + cls.p_2 = ARCSpecies(label='p2', smiles=smiles[3], xyz=p_2_xyz) cls.spc1 = ARCSpecies(label="Test_is_isomorphic_1",smiles="C(CO)CC") cls.spc2 = ARCSpecies(label="Test_is_isomorphic_2",smiles="OCCCC") @@ -129,15 +91,12 @@ def setUpClass(cls): (0.7709290243761263, -0.38422053705817527, 1.4054470816682596), (1.337910696892115, -1.3171321272490044, 0.001256204546378134), (2.595292874972368, 0.5254618254772234, 0.4066018700054956))}) - - cls.ch4_xyz = {'symbols': ('C', 'H', 'H', 'H', 'H'), 'isotopes': (12, 1, 1, 1, 1), 'coords': ((-5.45906343962835e-10, 4.233517924761169e-10, 2.9505240956083194e-10), (-0.6505520089868748, -0.7742801979689132, -0.4125187934483119), (-0.34927557824779626, 0.9815958255612931, -0.3276823191685369), (-0.022337921721882443, -0.04887374527620588, 1.0908766524267022), (1.0221655095024578, -0.15844188273952128, -0.350675540104908))} - cls.ch3cl_xyz = { "symbols": ("Cl", "C", "H", "H", "H"), "isotopes": (35, 12, 1, 1, 1), @@ -146,10 +105,7 @@ def setUpClass(cls): (-0.13902319255720844, -0.005942522549225256, 0.0025040105508790877), (-0.48751519000626353, -0.5782116958350602, -0.8600292104425608), (-0.45894137315516464, -0.4942789247294056, 0.9255869621295756), - (-0.540565317573775, 1.0089280316586664, -0.038774361328407066), - ), - } - + (-0.540565317573775, 1.0089280316586664, -0.038774361328407066))} cls.ch3_xyz_2 = { "symbols": ("C", "H", "H", "H"), "isotopes": (12, 1, 1, 1), @@ -157,18 +113,12 @@ def setUpClass(cls): (3.3746019998564553e-09, 5.828827384106545e-09, -4.859105107686622e-09), (1.0669051052331406, -0.17519582095514982, 0.05416492980439295), (-0.6853171627400634, -0.8375353626879753, -0.028085652887100996), - (-0.3815879458676787, 1.0127311778142964, -0.026079272058187608), - ), - } - + (-0.3815879458676787, 1.0127311778142964, -0.026079272058187608))} cls.h_rad_xyz = {"symbols": ("H",), "isotopes": (1,), "coords": ((0, 0, 0),)} - cls.hcl_xyz = { "symbols": ("H", "Cl"), "isotopes": (1, 35), - "coords": ((0.6878248644303301, 0.0, 0.0), (-0.6878248644303301, 0.0, 0.0)), - } - + "coords": ((0.6878248644303301, 0.0, 0.0), (-0.6878248644303301, 0.0, 0.0))} cls.ch4_xyz = { "symbols": ("C", "H", "H", "H", "H"), "isotopes": (12, 1, 1, 1, 1), @@ -177,10 +127,7 @@ def setUpClass(cls): (-0.6505520089868748, -0.7742801979689132, -0.4125187934483119), (-0.34927557824779626, 0.9815958255612931, -0.3276823191685369), (-0.022337921721882443, -0.04887374527620588, 1.0908766524267022), - (1.0221655095024578, -0.15844188273952128, -0.350675540104908), - ), - } - + (1.0221655095024578, -0.15844188273952128, -0.350675540104908))} cls.ch4_xyz_diff_order = """H -0.65055201 -0.77428020 -0.41251879 H -0.34927558 0.98159583 -0.32768232 C -0.00000000 0.00000000 0.00000000 @@ -562,14 +509,11 @@ def test_pair_reaction_products(self): product_pairs = pair_reaction_products(rxn_1, [prod_mol, prod_mol]) self.assertEqual(product_pairs, {0: 0, 1: 1}) - - def test_assign_labels_to_products(self): """Test assigning labels to products based on the atom map of the reaction""" rxn_1_test = ARCReaction(r_species=[self.r_1, self.r_2], p_species=[self.p_1, self.p_2]) assign_labels_to_products(rxn_1_test) print([atom.label for atom in rxn_1_test.p_species[0].mol.atoms]) - raise index = 0 for product in rxn_1_test.p_species: print(product.label, index) @@ -588,37 +532,31 @@ def test_inc_vals(self): def test_label_species_atoms(self): rxn_1_test = ARCReaction(r_species=[self.r_1, self.r_2], p_species=[self.p_1, self.p_2]) - assign_labels_to_products(rxn_1_test, self.p_label_dict_rxn_1) - + assign_labels_to_products(rxn_1_test) reactants, products = copy_species_list_for_mapping(rxn_1_test.r_species), copy_species_list_for_mapping(rxn_1_test.p_species) - label_species_atoms(reactants) label_species_atoms(products) - index = 0 for reactant in reactants: for atom in reactant.mol.atoms: self.assertEqual(atom.label,str(index)) - index +=1 - + index += 1 index = 0 for product in products: for atom in product.mol.atoms: - self.assertEqual(atom.label,str(index)) - index +=1 + self.assertEqual(atom.label, str(index)) + index += 1 def test_cut_species_based_on_atom_indices(self): """test the cut_species_for_mapping function""" - rxn_1_test = ARCReaction(r_species=[self.r_1, self.r_2], p_species=[self.p_1, self.p_2]) + rxn_1_test = ARCReaction(r_species=[self.r_1, self.r_2], p_species=[self.p_1, self.p_2], + rmg_family_set=['F_Abstraction']) reactants, products = copy_species_list_for_mapping(rxn_1_test.r_species), copy_species_list_for_mapping(rxn_1_test.p_species) label_species_atoms(reactants), label_species_atoms(products) - - r_bdes, p_bdes = find_all_bdes(rxn_1_test, self.r_label_dict_rxn_1, True), find_all_bdes(rxn_1_test, self.p_label_dict_rxn_1, False) - + r_bdes, p_bdes = find_all_bdes(rxn_1_test, True), find_all_bdes(rxn_1_test, False) r_cuts = cut_species_based_on_atom_indices(reactants, r_bdes) p_cuts = cut_species_based_on_atom_indices(products, p_bdes) - self.assertIn("C[CH]C", [r_cut.mol.copy(deep=True).smiles for r_cut in r_cuts]) self.assertIn("[F]", [r_cut.mol.copy(deep=True).smiles for r_cut in r_cuts]) self.assertIn("[CH3]", [r_cut.mol.copy(deep=True).smiles for r_cut in r_cuts]) @@ -626,9 +564,9 @@ def test_cut_species_based_on_atom_indices(self): self.assertIn("[F]", [p_cut.mol.copy(deep=True).smiles for p_cut in p_cuts]) self.assertIn("[CH3]", [p_cut.mol.copy(deep=True).smiles for p_cut in p_cuts]) - spc = ARCSpecies(label="test", smiles="CNC", bdes = [(1, 2), (2, 3)]) + spc = ARCSpecies(label="test", smiles="CNC", bdes=[(1, 2), (2, 3)]) for i, a in enumerate(spc.mol.atoms): - a.label=str(i) + a.label = str(i) cuts = cut_species_based_on_atom_indices([spc], [(1, 2), (2, 3)]) self.assertEqual(len(cuts), 3) for cut in cuts: @@ -638,12 +576,12 @@ def test_cut_species_based_on_atom_indices(self): h2 = ARCSpecies(label="H2", smiles="[H][H]") label_species_atoms([h2]) - cuts = cut_species_based_on_atom_indices([h2], [(1, 2)]) + cuts = cut_species_based_on_atom_indices([h2], [(1, 2)]) self.assertEqual(len(cuts), 2) for cut in cuts: self.assertEqual(cut.get_xyz()["symbols"], ('H',)) - - spcs = [ARCSpecies(label="r", smiles = 'O=C(O)CCF')] + + spcs = [ARCSpecies(label="r", smiles='O=C(O)CCF')] label_species_atoms(spcs) cuts = cut_species_based_on_atom_indices(spcs, [(6, 5), (4, 2), (3, 7)]) self.assertEqual(len(cuts), 4) @@ -652,9 +590,7 @@ def test_r_cut_p_cut_isomorphic(self): rxn_1_test = ARCReaction(r_species=[self.r_1, self.r_2], p_species=[self.p_1, self.p_2]) reactants, products = copy_species_list_for_mapping(rxn_1_test.r_species), copy_species_list_for_mapping(rxn_1_test.p_species) label_species_atoms(reactants), label_species_atoms(products) - - r_bdes, p_bdes = find_all_bdes(rxn_1_test, self.r_label_dict_rxn_1, True), find_all_bdes(rxn_1_test, self.p_label_dict_rxn_1, False) - + r_bdes, p_bdes = find_all_bdes(rxn_1_test, True), find_all_bdes(rxn_1_test, False) r_cuts = cut_species_based_on_atom_indices(reactants, r_bdes) p_cuts = cut_species_based_on_atom_indices(products, p_bdes) @@ -673,9 +609,7 @@ def test_pairing_reactants_and_products_for_mapping(self): rxn_1_test = ARCReaction(r_species=[self.r_1, self.r_2], p_species=[self.p_1, self.p_2]) reactants, products = copy_species_list_for_mapping(rxn_1_test.r_species), copy_species_list_for_mapping(rxn_1_test.p_species) label_species_atoms(reactants), label_species_atoms(products) - - r_bdes, p_bdes = find_all_bdes(rxn_1_test, self.r_label_dict_rxn_1, True), find_all_bdes(rxn_1_test, self.p_label_dict_rxn_1, False) - + r_bdes, p_bdes = find_all_bdes(rxn_1_test, True), find_all_bdes(rxn_1_test, False) r_cuts = cut_species_based_on_atom_indices(reactants, r_bdes) p_cuts = cut_species_based_on_atom_indices(products, p_bdes) @@ -686,16 +620,11 @@ def test_pairing_reactants_and_products_for_mapping(self): self.assertIn(str(pair[0].mol.copy(deep=True).to_smiles()), smiles) smiles.remove(pair[0].mol.copy(deep=True).to_smiles()) - def test_map_pairs(self): + def test_map_pairs(self): # H abs reverse 1st rxn_1_test = ARCReaction(r_species=[self.r_1, self.r_2], p_species=[self.p_1, self.p_2]) - rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=rxn_1_test, backend="ARC") - r_label_dict, p_label_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=rxn_1_test, - rmg_reaction=rmg_reactions[0]) reactants, products = copy_species_list_for_mapping(rxn_1_test.r_species), copy_species_list_for_mapping(rxn_1_test.p_species) label_species_atoms(reactants), label_species_atoms(products) - - r_bdes, p_bdes = find_all_bdes(rxn_1_test, r_label_dict, True), find_all_bdes(rxn_1_test, p_label_dict, False) - + r_bdes, p_bdes = find_all_bdes(rxn_1_test, True), find_all_bdes(rxn_1_test, False) r_cuts = cut_species_based_on_atom_indices(reactants, r_bdes) p_cuts = cut_species_based_on_atom_indices(products, p_bdes) pairs_of_reactant_and_products = pairing_reactants_and_products_for_mapping(r_cuts, p_cuts) @@ -717,15 +646,10 @@ def test_map_pairs(self): def test_glue_maps(self): rxn_1_test = ARCReaction(r_species=[self.r_1, self.r_2], p_species=[self.p_1, self.p_2]) - rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=rxn_1_test, backend="ARC") - r_label_dict, p_label_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=rxn_1_test, - rmg_reaction=rmg_reactions[0]) - assign_labels_to_products(rxn_1_test, p_label_dict) + assign_labels_to_products(rxn_1_test) reactants, products = copy_species_list_for_mapping(rxn_1_test.r_species), copy_species_list_for_mapping(rxn_1_test.p_species) label_species_atoms(reactants), label_species_atoms(products) - - r_bdes, p_bdes = find_all_bdes(rxn_1_test, r_label_dict, True), find_all_bdes(rxn_1_test, p_label_dict, False) - + r_bdes, p_bdes = find_all_bdes(rxn_1_test, True), find_all_bdes(rxn_1_test, False) r_cuts = cut_species_based_on_atom_indices(reactants, r_bdes) p_cuts = cut_species_based_on_atom_indices(products, p_bdes) pairs_of_reactant_and_products = pairing_reactants_and_products_for_mapping(r_cuts, p_cuts) @@ -737,206 +661,205 @@ def test_glue_maps(self): for index, value in enumerate(atom_map): self.assertEqual(atoms_r[index].symbol, atoms_p[value].symbol) - def test_get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(self): - """Test the get_atom_indices_of_labeled_atoms_in_an_rmg_reaction() function.""" - rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=self.arc_reaction_1) - r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=self.arc_reaction_1, - rmg_reaction=rmg_reactions[0]) - self.assertEqual(r_dict['*1'], 0) - self.assertIn(r_dict['*2'], [1, 2, 3, 4]) - self.assertEqual(r_dict['*3'], 5) - self.assertEqual(p_dict['*1'], 0) - self.assertIn(p_dict['*2'], [5, 6]) - self.assertEqual(p_dict['*3'], 4) - self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(self.arc_reaction_1, rmg_reactions)) - - rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=self.arc_reaction_2) - r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=self.arc_reaction_2, - rmg_reaction=rmg_reactions[0]) - self.assertIn(r_dict['*1'], [0, 2]) - self.assertIn(r_dict['*2'], [3, 4, 5, 8, 9, 10]) - self.assertEqual(r_dict['*3'], 11) - self.assertEqual(p_dict['*1'], 0) - self.assertIn(p_dict['*2'], [11, 12, 13]) - self.assertEqual(p_dict['*3'], 10) - self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(self.arc_reaction_2, rmg_reactions)) - - rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=self.arc_reaction_4) - r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=self.arc_reaction_4, - rmg_reaction=rmg_reactions[0]) - self.assertEqual(r_dict['*1'], 0) - self.assertEqual(r_dict['*2'], 2) - self.assertIn(r_dict['*3'], [7, 8]) - self.assertEqual(p_dict['*1'], 0) - self.assertEqual(p_dict['*2'], 2) - self.assertIn(p_dict['*3'], [3, 4, 5]) - self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(self.arc_reaction_4, rmg_reactions)) - - for atom, symbol in zip(self.rxn_2a.r_species[0].mol.atoms, ['C', 'C', 'C', 'H', 'H', 'H', 'H', 'H', 'H', 'H']): - self.assertEqual(atom.symbol, symbol) - self.assertEqual(self.rxn_2a.r_species[0].mol.atoms[0].radical_electrons, 0) - self.assertEqual(self.rxn_2a.r_species[0].mol.atoms[1].radical_electrons, 1) - self.assertEqual(self.rxn_2a.r_species[0].mol.atoms[2].radical_electrons, 0) - self.assertEqual(self.rxn_2a.p_species[0].mol.atoms[0].radical_electrons, 0) - self.assertEqual(self.rxn_2a.p_species[0].mol.atoms[1].radical_electrons, 0) - self.assertEqual(self.rxn_2a.p_species[0].mol.atoms[2].radical_electrons, 1) - rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=self.rxn_2a) - r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=self.rxn_2a, - rmg_reaction=rmg_reactions[0]) - self.assertEqual(r_dict['*1'], 1) - self.assertIn(r_dict['*2'], [0, 2]) - self.assertIn(r_dict['*3'], [4, 5, 6, 7, 8, 9]) - self.assertEqual(p_dict['*1'], 1) - self.assertEqual(p_dict['*2'], 2) - self.assertIn(p_dict['*3'], [3, 6]) - self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(self.rxn_2a, rmg_reactions)) - - for atom, symbol in zip(self.rxn_2b.r_species[0].mol.atoms, ['C', 'C', 'H', 'H', 'H', 'H', 'C', 'H', 'H', 'H']): - self.assertEqual(atom.symbol, symbol) - rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=self.rxn_2b) - r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=self.rxn_2b, - rmg_reaction=rmg_reactions[0]) - self.assertEqual(r_dict['*1'], 1) - self.assertIn(r_dict['*2'], [0, 6]) - self.assertIn(r_dict['*3'], [3, 4, 5, 7, 8, 9]) - self.assertEqual(p_dict['*1'], 1) - self.assertEqual(p_dict['*2'], 2) - self.assertIn(p_dict['*3'], [3, 6]) - self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(self.rxn_2b, rmg_reactions)) - - # C3H6O + C4H9O <=> C3H5O + C4H10O - r_1 = ARCSpecies(label='C3H6O', smiles='CCC=O') - r_2 = ARCSpecies(label='C4H9O', smiles='[CH2]C(C)CO') - p_1 = ARCSpecies(label='C3H5O', smiles='C[CH]C=O') - p_2 = ARCSpecies(label='C4H10O', smiles='CC(C)CO') - rxn_1 = ARCReaction(reactants=['C3H6O', 'C4H9O'], products=['C3H5O', 'C4H10O'], - r_species=[r_1, r_2], p_species=[p_1, p_2]) - rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=rxn_1) - for rmg_reaction in rmg_reactions: - r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=rxn_1, - rmg_reaction=rmg_reaction) - for d in [r_dict, p_dict]: - self.assertEqual(len(list(d.keys())), 3) - keys = list(d.keys()) - for label in ['*1', '*2', '*3']: - self.assertIn(label, keys) - self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(rxn_1, rmg_reactions)) - - p_1 = ARCSpecies(label='C3H5O', smiles='CC=C[O]') # Use a wrong resonance structure and repeat the above. - rxn_2 = ARCReaction(reactants=['C3H6O', 'C4H9O'], products=['C3H5O', 'C4H10O'], - r_species=[r_1, r_2], p_species=[p_1, p_2]) - rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=rxn_2) - for rmg_reaction in rmg_reactions: - r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=rxn_2, - rmg_reaction=rmg_reaction) - for d in [r_dict, p_dict]: - self.assertEqual(len(list(d.keys())), 3) - keys = list(d.keys()) - for label in ['*1', '*2', '*3']: - self.assertIn(label, keys) - self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(rxn_2, rmg_reactions)) - - # C3H6O + C4H9O <=> C3H5O + C4H10O - r_1 = ARCSpecies(label='C3H6O', smiles='CCC=O', xyz=self.c3h6o_xyz) - r_2 = ARCSpecies(label='C4H9O', smiles='[CH2]C(C)CO', xyz=self.c4h9o_xyz) - p_1 = ARCSpecies(label='C3H5O', smiles='C[CH]C=O', xyz=self.c3h5o_xyz) - p_2 = ARCSpecies(label='C4H10O', smiles='CC(C)CO', xyz=self.c4h10o_xyz) - rxn_3 = ARCReaction(reactants=['C3H6O', 'C4H9O'], products=['C3H5O', 'C4H10O'], - r_species=[r_1, r_2], p_species=[p_1, p_2]) - rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=rxn_3) - r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=rxn_3, - rmg_reaction=rmg_reactions[0]) - self.assertEqual(r_dict, {'*3': 10, '*1': 1, '*2': 7}) - self.assertEqual(p_dict, {'*1': 1, '*3': 9, '*2': 16}) - self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(rxn_3, rmg_reactions)) - - def test_map_arc_rmg_species(self): - """Test the map_arc_rmg_species() function.""" - r_map, p_map = map_arc_rmg_species( - arc_reaction=ARCReaction(r_species=[ARCSpecies(label='CCjC', smiles='C[CH]C')], - p_species=[ARCSpecies(label='CjCC', smiles='[CH2]CC')]), - rmg_reaction=Reaction(reactants=[Species(smiles='C[CH]C')], - products=[Species(smiles='[CH2]CC')]), - concatenate=False) - self.assertEqual(r_map, {0: 0}) - self.assertEqual(p_map, {0: 0}) - - r_map, p_map = map_arc_rmg_species( - arc_reaction=ARCReaction(r_species=[ARCSpecies(label='CCjC', smiles='C[CH]C')], - p_species=[ARCSpecies(label='CjCC', smiles='[CH2]CC')]), - rmg_reaction=Reaction(reactants=[Species(smiles='C[CH]C')], - products=[Species(smiles='[CH2]CC')])) - self.assertEqual(r_map, {0: [0]}) - self.assertEqual(p_map, {0: [0]}) - - r_map, p_map = map_arc_rmg_species(arc_reaction=self.arc_reaction_1, rmg_reaction=self.rmg_reaction_1) - self.assertEqual(r_map, {0: [0], 1: [1]}) - self.assertEqual(p_map, {0: [0], 1: [1]}) - - r_map, p_map = map_arc_rmg_species(arc_reaction=self.arc_reaction_1, rmg_reaction=self.rmg_reaction_2) - self.assertEqual(r_map, {0: [1], 1: [0]}) - self.assertEqual(p_map, {0: [0], 1: [1]}) - - r_map, p_map = map_arc_rmg_species(arc_reaction=self.arc_reaction_3, rmg_reaction=self.rmg_reaction_3) - self.assertEqual(r_map, {0: [0, 1], 1: [0, 1]}) - self.assertEqual(p_map, {0: [0]}) - - rmg_reaction_1 = Reaction(reactants=[Species(smiles='N'), Species(smiles='[H]')], - products=[Species(smiles='[NH2]'), Species(smiles='[H][H]')]) - rmg_reaction_2 = Reaction(reactants=[Species(smiles='[H]'), Species(smiles='N')], - products=[Species(smiles='[H][H]'), Species(smiles='[NH2]')]) - rmg_reaction_3 = Reaction(reactants=[Species(smiles='N'), Species(smiles='[H]')], - products=[Species(smiles='[H][H]'), Species(smiles='[NH2]')]) - arc_reaction = ARCReaction(r_species=[ARCSpecies(label='NH3', smiles='N'), ARCSpecies(label='H', smiles='[H]')], - p_species=[ARCSpecies(label='NH2', smiles='[NH2]'), ARCSpecies(label='H2', smiles='[H][H]')]) - - r_map, p_map = map_arc_rmg_species(rmg_reaction=rmg_reaction_1, arc_reaction=arc_reaction) - self.assertEqual(r_map, {0: [0], 1: [1]}) - self.assertEqual(p_map, {0: [0], 1: [1]}) - - r_map, p_map = map_arc_rmg_species(rmg_reaction=rmg_reaction_2, arc_reaction=arc_reaction) - self.assertEqual(r_map, {0: [1], 1: [0]}) - self.assertEqual(p_map, {0: [1], 1: [0]}) - - r_map, p_map = map_arc_rmg_species(rmg_reaction=rmg_reaction_3, arc_reaction=arc_reaction) - self.assertEqual(r_map, {0: [0], 1: [1]}) - self.assertEqual(p_map, {0: [1], 1: [0]}) - rmg_reaction = Reaction(reactants=[Species(smiles='[CH3]'), Species(smiles='[CH3]')], - products=[Species(smiles='CC')]) - arc_reaction = ARCReaction(r_species=[ARCSpecies(label='CH3', smiles='[CH3]'), - ARCSpecies(label='CH3', smiles='[CH3]')], - p_species=[ARCSpecies(label='C2H6', smiles='CC')]) - - r_map, p_map = map_arc_rmg_species(rmg_reaction=rmg_reaction, arc_reaction=arc_reaction) - self.assertEqual(r_map, {0: [0, 1], 1: [0, 1]}) - self.assertEqual(p_map, {0: [0]}) - - rmg_reaction = Reaction(reactants=[Species(smiles='c1ccccc1CCC=C'), Species(smiles='CCO[O]')], - products=[Species(smiles='c1ccccc1[CH]CC=C'), Species(smiles='CCOO')]) - arc_reaction = ARCReaction(r_species=[ARCSpecies(label='butenylbenzene', smiles='c1ccccc1CCC=C'), - ARCSpecies(label='peroxyl', smiles='CCO[O]')], - p_species=[ARCSpecies(label='peroxide', smiles='CCOO'), - ARCSpecies(label='butenylbenzene_rad', smiles='c1ccccc1[CH]CC=C')]) - r_map, p_map = map_arc_rmg_species(rmg_reaction=rmg_reaction, - arc_reaction=arc_reaction, - concatenate=False, - ) - self.assertEqual(r_map, {0: 0, 1: 1}) - self.assertEqual(p_map, {0: 1, 1: 0}) - - rmg_reaction = Reaction(reactants=[Species(smiles='[N]N'), Species(smiles='NN')], - products=[Species(smiles='[NH]N'), Species(smiles='[NH]N')]) - arc_reaction = ARCReaction(r_species=[ARCSpecies(label='H2NN(T)', smiles='[N]N'), - ARCSpecies(label='N2H4', smiles='NN')], - p_species=[ARCSpecies(label='N2H3', smiles='[NH]N'), - ARCSpecies(label='N2H3', smiles='[NH]N')]) - r_map, p_map = map_arc_rmg_species(rmg_reaction=rmg_reaction, - arc_reaction=arc_reaction, - concatenate=False, - ) - self.assertEqual(r_map, {0: 0, 1: 1}) - self.assertEqual(p_map, {0: 0, 1: 1}) - + # def test_get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(self): + # """Test the get_atom_indices_of_labeled_atoms_in_an_rmg_reaction() function.""" + # rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=self.arc_reaction_1) + # r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=self.arc_reaction_1, + # rmg_reaction=rmg_reactions[0]) + # self.assertEqual(r_dict['*1'], 0) + # self.assertIn(r_dict['*2'], [1, 2, 3, 4]) + # self.assertEqual(r_dict['*3'], 5) + # self.assertEqual(p_dict['*1'], 0) + # self.assertIn(p_dict['*2'], [5, 6]) + # self.assertEqual(p_dict['*3'], 4) + # self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(self.arc_reaction_1, rmg_reactions)) + # + # rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=self.arc_reaction_2) + # r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=self.arc_reaction_2, + # rmg_reaction=rmg_reactions[0]) + # self.assertIn(r_dict['*1'], [0, 2]) + # self.assertIn(r_dict['*2'], [3, 4, 5, 8, 9, 10]) + # self.assertEqual(r_dict['*3'], 11) + # self.assertEqual(p_dict['*1'], 0) + # self.assertIn(p_dict['*2'], [11, 12, 13]) + # self.assertEqual(p_dict['*3'], 10) + # self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(self.arc_reaction_2, rmg_reactions)) + # + # rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=self.arc_reaction_4) + # r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=self.arc_reaction_4, + # rmg_reaction=rmg_reactions[0]) + # self.assertEqual(r_dict['*1'], 0) + # self.assertEqual(r_dict['*2'], 2) + # self.assertIn(r_dict['*3'], [7, 8]) + # self.assertEqual(p_dict['*1'], 0) + # self.assertEqual(p_dict['*2'], 2) + # self.assertIn(p_dict['*3'], [3, 4, 5]) + # self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(self.arc_reaction_4, rmg_reactions)) + # + # for atom, symbol in zip(self.rxn_2a.r_species[0].mol.atoms, ['C', 'C', 'C', 'H', 'H', 'H', 'H', 'H', 'H', 'H']): + # self.assertEqual(atom.symbol, symbol) + # self.assertEqual(self.rxn_2a.r_species[0].mol.atoms[0].radical_electrons, 0) + # self.assertEqual(self.rxn_2a.r_species[0].mol.atoms[1].radical_electrons, 1) + # self.assertEqual(self.rxn_2a.r_species[0].mol.atoms[2].radical_electrons, 0) + # self.assertEqual(self.rxn_2a.p_species[0].mol.atoms[0].radical_electrons, 0) + # self.assertEqual(self.rxn_2a.p_species[0].mol.atoms[1].radical_electrons, 0) + # self.assertEqual(self.rxn_2a.p_species[0].mol.atoms[2].radical_electrons, 1) + # rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=self.rxn_2a) + # r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=self.rxn_2a, + # rmg_reaction=rmg_reactions[0]) + # self.assertEqual(r_dict['*1'], 1) + # self.assertIn(r_dict['*2'], [0, 2]) + # self.assertIn(r_dict['*3'], [4, 5, 6, 7, 8, 9]) + # self.assertEqual(p_dict['*1'], 1) + # self.assertEqual(p_dict['*2'], 2) + # self.assertIn(p_dict['*3'], [3, 6]) + # self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(self.rxn_2a, rmg_reactions)) + # + # for atom, symbol in zip(self.rxn_2b.r_species[0].mol.atoms, ['C', 'C', 'H', 'H', 'H', 'H', 'C', 'H', 'H', 'H']): + # self.assertEqual(atom.symbol, symbol) + # rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=self.rxn_2b) + # r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=self.rxn_2b, + # rmg_reaction=rmg_reactions[0]) + # self.assertEqual(r_dict['*1'], 1) + # self.assertIn(r_dict['*2'], [0, 6]) + # self.assertIn(r_dict['*3'], [3, 4, 5, 7, 8, 9]) + # self.assertEqual(p_dict['*1'], 1) + # self.assertEqual(p_dict['*2'], 2) + # self.assertIn(p_dict['*3'], [3, 6]) + # self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(self.rxn_2b, rmg_reactions)) + # + # # C3H6O + C4H9O <=> C3H5O + C4H10O + # r_1 = ARCSpecies(label='C3H6O', smiles='CCC=O') + # r_2 = ARCSpecies(label='C4H9O', smiles='[CH2]C(C)CO') + # p_1 = ARCSpecies(label='C3H5O', smiles='C[CH]C=O') + # p_2 = ARCSpecies(label='C4H10O', smiles='CC(C)CO') + # rxn_1 = ARCReaction(reactants=['C3H6O', 'C4H9O'], products=['C3H5O', 'C4H10O'], + # r_species=[r_1, r_2], p_species=[p_1, p_2]) + # rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=rxn_1) + # for rmg_reaction in rmg_reactions: + # r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=rxn_1, + # rmg_reaction=rmg_reaction) + # for d in [r_dict, p_dict]: + # self.assertEqual(len(list(d.keys())), 3) + # keys = list(d.keys()) + # for label in ['*1', '*2', '*3']: + # self.assertIn(label, keys) + # self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(rxn_1, rmg_reactions)) + # + # p_1 = ARCSpecies(label='C3H5O', smiles='CC=C[O]') # Use a wrong resonance structure and repeat the above. + # rxn_2 = ARCReaction(reactants=['C3H6O', 'C4H9O'], products=['C3H5O', 'C4H10O'], + # r_species=[r_1, r_2], p_species=[p_1, p_2]) + # rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=rxn_2) + # for rmg_reaction in rmg_reactions: + # r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=rxn_2, + # rmg_reaction=rmg_reaction) + # for d in [r_dict, p_dict]: + # self.assertEqual(len(list(d.keys())), 3) + # keys = list(d.keys()) + # for label in ['*1', '*2', '*3']: + # self.assertIn(label, keys) + # self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(rxn_2, rmg_reactions)) + # + # # C3H6O + C4H9O <=> C3H5O + C4H10O + # r_1 = ARCSpecies(label='C3H6O', smiles='CCC=O', xyz=self.c3h6o_xyz) + # r_2 = ARCSpecies(label='C4H9O', smiles='[CH2]C(C)CO', xyz=self.c4h9o_xyz) + # p_1 = ARCSpecies(label='C3H5O', smiles='C[CH]C=O', xyz=self.c3h5o_xyz) + # p_2 = ARCSpecies(label='C4H10O', smiles='CC(C)CO', xyz=self.c4h10o_xyz) + # rxn_3 = ARCReaction(reactants=['C3H6O', 'C4H9O'], products=['C3H5O', 'C4H10O'], + # r_species=[r_1, r_2], p_species=[p_1, p_2]) + # rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=rxn_3) + # r_dict, p_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=rxn_3, + # rmg_reaction=rmg_reactions[0]) + # self.assertEqual(r_dict, {'*3': 10, '*1': 1, '*2': 7}) + # self.assertEqual(p_dict, {'*1': 1, '*3': 9, '*2': 16}) + # self.assertTrue(_check_r_n_p_symbols_between_rmg_and_arc_rxns(rxn_3, rmg_reactions)) + # + # def test_map_arc_rmg_species(self): + # """Test the map_arc_rmg_species() function.""" + # r_map, p_map = map_arc_rmg_species( + # arc_reaction=ARCReaction(r_species=[ARCSpecies(label='CCjC', smiles='C[CH]C')], + # p_species=[ARCSpecies(label='CjCC', smiles='[CH2]CC')]), + # rmg_reaction=Reaction(reactants=[Species(smiles='C[CH]C')], + # products=[Species(smiles='[CH2]CC')]), + # concatenate=False) + # self.assertEqual(r_map, {0: 0}) + # self.assertEqual(p_map, {0: 0}) + # + # r_map, p_map = map_arc_rmg_species( + # arc_reaction=ARCReaction(r_species=[ARCSpecies(label='CCjC', smiles='C[CH]C')], + # p_species=[ARCSpecies(label='CjCC', smiles='[CH2]CC')]), + # rmg_reaction=Reaction(reactants=[Species(smiles='C[CH]C')], + # products=[Species(smiles='[CH2]CC')])) + # self.assertEqual(r_map, {0: [0]}) + # self.assertEqual(p_map, {0: [0]}) + # + # r_map, p_map = map_arc_rmg_species(arc_reaction=self.arc_reaction_1, rmg_reaction=self.rmg_reaction_1) + # self.assertEqual(r_map, {0: [0], 1: [1]}) + # self.assertEqual(p_map, {0: [0], 1: [1]}) + # + # r_map, p_map = map_arc_rmg_species(arc_reaction=self.arc_reaction_1, rmg_reaction=self.rmg_reaction_2) + # self.assertEqual(r_map, {0: [1], 1: [0]}) + # self.assertEqual(p_map, {0: [0], 1: [1]}) + # + # r_map, p_map = map_arc_rmg_species(arc_reaction=self.arc_reaction_3, rmg_reaction=self.rmg_reaction_3) + # self.assertEqual(r_map, {0: [0, 1], 1: [0, 1]}) + # self.assertEqual(p_map, {0: [0]}) + # + # rmg_reaction_1 = Reaction(reactants=[Species(smiles='N'), Species(smiles='[H]')], + # products=[Species(smiles='[NH2]'), Species(smiles='[H][H]')]) + # rmg_reaction_2 = Reaction(reactants=[Species(smiles='[H]'), Species(smiles='N')], + # products=[Species(smiles='[H][H]'), Species(smiles='[NH2]')]) + # rmg_reaction_3 = Reaction(reactants=[Species(smiles='N'), Species(smiles='[H]')], + # products=[Species(smiles='[H][H]'), Species(smiles='[NH2]')]) + # arc_reaction = ARCReaction(r_species=[ARCSpecies(label='NH3', smiles='N'), ARCSpecies(label='H', smiles='[H]')], + # p_species=[ARCSpecies(label='NH2', smiles='[NH2]'), ARCSpecies(label='H2', smiles='[H][H]')]) + # + # r_map, p_map = map_arc_rmg_species(rmg_reaction=rmg_reaction_1, arc_reaction=arc_reaction) + # self.assertEqual(r_map, {0: [0], 1: [1]}) + # self.assertEqual(p_map, {0: [0], 1: [1]}) + # + # r_map, p_map = map_arc_rmg_species(rmg_reaction=rmg_reaction_2, arc_reaction=arc_reaction) + # self.assertEqual(r_map, {0: [1], 1: [0]}) + # self.assertEqual(p_map, {0: [1], 1: [0]}) + # + # r_map, p_map = map_arc_rmg_species(rmg_reaction=rmg_reaction_3, arc_reaction=arc_reaction) + # self.assertEqual(r_map, {0: [0], 1: [1]}) + # self.assertEqual(p_map, {0: [1], 1: [0]}) + # rmg_reaction = Reaction(reactants=[Species(smiles='[CH3]'), Species(smiles='[CH3]')], + # products=[Species(smiles='CC')]) + # arc_reaction = ARCReaction(r_species=[ARCSpecies(label='CH3', smiles='[CH3]'), + # ARCSpecies(label='CH3', smiles='[CH3]')], + # p_species=[ARCSpecies(label='C2H6', smiles='CC')]) + # + # r_map, p_map = map_arc_rmg_species(rmg_reaction=rmg_reaction, arc_reaction=arc_reaction) + # self.assertEqual(r_map, {0: [0, 1], 1: [0, 1]}) + # self.assertEqual(p_map, {0: [0]}) + # + # rmg_reaction = Reaction(reactants=[Species(smiles='c1ccccc1CCC=C'), Species(smiles='CCO[O]')], + # products=[Species(smiles='c1ccccc1[CH]CC=C'), Species(smiles='CCOO')]) + # arc_reaction = ARCReaction(r_species=[ARCSpecies(label='butenylbenzene', smiles='c1ccccc1CCC=C'), + # ARCSpecies(label='peroxyl', smiles='CCO[O]')], + # p_species=[ARCSpecies(label='peroxide', smiles='CCOO'), + # ARCSpecies(label='butenylbenzene_rad', smiles='c1ccccc1[CH]CC=C')]) + # r_map, p_map = map_arc_rmg_species(rmg_reaction=rmg_reaction, + # arc_reaction=arc_reaction, + # concatenate=False, + # ) + # self.assertEqual(r_map, {0: 0, 1: 1}) + # self.assertEqual(p_map, {0: 1, 1: 0}) + # + # rmg_reaction = Reaction(reactants=[Species(smiles='[N]N'), Species(smiles='NN')], + # products=[Species(smiles='[NH]N'), Species(smiles='[NH]N')]) + # arc_reaction = ARCReaction(r_species=[ARCSpecies(label='H2NN(T)', smiles='[N]N'), + # ARCSpecies(label='N2H4', smiles='NN')], + # p_species=[ARCSpecies(label='N2H3', smiles='[NH]N'), + # ARCSpecies(label='N2H3', smiles='[NH]N')]) + # r_map, p_map = map_arc_rmg_species(rmg_reaction=rmg_reaction, + # arc_reaction=arc_reaction, + # concatenate=False, + # ) + # self.assertEqual(r_map, {0: 0, 1: 1}) + # self.assertEqual(p_map, {0: 0, 1: 1}) def test_find_equivalent_atoms_in_reactants_and_products(self): """Test the find_equivalent_atoms_in_reactants_and_products() function""" @@ -1155,7 +1078,6 @@ def test_check_species_before_mapping(self): verbose=True)) self.assertTrue(check_species_before_mapping(spc_1=self.spc1, spc_2=self.spc2, verbose=True)) - def test_get_bonds_dict(self): """Test the get_bonds_dict function.""" bond_dict = get_bonds_dict(spc=ARCSpecies(label='CH4', smiles='C')) @@ -1408,13 +1330,13 @@ def test_make_bond_changes(self): rxn = ARCReaction(r_species=[ARCSpecies(label="N4", smiles="NNNN")], p_species=[ARCSpecies(label="NH3", smiles="N"), ARCSpecies(label="NH2NHN_p", smiles="[N-]=[NH+]N")]) - rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=rxn, backend="arc") - r_label_dict, p_label_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=rxn, - rmg_reaction=rmg_reactions[0]) - assign_labels_to_products(rxn, p_label_dict) + # rmg_reactions = get_rmg_reactions_from_arc_reaction(arc_reaction=rxn, backend="arc") + # r_label_dict, p_label_dict = get_atom_indices_of_labeled_atoms_in_an_rmg_reaction(arc_reaction=rxn, + # rmg_reaction=rmg_reactions[0]) + assign_labels_to_products(rxn) reactants, products = copy_species_list_for_mapping(rxn.r_species), copy_species_list_for_mapping(rxn.p_species) label_species_atoms(reactants), label_species_atoms(products) - r_bdes, _ = find_all_bdes(rxn, r_label_dict, True), find_all_bdes(rxn, p_label_dict, False) + r_bdes, _ = find_all_bdes(rxn, True), find_all_bdes(rxn, False) r_cuts = cut_species_based_on_atom_indices(reactants, r_bdes) self.assertFalse(r_cuts[1].mol.is_isomorphic(rxn.p_species[1].mol)) make_bond_changes(rxn=rxn, @@ -1424,15 +1346,15 @@ def test_make_bond_changes(self): def test_update_xyz(self): """tests the update_xyz function""" - spc = ARCSpecies(label="test_UX",smiles = "BrC(F)Cl") + spc = ARCSpecies(label="test_UX", smiles="BrC(F)Cl") shuffle(spc.mol.atoms) update_xyz([spc]) xyz = spc.get_xyz()["symbols"] atoms = [atom.element.symbol for atom in spc.mol.atoms] for label1, label2 in zip(atoms, xyz): self.assertEqual(label1, label2) - - spc = ARCSpecies(label="test_UX", smiles = "OCl") + + spc = ARCSpecies(label="test_UX", smiles="OCl") shuffle(spc.mol.atoms) update_xyz([spc]) xyz = spc.get_xyz()["symbols"] @@ -1440,36 +1362,36 @@ def test_update_xyz(self): for label1, label2 in zip(atoms, xyz): self.assertEqual(label1, label2) - spc = ARCSpecies(label="test_UX", smiles = "BrOCl") + spc = ARCSpecies(label="test_UX", smiles="BrOCl") shuffle(spc.mol.atoms) update_xyz([spc]) xyz = spc.get_xyz()["symbols"] atoms = [atom.element.symbol for atom in spc.mol.atoms] - for label1,label2 in zip(atoms, xyz): + for label1, label2 in zip(atoms, xyz): self.assertEqual(label1, label2) def test_add_adjacent_hydrogen_atoms_to_map_based_on_a_specific_torsion(self): "test the add_adjacent_hydrogen_atoms_to_map_based_on_a_specific_torsion function" - spc1 = ARCSpecies(label="c4", smiles = "CCCC") + spc1 = ARCSpecies(label="c4", smiles="CCCC") spc2 = spc1.copy() out_dict = add_adjacent_hydrogen_atoms_to_map_based_on_a_specific_torsion(spc1, spc2, spc1.mol.atoms[0], spc2.mol.atoms[0], - [0,1,2,3], - {0: 0,1: 1,2: 2,3: 3}, + [0, 1, 2, 3], + {0: 0, 1: 1, 2: 2, 3: 3}, True) for key in range(7): self.assertEqual(out_dict[key], key) - - spc1 = ARCSpecies(label="c4", smiles = "cccc") + + spc1 = ARCSpecies(label="c4", smiles="cccc") spc2 = spc1.copy() out_dict = add_adjacent_hydrogen_atoms_to_map_based_on_a_specific_torsion(spc1, spc2, spc1.mol.atoms[0], spc2.mol.atoms[0], - [0,1,2,3], - {0: 0,1: 1,2: 2,3: 3}, + [0, 1, 2, 3], + {0: 0, 1: 1, 2: 2, 3: 3}, True) for key in range(4): self.assertEqual(out_dict[key], key) @@ -1480,8 +1402,8 @@ def test_add_adjacent_hydrogen_atoms_to_map_based_on_a_specific_torsion(self): spc2, spc1.mol.atoms[0], spc2.mol.atoms[0], - [0,1,2,3], - {0: 0,1: 1,2: 2,3: 3}, + [0, 1, 2, 3], + {0: 0, 1: 1, 2: 2, 3: 3}, True) for key in range(4): self.assertEqual(out_dict[key], key) @@ -1489,8 +1411,8 @@ def test_add_adjacent_hydrogen_atoms_to_map_based_on_a_specific_torsion(self): def test_find_all_bdes(self): """tests the find_all_bdes function""" rxn = ARCReaction(r_species=[self.r_1, self.r_2], p_species=[self.p_1, self.p_2]) - bdes = find_all_bdes(rxn, self.r_label_dict_rxn_1, True) - self.assertEqual(bdes, [(2, 4)]) + bdes = find_all_bdes(rxn, True) + self.assertEqual(bdes, [(4, 6)]) def test_copy_species_list_for_mapping(self): """tests the copy_species_list_for_mapping function""" @@ -1511,5 +1433,6 @@ def test_determine_bdes_on_spc_based_on_atom_labels(self): self.assertTrue(determine_bdes_on_spc_based_on_atom_labels(spc, (1, 2))) self.assertFalse(determine_bdes_on_spc_based_on_atom_labels(spc, (2, 3))) + if __name__ == '__main__': unittest.main(testRunner=unittest.TextTestRunner(verbosity=2)) \ No newline at end of file diff --git a/arc/reaction/family.py b/arc/reaction/family.py index 31d8584284..8d840d99f1 100644 --- a/arc/reaction/family.py +++ b/arc/reaction/family.py @@ -40,6 +40,8 @@ def __init__(self, label: str, consider_arc_families: bool = True, ): + if label is None: + raise ValueError('Cannot initialize a ReactionFamily object without a label') self.label = label self.groups_as_lines = self.get_groups_file_as_lines(consider_arc_families=consider_arc_families) self.reversible = is_reversible(self.groups_as_lines) @@ -370,7 +372,7 @@ def get_reactant_num(self) -> int: def get_reaction_family_products(rxn: 'ARCReaction', - rmg_family_set: str = 'default', # todo: consider giving a particular set of families + rmg_family_set: Optional[Union[List[str], str]] = None, consider_rmg_families: bool = True, consider_arc_families: bool = True, discover_own_reverse_rxns_in_reverse: bool = False, @@ -380,8 +382,11 @@ def get_reaction_family_products(rxn: 'ARCReaction', Args: rxn ('ARCReaction'): The ARC reaction object. - rmg_family_set (str, optional): The RMG family set to use from RMG-database/input/kinetics/families/recommended.py. - Note that surface families are excluded if 'all' is used. + rmg_family_set (Union[List[str], str], optional): The RMG family set to use from + RMG-database/input/kinetics/families/recommended.py. + Can be a name of a defined set, or a list + of explicit family lables to consider. + Note that surface families are excluded if 'all' is used. consider_rmg_families (bool, optional): Whether to consider RMG's families. consider_arc_families (bool, optional): Whether to consider ARC's custom families. discover_own_reverse_rxns_in_reverse (bool, optional): Whether to discover reactions belonging to a family @@ -395,7 +400,7 @@ def get_reaction_family_products(rxn: 'ARCReaction', List[dict]: The list of product dictionaries with the reaction family label. Keys are: 'family', 'group_labels', 'products', 'own_reverse', 'discovered_in_reverse', 'actions'. """ - family_labels = get_all_families(rmg_family_set=rmg_family_set, + family_labels = get_all_families(rmg_family_set=rmg_family_set or rxn.rmg_family_set, consider_rmg_families=consider_rmg_families, consider_arc_families=consider_arc_families) product_dicts = list() @@ -536,7 +541,7 @@ def check_product_isomorphism(products: List['Molecule'], # return all_reactions -def get_all_families(rmg_family_set: str = 'default', +def get_all_families(rmg_family_set: Union[List[str], str] = 'default', consider_rmg_families: bool = True, consider_arc_families: bool = True, ) -> List[str]: @@ -544,8 +549,7 @@ def get_all_families(rmg_family_set: str = 'default', Get all available RMG and ARC families. Args: - rmg_family_set (str, optional): The RMG family set to use from RMG-database/input/kinetics/families/recommended.py. - Note that surface families are excluded if 'all' is used. + rmg_family_set (Union[List[str], str], optional): The RMG family set to use. consider_rmg_families (bool, optional): Whether to consider RMG's families. consider_arc_families (bool, optional): Whether to consider ARC's custom families. @@ -555,14 +559,15 @@ def get_all_families(rmg_family_set: str = 'default', rmg_families, arc_families = list(), list() if consider_rmg_families: family_sets = get_rmg_recommended_family_sets() - if rmg_family_set not in list(family_sets.keys()) + ['all']: + if not isinstance(rmg_families, list) and rmg_family_set not in list(family_sets) + ['all']: raise ValueError(f'Invalid RMG family set: {rmg_family_set}') if rmg_family_set == 'all': for family_set_label, families in family_sets.items(): if 'surface' not in family_set_label: rmg_families.extend(list(families)) else: - rmg_families = list(family_sets[rmg_family_set]) + rmg_families = list(family_sets[rmg_family_set]) \ + if isinstance(rmg_family_set, str) and rmg_family_set in family_sets else rmg_family_set if consider_arc_families: arc_families = [os.path.splitext(family)[0] for family in os.listdir(ARC_FAMILIES_PATH)] return rmg_families + arc_families diff --git a/arc/reaction/family_test.py b/arc/reaction/family_test.py index bd64bcbff4..f80c3d6efa 100644 --- a/arc/reaction/family_test.py +++ b/arc/reaction/family_test.py @@ -59,7 +59,7 @@ def test_arc_families_path(self): self.assertIn('families', ARC_FAMILIES_PATH) self.assertTrue(os.path.isdir(ARC_FAMILIES_PATH)) - def test_get_reaction_family_products(self): + def test_get_reaction_family_products_abstraction(self): """Test determining the reaction family using product dicts""" rxn_0a = ARCReaction(r_species=[ARCSpecies(label='CH4', smiles='C'), ARCSpecies(label='O2', smiles='[O][O]')], p_species=[ARCSpecies(label='CH3', smiles='[CH3]'), @@ -130,7 +130,8 @@ def test_get_reaction_family_products(self): H -0.83821148 -0.26602407 0.00000000""" rxn_0b = ARCReaction(r_species=[ARCSpecies(label='CH4', xyz=ch4_xyz), ARCSpecies(label='O2', xyz=o2_xyz, multiplicity=3)], - p_species=[ARCSpecies(label='CH3', xyz=ch3_xyz), ARCSpecies(label='HO2', xyz=ho2_xyz)]) + p_species=[ARCSpecies(label='CH3', xyz=ch3_xyz), ARCSpecies(label='HO2', xyz=ho2_xyz)], + rmg_family_set=['H_Abstraction']) products = get_reaction_family_products(rxn_0b) expected_products = [{'family': 'H_Abstraction', 'group_labels': ('X_H', 'Y_rad'), 'products': [Molecule(smiles="[O]O"), Molecule(smiles="[CH3]")], @@ -168,7 +169,8 @@ def test_get_reaction_family_products(self): rxn_1 = ARCReaction(reactants=['NH2', 'NH2'], products=['NH', 'NH3'], r_species=[ARCSpecies(label='NH2', smiles='[NH2]')], - p_species=[ARCSpecies(label='NH', smiles='[NH]'), ARCSpecies(label='NH3', smiles='N')]) + p_species=[ARCSpecies(label='NH', smiles='[NH]'), ARCSpecies(label='NH3', smiles='N')], + rmg_family_set=['H_Abstraction']) products = get_reaction_family_products(rxn_1) expected_products = [{'discovered_in_reverse': False, 'family': 'H_Abstraction', @@ -196,6 +198,8 @@ def test_get_reaction_family_products(self): 'products': [Molecule(smiles="N"), Molecule(smiles="[NH]")]}] self.assertEqual(products, expected_products) + def test_get_reaction_family_products_cyclic_ether_formation(self): + """Test determining the reaction family using product dicts""" rxn_2f = ARCReaction(r_species=[ARCSpecies(label='C2H3O3', smiles='[CH2]C(=O)OO')], p_species=[ARCSpecies(label='C2H2O2', smiles='O=C1CO1'), ARCSpecies(label='OH', smiles='[OH]')]) @@ -216,7 +220,8 @@ def test_get_reaction_family_products(self): rxn_2_r = ARCReaction(p_species=[ARCSpecies(label='C2H3O3', smiles='[CH2]C(=O)OO')], r_species=[ARCSpecies(label='C2H2O2', smiles='O=C1CO1'), - ARCSpecies(label='OH', smiles='[OH]')]) + ARCSpecies(label='OH', smiles='[OH]')], + rmg_family_set=['Cyclic_Ether_Formation']) products = get_reaction_family_products(rxn_2_r) expected_products = [{'discovered_in_reverse': True, 'family': 'Cyclic_Ether_Formation', @@ -232,8 +237,11 @@ def test_get_reaction_family_products(self): 'products': [Molecule(smiles="[OH]"), Molecule(smiles="O=C1CO1")]}] self.assertEqual(products, expected_products) + def test_get_reaction_family_migration(self): + """Test determining the reaction family using product dicts""" rxn_3 = ARCReaction(r_species=[ARCSpecies(label='CCC[O]', smiles='CCC[O]')], - p_species=[ARCSpecies(label='C[CH]CO', smiles='C[CH]CO')]) + p_species=[ARCSpecies(label='C[CH]CO', smiles='C[CH]CO')], + rmg_family_set=['intra_H_migration']) products = get_reaction_family_products(rxn_3) expected_products = [{'discovered_in_reverse': False, 'family': 'intra_H_migration', @@ -249,8 +257,65 @@ def test_get_reaction_family_products(self): 'products': [Molecule(smiles="C[CH]CO")]}] self.assertEqual(products, expected_products) + rxn_8 = ARCReaction(r_species=[ARCSpecies(label='C3H7O2a', smiles='C[CH]COO')], + p_species=[ARCSpecies(label='C3H7O2b', smiles='CC(O)C[O]')], + rmg_family_set=['intra_OH_migration']) + products = get_reaction_family_products(rxn_8) + expected_products = [{'discovered_in_reverse': False, + 'family': 'intra_OH_migration', + 'group_labels': 'R2OOH', + 'label_map': {'*1': 1, '*2': 3, '*3': 4, '*4': 2}, + 'own_reverse': False, + 'products': [Molecule(smiles="CC(O)C[O]")]}] + self.assertEqual(products, expected_products) + + rxn_9 = ARCReaction(r_species=[ARCSpecies(label='R', smiles='C[CH]CCCC')], # iso-teleological paths + p_species=[ARCSpecies(label='P', smiles='[CH2]CCCCC')], + rmg_family_set=['intra_H_migration']) + products = get_reaction_family_products(rxn_9) + expected_products = [{'discovered_in_reverse': False, + 'family': 'intra_H_migration', + 'group_labels': 'R5Hall', + 'label_map': {'*1': 1, '*2': 5, '*3': 16, '*4': 2, '*5': 4, '*6': 3}, + 'own_reverse': True, + 'products': [Molecule(smiles="[CH2]CCCCC")]}, + {'discovered_in_reverse': False, + 'family': 'intra_H_migration', + 'group_labels': 'R5Hall', + 'label_map': {'*1': 1, '*2': 5, '*3': 17, '*4': 2, '*5': 4, '*6': 3}, + 'own_reverse': True, + 'products': [Molecule(smiles="[CH2]CCCCC")]}, + {'discovered_in_reverse': False, + 'family': 'intra_H_migration', + 'group_labels': 'R5Hall', + 'label_map': {'*1': 1, '*2': 5, '*3': 18, '*4': 2, '*5': 4, '*6': 3}, + 'own_reverse': True, + 'products': [Molecule(smiles="[CH2]CCCCC")]}, + {'discovered_in_reverse': False, + 'family': 'intra_H_migration', + 'group_labels': 'R2H', + 'label_map': {'*1': 1, '*2': 0, '*3': 6}, + 'own_reverse': True, + 'products': [Molecule(smiles="[CH2]CCCCC")]}, + {'discovered_in_reverse': False, + 'family': 'intra_H_migration', + 'group_labels': 'R2H', + 'label_map': {'*1': 1, '*2': 0, '*3': 7}, + 'own_reverse': True, + 'products': [Molecule(smiles="[CH2]CCCCC")]}, + {'discovered_in_reverse': False, + 'family': 'intra_H_migration', + 'group_labels': 'R2H', + 'label_map': {'*1': 1, '*2': 0, '*3': 8}, + 'own_reverse': True, + 'products': [Molecule(smiles="[CH2]CCCCC")]}] + self.assertEqual(products, expected_products) + + def test_get_reaction_family_ketoenol(self): + """Test determining the reaction family using product dicts""" rxn_4f = ARCReaction(r_species=[ARCSpecies(label='HOCCOH', smiles='OC=CO')], - p_species=[ARCSpecies(label='HOCCO', smiles='OCC=O')]) + p_species=[ARCSpecies(label='HOCCO', smiles='OCC=O')], + rmg_family_set=['Ketoenol']) products = get_reaction_family_products(rxn_4f) expected_products = [{'discovered_in_reverse': False, 'family': 'Ketoenol', @@ -267,7 +332,8 @@ def test_get_reaction_family_products(self): self.assertEqual(products, expected_products) rxn_4r = ARCReaction(r_species=[ARCSpecies(label='HOCCO', smiles='OCC=O')], - p_species=[ARCSpecies(label='HOCCOH', smiles='OC=CO')]) + p_species=[ARCSpecies(label='HOCCOH', smiles='OC=CO')], + rmg_family_set=['Ketoenol']) products = get_reaction_family_products(rxn_4r) expected_products = [{'discovered_in_reverse': True, 'family': 'Ketoenol', @@ -283,10 +349,13 @@ def test_get_reaction_family_products(self): 'products': [Molecule(smiles="O=CCO")]}] self.assertEqual(products, expected_products) + def test_get_reaction_family_insertion(self): + """Test determining the reaction family using product dicts""" # Test a family that is bimolecular, but the two reactants are defined as a single group in the recipe: rxn_5f = ARCReaction(r_species=[ARCSpecies(label='CO2', smiles='O=C=O'), ARCSpecies(label='H2', smiles='[H][H]')], - p_species=[ARCSpecies(label='CHOOH', smiles='O=CO')]) + p_species=[ARCSpecies(label='CHOOH', smiles='O=CO')], + rmg_family_set=['1,3_Insertion_CO2']) products = get_reaction_family_products(rxn_5f) expected_products = [{'discovered_in_reverse': False, 'family': '1,3_Insertion_CO2', @@ -316,7 +385,8 @@ def test_get_reaction_family_products(self): rxn_5r = ARCReaction(r_species=[ARCSpecies(label='CHOOH', smiles='O=CO')], p_species=[ARCSpecies(label='CO2', smiles='O=C=O'), - ARCSpecies(label='H2', smiles='[H][H]')]) + ARCSpecies(label='H2', smiles='[H][H]')], + rmg_family_set=['1,3_Insertion_CO2']) products = get_reaction_family_products(rxn_5r) expected_products = [{'discovered_in_reverse': True, 'family': '1,3_Insertion_CO2', @@ -344,9 +414,12 @@ def test_get_reaction_family_products(self): 'products': [Molecule(smiles="O=CO")]}] self.assertEqual(products, expected_products) + def test_get_reaction_family_elimination(self): + """Test determining the reaction family using product dicts""" rxn_6 = ARCReaction(r_species=[ARCSpecies(label='C2H3O3', smiles='CC(=O)O[O]')], p_species=[ARCSpecies(label='C2H2O', smiles='C=C=O'), - ARCSpecies(label='HO2', smiles='O[O]')]) + ARCSpecies(label='HO2', smiles='O[O]')], + rmg_family_set=['HO2_Elimination_from_PeroxyRadical']) products = get_reaction_family_products(rxn_6) expected_products = [{'discovered_in_reverse': False, 'family': 'HO2_Elimination_from_PeroxyRadical', @@ -368,11 +441,14 @@ def test_get_reaction_family_products(self): 'products': [Molecule(smiles="[O]O"), Molecule(smiles="C=C=O")]}] self.assertEqual(products, expected_products) + def test_get_reaction_family_r_add_mult_bond(self): + """Test determining the reaction family using product dicts""" spc_1 = ARCSpecies(label='C2H3O3', smiles='[CH2]C(=O)OO') spc_1.mol_list = generate_resonance_structures(object_=spc_1.mol, keep_isomorphic=True) rxn_7 = ARCReaction(r_species=[spc_1], p_species=[ARCSpecies(label='C2H2O', smiles='C=C=O'), - ARCSpecies(label='HO2', smiles='O[O]')]) + ARCSpecies(label='HO2', smiles='O[O]')], + rmg_family_set=['R_Addition_MultipleBond']) products = get_reaction_family_products(rxn_7) expected_products = [{'discovered_in_reverse': True, 'family': 'R_Addition_MultipleBond', @@ -388,64 +464,13 @@ def test_get_reaction_family_products(self): 'products': [Molecule(smiles="C=C([O])OO")]}] self.assertEqual(products, expected_products) - rxn_8 = ARCReaction(r_species=[ARCSpecies(label='C3H7O2a', smiles='C[CH]COO')], - p_species=[ARCSpecies(label='C3H7O2b', smiles='CC(O)C[O]')]) - products = get_reaction_family_products(rxn_8) - expected_products = [{'discovered_in_reverse': False, - 'family': 'intra_OH_migration', - 'group_labels': 'R2OOH', - 'label_map': {'*1': 1, '*2': 3, '*3': 4, '*4': 2}, - 'own_reverse': False, - 'products': [Molecule(smiles="CC(O)C[O]")]}] - self.assertEqual(products, expected_products) - - rxn_9 = ARCReaction(r_species=[ARCSpecies(label='R', smiles='C[CH]CCCC')], # iso-teleological paths - p_species=[ARCSpecies(label='P', smiles='[CH2]CCCCC')]) - products = get_reaction_family_products(rxn_9) - expected_products = [{'discovered_in_reverse': False, - 'family': 'intra_H_migration', - 'group_labels': 'R5Hall', - 'label_map': {'*1': 1, '*2': 5, '*3': 16, '*4': 2, '*5': 4, '*6': 3}, - 'own_reverse': True, - 'products': [Molecule(smiles="[CH2]CCCCC")]}, - {'discovered_in_reverse': False, - 'family': 'intra_H_migration', - 'group_labels': 'R5Hall', - 'label_map': {'*1': 1, '*2': 5, '*3': 17, '*4': 2, '*5': 4, '*6': 3}, - 'own_reverse': True, - 'products': [Molecule(smiles="[CH2]CCCCC")]}, - {'discovered_in_reverse': False, - 'family': 'intra_H_migration', - 'group_labels': 'R5Hall', - 'label_map': {'*1': 1, '*2': 5, '*3': 18, '*4': 2, '*5': 4, '*6': 3}, - 'own_reverse': True, - 'products': [Molecule(smiles="[CH2]CCCCC")]}, - {'discovered_in_reverse': False, - 'family': 'intra_H_migration', - 'group_labels': 'R2H', - 'label_map': {'*1': 1, '*2': 0, '*3': 6}, - 'own_reverse': True, - 'products': [Molecule(smiles="[CH2]CCCCC")]}, - {'discovered_in_reverse': False, - 'family': 'intra_H_migration', - 'group_labels': 'R2H', - 'label_map': {'*1': 1, '*2': 0, '*3': 7}, - 'own_reverse': True, - 'products': [Molecule(smiles="[CH2]CCCCC")]}, - {'discovered_in_reverse': False, - 'family': 'intra_H_migration', - 'group_labels': 'R2H', - 'label_map': {'*1': 1, '*2': 0, '*3': 8}, - 'own_reverse': True, - 'products': [Molecule(smiles="[CH2]CCCCC")]}] - self.assertEqual(products, expected_products) - def test_get_reaction_family_products_considering_resonance(self): """Test determining the reaction family using product dicts for the non-trivial resonance structure""" rxn_0a = ARCReaction(r_species=[ARCSpecies(label='R1', smiles='CC(C=S)S'), ARCSpecies(label='R2', smiles='C=C[CH]C')], p_species=[ARCSpecies(label='P1', smiles='[S]C=C(S)C'), - ARCSpecies(label='P2', smiles='C=CCC')]) + ARCSpecies(label='P2', smiles='C=CCC')], + rmg_family_set=['H_Abstraction']) products = get_reaction_family_products(rxn_0a) expected_products = [{'family': 'H_Abstraction', 'group_labels': ('X_H', 'Y_rad'), 'products': [Molecule(smiles="C=CCC"), Molecule(smiles="C[C](C=S)S")], @@ -460,7 +485,8 @@ def test_get_reaction_family_products_considering_resonance(self): def test_determine_possible_reaction_products_from_family(self): """Test determining the possible reaction products from a family""" rxn_1 = ARCReaction(r_species=[ARCSpecies(label='nC3H7', smiles='[CH2]CC')], - p_species=[ARCSpecies(label='iC3H7', smiles='C[CH]C')]) + p_species=[ARCSpecies(label='iC3H7', smiles='C[CH]C')], + rmg_family_set=['intra_H_migration']) product_dicts = determine_possible_reaction_products_from_family(rxn_1, family_label='intra_H_migration') expected_product_dicts = [{'discovered_in_reverse': False, 'family': 'intra_H_migration', diff --git a/arc/reaction/reaction.py b/arc/reaction/reaction.py index 85331a96f3..d019ba6da7 100644 --- a/arc/reaction/reaction.py +++ b/arc/reaction/reaction.py @@ -43,7 +43,10 @@ class ARCReaction(object): r_species (List[ARCSpecies], optional): A list of reactants :ref:`ARCSpecies ` objects. p_species (List[ARCSpecies], optional): A list of products :ref:`ARCSpecies ` objects. ts_label (str, optional): The :ref:`ARCSpecies ` label of the respective TS. - rmg_reaction (Reaction, optional): An RMG Reaction class. + rmg_family_set (Union[List[str], str], optional): The RMG family set to which the reaction belongs + (corresponds to sets in db families recommended.py) + or an explicit names of families to consider + given as a list of strings. ts_xyz_guess (list, optional): A list of TS XYZ user guesses, each in a string format. xyz (list, optional): Identical to `ts_xyz_guess`, used as a shortcut. multiplicity (int, optional): The reaction surface multiplicity. A trivial guess will be made unless provided. @@ -68,7 +71,7 @@ class ARCReaction(object): dh_rxn298 (float): The heat of reaction at 298K. kinetics (Arrhenius): The high pressure limit rate coefficient calculated by ARC. rmg_kinetics (Arrhenius): The kinetics generated by RMG, for reality-check. - rmg_reaction (Reaction): An RMG Reaction class. + rmg_family_set (Union[List[str]]): The RMG family set or a list of specific families to consider. rmg_reactions (list): A list of RMG Reaction objects with RMG rates for comparisons. long_kinetic_description (str): A description for the species entry in the thermo library outputted. ts_xyz_guess (list): A list of TS XYZ user guesses, each in a string format. @@ -91,7 +94,7 @@ def __init__(self, r_species: Optional[List[ARCSpecies]] = None, p_species: Optional[List[ARCSpecies]] = None, ts_label: Optional[str] = None, - rmg_reaction: Optional[Reaction] = None, + rmg_family_set: Union[List[str]] = 'default', ts_xyz_guess: Optional[list] = None, xyz: Optional[list] = None, multiplicity: Optional[int] = None, @@ -106,6 +109,7 @@ def __init__(self, self.p_species = p_species or list() self.kinetics = None self.rmg_kinetics = None + self.rmg_family_set = rmg_family_set self.long_kinetic_description = '' self._family = None self._family_own_reverse = False @@ -129,10 +133,8 @@ def __init__(self, self.reactants = [check_label(reactant)[0] for reactant in reactants] if reactants else list() products = products or [spc.label for spc in self.p_species] or None self.products = [check_label(product)[0] for product in products] if products else list() - self.rmg_reaction = rmg_reaction - if self.rmg_reaction is None and not self.label \ - and not (len(self.reactants) * len(self.products)) \ - and not (len(self.r_species) * len(self.p_species)): + if (not self.label and not (len(self.reactants) * len(self.products)) + and not (len(self.r_species) * len(self.p_species))): raise InputError(f'Cannot determine reactants and/or products labels for reaction {self.label}') self.set_label_reactants_products() self.ts_xyz_guess = ts_xyz_guess or list() @@ -142,7 +144,6 @@ def __init__(self, f'reactants and {len(self.products)} products for reaction {self.label}.') if not isinstance(self.ts_xyz_guess, list): self.ts_xyz_guess = [self.ts_xyz_guess] - self.arc_species_from_rmg_reaction() self.remove_dup_species() self.check_atom_balance() @@ -172,9 +173,6 @@ def charge(self): if self._charge is None: if len(self.r_species): self._charge = self.get_rxn_charge() - elif self.rmg_reaction is not None: - self._charge = sum([spc.molecule[0].get_net_charge() for spc in - self.rmg_reaction.reactants + self.rmg_reaction.products]) return self._charge @charge.setter @@ -232,7 +230,8 @@ def __str__(self) -> str: """Return a string representation of the object""" str_representation = f'ARCReaction(' str_representation += f'label="{self.label}", ' - str_representation += f'rmg_reaction="{self.rmg_reaction}", ' + if self.family is not None: + str_representation += f'family="{self.family}", ' if self.preserve_param_in_scan is not None: str_representation += f'preserve_param_in_scan="{self.preserve_param_in_scan}", ' str_representation += f'multiplicity={self.multiplicity}, ' @@ -274,8 +273,8 @@ def as_dict(self, reaction_dict['done_opt_r_n_p'] = self.done_opt_r_n_p if self.preserve_param_in_scan is not None: reaction_dict['preserve_param_in_scan'] = self.preserve_param_in_scan - if 'rmg_reaction' in reaction_dict: - reaction_dict['rmg_reaction'] = self.rmg_reaction_to_str() + if 'rmg_family_set' in reaction_dict and reaction_dict['rmg_family_set'] != 'default': + reaction_dict['rmg_family_set'] = self.rmg_family_set if report_family: if self.family is not None: reaction_dict['family'] = self.family @@ -306,21 +305,10 @@ def from_dict(self, if 'family' in reaction_dict and reaction_dict['family'] is not None: self.family = reaction_dict['family'] self.family_own_reverse = reaction_dict['family_own_reverse'] if 'family_own_reverse' in reaction_dict else False - if 'rmg_reaction' in reaction_dict: - self.rmg_reaction_from_str(reaction_string=reaction_dict['rmg_reaction']) - else: - self.rmg_reaction = None + self.rmg_family_set = reaction_dict['rmg_family_set'] if 'rmg_family_set' in reaction_dict else 'default' self.set_label_reactants_products(species_list) - if self.rmg_reaction is None and not (len(self.reactants) * len(self.products)): - raise InputError(f'Cannot determine reactants and/or products labels for reaction {self.label}') if not (len(self.reactants) * len(self.products)): - if not all([spc.label for spc in self.rmg_reaction.reactants + self.rmg_reaction.products]): - raise InputError(f'All species in a reaction must be labeled (and the labels must correspond ' - f'to respective Species in ARC). If an RMG Reaction object was passes, make ' - f'sure that all species in the reactants and products are correctly labeled. ' - f'Problematic reaction: {self.label}') - self.reactants = [check_label(spc.label)[0] for spc in self.rmg_reaction.reactants] - self.products = [check_label(spc.label)[0] for spc in self.rmg_reaction.products] + raise InputError(f'Cannot determine reactants and/or products labels for reaction {self.label}') if self.ts_label is None: self.ts_label = reaction_dict['ts_label'] if 'ts_label' in reaction_dict else None if species_list is not None and 'r_species' in reaction_dict and len(reaction_dict['r_species']) \ @@ -371,8 +359,7 @@ def flip_reaction(self, report_family: bool = True): ARCReaction: A copy of this object instance with flipped reactants and products. """ reaction_dict = self.as_dict(reset_atom_ids=True, report_family=report_family) - reset_keys = ['label', 'index', 'atom_map', 'rmg_reaction', - 'family', 'family_own_reverse', 'long_kinetic_description'] + reset_keys = ['label', 'index', 'atom_map', 'family', 'family_own_reverse', 'long_kinetic_description'] if 'r_species' in reaction_dict.keys() and 'p_species' in reaction_dict.keys(): reaction_dict['r_species'], reaction_dict['p_species'] = reaction_dict['p_species'], reaction_dict['r_species'] else: @@ -388,7 +375,6 @@ def flip_reaction(self, report_family: bool = True): reaction_dict['label'] = self.arrow.join([label_splits[1], label_splits[0]]) flipped_rxn = ARCReaction(reaction_dict=reaction_dict) flipped_rxn.set_label_reactants_products() - flipped_rxn.rmg_reaction_from_arc_species() return flipped_rxn def is_isomerization(self): @@ -433,9 +419,6 @@ def set_label_reactants_products(self, species_list: Optional[List[ARCSpecies]] if not len(self.p_species) and len(self.products): self.p_species = [spc for spc in species_list if spc.label in self.products or spc.original_label in self.products] - elif self.rmg_reaction is not None: - self.reactants = [r.label for r in self.rmg_reaction.reactants] - self.products = [p.label for p in self.rmg_reaction.products] elif len(self.r_species) and len(self.p_species): self.reactants = [r.label for r in self.r_species] self.products = [p.label for p in self.p_species] @@ -446,12 +429,6 @@ def set_label_reactants_products(self, species_list: Optional[List[ARCSpecies]] elif len(self.r_species) and len(self.p_species): self.label = self.arrow.join([self.plus.join(r.label for r in self.r_species), self.plus.join(p.label for p in self.p_species)]) - elif self.rmg_reaction is not None: - # This will probably never be executed, but OK to keep. - self.label = self.arrow.join([self.plus.join(r.label for r in self.rmg_reaction.reactants), - self.plus.join(p.label for p in self.rmg_reaction.products)]) - if self.rmg_reaction is None: - self.rmg_reaction_from_arc_species() elif not self.label and not (len(self.reactants) * len(self.products)): raise ReactionError('Either a label or reactants and products lists must be specified') self.reactants = [check_label(reactant)[0] for reactant in self.reactants] @@ -460,51 +437,6 @@ def set_label_reactants_products(self, species_list: Optional[List[ARCSpecies]] raise ReactionError(f'Both the reactants and products must be specified for a reaction, ' f'got: reactants = {self.reactants}, products = {self.products}.') - def rmg_reaction_to_str(self) -> str: - """A helper function for dumping the RMG Reaction object as a string for the YAML restart dictionary""" - return self.arrow.join([self.plus.join(r.molecule[0].copy(deep=True).to_smiles() - for r in self.rmg_reaction.reactants), - self.plus.join(p.molecule[0].copy(deep=True).to_smiles() - for p in self.rmg_reaction.products)]) - - def rmg_reaction_from_str(self, reaction_string: str): - """A helper function for regenerating the RMG Reaction object from a string representation""" - reactants, products = reaction_string.split(self.arrow) - reactants = [Species(smiles=smiles) for smiles in reactants.split(self.plus)] - products = [Species(smiles=smiles) for smiles in products.split(self.plus)] - self.rmg_reaction = Reaction(reactants=reactants, products=products) - - def rmg_reaction_from_arc_species(self): - """ - A helper function for generating the RMG Reaction object from ARCSpecies - Used for determining the family. - """ - if self.rmg_reaction is None and len(self.r_species) and len(self.p_species) and \ - all([arc_spc.mol is not None for arc_spc in self.r_species + self.p_species]): - reactants, products = self.get_reactants_and_products(arc=False) # Returns RMG Species. - self.rmg_reaction = Reaction(reactants=reactants, products=products) - - def arc_species_from_rmg_reaction(self): - """ - A helper function for generating the ARC Species (.r_species and .p_species) from the RMG Reaction object - """ - if self.rmg_reaction is not None and not len(self.r_species) and not len(self.p_species): - self.r_species, self.p_species = list(), list() - for i, rmg_reactant in enumerate(self.rmg_reaction.reactants): - if len(self.reactants) > i: - label = self.reactants[i] - else: - label = rmg_reactant.label or rmg_reactant.molecule[0].copy(deep=True).to_smiles() - label = check_label(label)[0] - self.r_species.append(ARCSpecies(label=label, mol=rmg_reactant.molecule[0])) - for i, rmg_product in enumerate(self.rmg_reaction.products): - if len(self.products) > i: - label = self.products[i] - else: - label = rmg_product.label or rmg_product.molecule[0].copy(deep=True).to_smiles() - label = check_label(label)[0] - self.p_species.append(ARCSpecies(label=label, mol=rmg_product.molecule[0])) - def get_rxn_charge(self): """A helper function for determining the surface charge""" if len(self.r_species): @@ -577,7 +509,6 @@ def get_rxn_multiplicity(self): return multiplicity def determine_family(self, - rmg_family_set: str = 'default', consider_rmg_families: bool = True, consider_arc_families: bool = True, discover_own_reverse_rxns_in_reverse: bool = False, @@ -587,22 +518,20 @@ def determine_family(self, Populates the .family, and .family_own_reverse attributes. Args: - rmg_family_set (str, optional): The RMG family set to use. consider_rmg_families (bool, optional): Whether to consider RMG's families in addition to ARC's. consider_arc_families (bool, optional): Whether to consider ARC's families in addition to RMG's. discover_own_reverse_rxns_in_reverse (bool, optional): Whether to discover own reverse reactions in reverse. """ - if self.rmg_reaction is not None: - product_dicts = get_reaction_family_products(rxn=self, - rmg_family_set=rmg_family_set, - consider_rmg_families=consider_rmg_families, - consider_arc_families=consider_arc_families, - discover_own_reverse_rxns_in_reverse=discover_own_reverse_rxns_in_reverse, - ) - if len(product_dicts): - family, family_own_reverse = product_dicts[0]['family'], product_dicts[0]['own_reverse'] - return family, family_own_reverse - return None, None + product_dicts = get_reaction_family_products(rxn=self, + rmg_family_set=self.rmg_family_set, + consider_rmg_families=consider_rmg_families, + consider_arc_families=consider_arc_families, + discover_own_reverse_rxns_in_reverse=discover_own_reverse_rxns_in_reverse, + ) + if len(product_dicts): + family, family_own_reverse = product_dicts[0]['family'], product_dicts[0]['own_reverse'] + return family, family_own_reverse + return None, None def check_attributes(self): """Check that the Reaction object is defined correctly""" diff --git a/arc/reaction/reaction_test.py b/arc/reaction/reaction_test.py index 3b3a00143f..73dbc01b51 100644 --- a/arc/reaction/reaction_test.py +++ b/arc/reaction/reaction_test.py @@ -95,24 +95,14 @@ def setUpClass(cls): H -0.36439762 -1.57761595 -0.41622918 H -1.43807526 1.62776079 0.73816131 H -1.28677889 1.04716138 -1.01532486""" - cls.rxn1 = ARCReaction(reactants=['CH4', 'OH'], products=['CH3', 'H2O'], - rmg_reaction=Reaction(reactants=[Species(label='CH4', smiles='C'), - Species(label='OH', smiles='[OH]')], - products=[Species(label='CH3', smiles='[CH3]'), - Species(label='H2O', smiles='O')])) - cls.rxn2 = ARCReaction(reactants=['C2H5', 'OH'], products=['C2H4', 'H2O'], - rmg_reaction=Reaction(reactants=[Species(label='C2H5', smiles='C[CH2]'), - Species(label='OH', smiles='[OH]')], - products=[Species(label='C2H4', smiles='C=C'), - Species(label='H2O', smiles='O')])) - cls.rxn3 = ARCReaction(reactants=['CH3CH2NH'], products=['CH2CH2NH2'], - rmg_reaction=Reaction(reactants=[Species(label='CH3CH2NH', smiles='CC[NH]')], - products=[Species(label='CH2CH2NH2', smiles='[CH2]CN')])) - cls.rxn4 = ARCReaction(reactants=['NH2', 'NH2NH'], products=['N', 'NH2N'], - rmg_reaction=Reaction(reactants=[Species(label='NH2', smiles='[NH2]'), - Species(label='NH2NH', smiles='N[NH]')], - products=[Species(label='N', smiles='N'), - Species(label='NH2N', smiles='N[N]')])) + cls.rxn1 = ARCReaction(r_species=[ARCSpecies(label='CH4', smiles='C'), ARCSpecies(label='OH', smiles='[OH]')], + p_species=[ARCSpecies(label='CH3', smiles='[CH3]'), ARCSpecies(label='H2O', smiles='O')]) + cls.rxn2 = ARCReaction(r_species=[ARCSpecies(label='C2H5', smiles='C[CH2]'), ARCSpecies(label='OH', smiles='[OH]')], + p_species=[ARCSpecies(label='C=C', smiles='C2H4'), ARCSpecies(label='H2O', smiles='O')]) + cls.rxn3 = ARCReaction(r_species=[ARCSpecies(label='CC[NH]', smiles='CH3CH2NH')], + p_species=[ARCSpecies(label='[CH2]CH', smiles='CH2CH2NH2')]) + cls.rxn4 = ARCReaction(r_species=[ARCSpecies(label='NH2', smiles='[NH2]'), ARCSpecies(label='N[NH]', smiles='NH2NH')], + p_species=[ARCSpecies(label='NH3', smiles='N'), ARCSpecies(label='N[N]', smiles='NH2N')]) cls.rxn5 = ARCReaction(reactants=['NH2', 'NH2'], products=['NH', 'NH3'], r_species=[ARCSpecies(label='NH2', smiles='[NH2]')], p_species=[ARCSpecies(label='NH', smiles='[NH]'), @@ -142,10 +132,10 @@ def setUpClass(cls): cls.rxn9 = ARCReaction(r_species=[ARCSpecies(label='NCOO', smiles='NCO[O]')], p_species=[ARCSpecies(label='NHCH2', smiles='N=C'), ARCSpecies(label='HO2', smiles='O[O]')]) - cls.rxn10 = ARCReaction(rmg_reaction=Reaction(reactants=[Species(label='HNO', smiles='N=O'), - Species(label='NO2', smiles='O=[N+][O-]')], - products=[Species(label='HNO2', smiles='[O-][NH+]=O'), - Species(label='NO', smiles='[N]=O')])) + cls.rxn10 = ARCReaction(r_species=[ARCSpecies(label='HNO', smiles='N=O'), + ARCSpecies(label='NO2', smiles='O=[N+][O-]')], + p_species=[ARCSpecies(label='HNO2', smiles='[O-][NH+]=O'), + ARCSpecies(label='NO', smiles='[N]=O')]) cls.rxn11 = ARCReaction(r_species=[ARCSpecies(label='C[CH]C', smiles='C[CH]C', xyz=cls.ch3chch3_xyz)], p_species=[ARCSpecies(label='[CH2]CC', smiles='[CH2]CC', xyz=cls.ch3ch2ch2_xyz)]) @@ -368,30 +358,6 @@ def test_is_isomerization(self): self.assertFalse(self.rxn8.is_isomerization()) self.assertFalse(self.rxn9.is_isomerization()) - def test_from_rmg_reaction(self): - """Test setting up an ARCReaction from an RMG Reaction""" - rmg_rxn_1 = Reaction(reactants=[Species(label='nC3H7', smiles='[CH2]CC')], - products=[Species(label='iC3H7', smiles='C[CH]C')]) - rxn_1 = ARCReaction(rmg_reaction=rmg_rxn_1) - self.assertEqual(rxn_1.label, 'nC3H7 <=> iC3H7') - - rmg_rxn_2 = Reaction(reactants=[Species(label='OH', smiles='[OH]'), Species(label='OH', smiles='[OH]')], - products=[Species(label='O', smiles='[O]'), Species(label='H2O', smiles='O')]) - rxn_2 = ARCReaction(rmg_reaction=rmg_rxn_2) - self.assertEqual(rxn_2.label, 'OH + OH <=> O + H2O') - - def test_rmg_reaction_to_str(self): - """Test the rmg_reaction_to_str() method and the reaction label generated""" - spc1 = Species().from_smiles('CON=O') - spc1.label = 'CONO' - spc2 = Species().from_smiles('C[N+](=O)[O-]') - spc2.label = 'CNO2' - rmg_reaction = Reaction(reactants=[spc1], products=[spc2]) - rxn = ARCReaction(rmg_reaction=rmg_reaction) - rxn_str = rxn.rmg_reaction_to_str() - self.assertEqual(rxn_str, 'CON=O <=> [O-][N+](=O)C') - self.assertEqual(rxn.label, 'CONO <=> CNO2') - def test_rxn_family(self): """Test that ARC gets the correct RMG family for different reactions""" self.assertEqual(self.rxn1.family, 'H_Abstraction') @@ -434,7 +400,6 @@ def test_charge_property(self): def test_multiplicity_property(self): """Test determining multiplicity""" self.assertEqual(self.rxn1.multiplicity, 2) - self.rxn2.arc_species_from_rmg_reaction() self.assertEqual(self.rxn2.multiplicity, 1) self.assertEqual(self.rxn3.multiplicity, 2) self.assertEqual(self.rxn4.multiplicity, 3) @@ -564,7 +529,6 @@ def test_get_species_count(self): def test_get_reactants_and_products(self): """Test getting reactants and products""" - self.rxn1.arc_species_from_rmg_reaction() self.rxn1.remove_dup_species() reactants, products = self.rxn1.get_reactants_and_products(arc=True) for spc in reactants + products: diff --git a/arc/scheduler.py b/arc/scheduler.py index f5c127a1be..ce33dd89f5 100644 --- a/arc/scheduler.py +++ b/arc/scheduler.py @@ -324,7 +324,6 @@ def __init__(self, rxn.r_species.append(spc) if spc.label in rxn.products: rxn.p_species.append(spc) - rxn.rmg_reaction_from_arc_species() rxn.check_attributes() family_text = '' if rxn.family is not None: