Merge pull request ReactionMechanismGenerator#2632 from ReactionMecha…

…nismGenerator/fix_generate_reactions

Clean up some code for calculating degeneracy in kinetics codebase

rmgpy/data/kinetics/common.py

@@ -222,9 +222,55 @@ def independent_ids():
         for species in input_species:
             species.generate_resonance_structures(keep_isomorphic=True)
 
+def check_for_same_reactants(reactants):
+    """
+    Given a list reactants, check if the reactants are the same.
+    If they refer to the same memory address, then make a deep copy so they can be manipulated independently.
+    
+    Returns a tuple containing the modified reactants list, and an integer containing the number of identical reactants in the reactants list. 
+    
+    """
+
+    same_reactants = 0
+    if len(reactants) == 2:
+        if reactants[0] is reactants[1]:
+            reactants[1] = reactants[1].copy(deep=True)
+            same_reactants = 2
+        elif reactants[0].is_isomorphic(reactants[1]):
+            same_reactants = 2
+    elif len(reactants) == 3:
+        same_01 = reactants[0] is reactants[1]
+        same_02 = reactants[0] is reactants[2]
+        if same_01 and same_02:
+            same_reactants = 3
+            reactants[1] = reactants[1].copy(deep=True)
+            reactants[2] = reactants[2].copy(deep=True)
+        elif same_01:
+            same_reactants = 2
+            reactants[1] = reactants[1].copy(deep=True)
+        elif same_02:
+            same_reactants = 2
+            reactants[2] = reactants[2].copy(deep=True)
+        elif reactants[1] is reactants[2]:
+            same_reactants = 2
+            reactants[2] = reactants[2].copy(deep=True)
+        else:
+            same_01 = reactants[0].is_isomorphic(reactants[1])
+            same_02 = reactants[0].is_isomorphic(reactants[2])
+            if same_01 and same_02:
+                same_reactants = 3
+            elif same_01 or same_02:
+                same_reactants = 2
+            elif reactants[1].is_isomorphic(reactants[2]):
+                same_reactants = 2
+    elif len(reactants) > 3:
+        raise ValueError('Cannot check for duplicate reactants if provided number of reactants is greater than 3. ' 
+                         'Got: {} reactants'.format(len(reactants))) 
+
+    return reactants, same_reactants
 
 def find_degenerate_reactions(rxn_list, same_reactants=None, template=None, kinetics_database=None,
-                              kinetics_family=None, save_order=False):
+                              kinetics_family=None, save_order=False, resonance=True):
     """
     Given a list of Reaction objects, this method combines degenerate
     reactions and increments the reaction degeneracy value. For multiple
@@ -250,6 +296,7 @@ def find_degenerate_reactions(rxn_list, same_reactants=None, template=None, kine
         kinetics_database (KineticsDatabase, optional): provide a KineticsDatabase instance for calculating degeneracy
         kinetics_family (KineticsFamily, optional):     provide a KineticsFamily instance for calculating degeneracy
         save_order (bool, optional):                    reset atom order after performing atom isomorphism
+        resonance (bool, optional):                     whether to consider resonance when computing degeneracy 
 
     Returns:
         Reaction list with degenerate reactions combined with proper degeneracy values
@@ -340,7 +387,7 @@ def find_degenerate_reactions(rxn_list, same_reactants=None, template=None, kine
                 from rmgpy.data.rmg import get_db
                 family = get_db('kinetics').families[rxn.family]
             if not family.own_reverse:
-                rxn.degeneracy = family.calculate_degeneracy(rxn)
+                rxn.degeneracy = family.calculate_degeneracy(rxn, resonance=resonance)
 
     return rxn_list
 

rmgpy/data/kinetics/database.py

@@ -37,7 +37,8 @@
 import rmgpy.constants as constants
 from rmgpy.data.base import LogicNode
 from rmgpy.data.kinetics.common import ensure_species, generate_molecule_combos, \
-                                       find_degenerate_reactions, ensure_independent_atom_ids
+                                       find_degenerate_reactions, ensure_independent_atom_ids, \
+                                       check_for_same_reactants
 from rmgpy.data.kinetics.family import KineticsFamily
 from rmgpy.data.kinetics.library import LibraryReaction, KineticsLibrary
 from rmgpy.exceptions import DatabaseError
@@ -427,7 +428,7 @@ def generate_reactions(self, reactants, products=None, only_families=None, reson
         if only_families is None:
             reaction_list.extend(self.generate_reactions_from_libraries(reactants, products))
         reaction_list.extend(self.generate_reactions_from_families(reactants, products,
-                                                                   only_families=None, resonance=resonance))
+                                                                   only_families=only_families, resonance=resonance))
         return reaction_list
 
     def generate_reactions_from_libraries(self, reactants, products=None):
@@ -487,43 +488,10 @@ def generate_reactions_from_families(self, reactants, products=None, only_famili
         Returns:
             List of reactions containing Species objects with the specified reactants and products.
         """
-        # Check if the reactants are the same
-        # If they refer to the same memory address, then make a deep copy so
-        # they can be manipulated independently
         if isinstance(reactants, tuple):
             reactants = list(reactants)
-        same_reactants = 0
-        if len(reactants) == 2:
-            if reactants[0] is reactants[1]:
-                reactants[1] = reactants[1].copy(deep=True)
-                same_reactants = 2
-            elif reactants[0].is_isomorphic(reactants[1]):
-                same_reactants = 2
-        elif len(reactants) == 3:
-            same_01 = reactants[0] is reactants[1]
-            same_02 = reactants[0] is reactants[2]
-            if same_01 and same_02:
-                same_reactants = 3
-                reactants[1] = reactants[1].copy(deep=True)
-                reactants[2] = reactants[2].copy(deep=True)
-            elif same_01:
-                same_reactants = 2
-                reactants[1] = reactants[1].copy(deep=True)
-            elif same_02:
-                same_reactants = 2
-                reactants[2] = reactants[2].copy(deep=True)
-            elif reactants[1] is reactants[2]:
-                same_reactants = 2
-                reactants[2] = reactants[2].copy(deep=True)
-            else:
-                same_01 = reactants[0].is_isomorphic(reactants[1])
-                same_02 = reactants[0].is_isomorphic(reactants[2])
-                if same_01 and same_02:
-                    same_reactants = 3
-                elif same_01 or same_02:
-                    same_reactants = 2
-                elif reactants[1].is_isomorphic(reactants[2]):
-                    same_reactants = 2
+
+        reactants, same_reactants = check_for_same_reactants(reactants)
 
         # Label reactant atoms for proper degeneracy calculation (cannot be in tuple)
         ensure_independent_atom_ids(reactants, resonance=resonance)
@@ -536,7 +504,8 @@ def generate_reactions_from_families(self, reactants, products=None, only_famili
                                                       prod_resonance=resonance))
 
         # Calculate reaction degeneracy
-        reaction_list = find_degenerate_reactions(reaction_list, same_reactants, kinetics_database=self)
+        reaction_list = find_degenerate_reactions(reaction_list, same_reactants, kinetics_database=self,
+                                                  resonance=resonance)
         # Add reverse attribute to families with ownReverse
         to_delete = []
         for i, rxn in enumerate(reaction_list):

rmgpy/data/kinetics/family.py

@@ -48,7 +48,7 @@
 from rmgpy.constraints import fails_species_constraints
 from rmgpy.data.base import Database, Entry, LogicNode, LogicOr, ForbiddenStructures, get_all_combinations
 from rmgpy.data.kinetics.common import save_entry, find_degenerate_reactions, generate_molecule_combos, \
-                                       ensure_independent_atom_ids
+                                       ensure_independent_atom_ids, check_for_same_reactants
 from rmgpy.data.kinetics.depository import KineticsDepository
 from rmgpy.data.kinetics.groups import KineticsGroups
 from rmgpy.data.kinetics.rules import KineticsRules
@@ -1741,65 +1741,29 @@ def add_reverse_attribute(self, rxn, react_non_reactive=True):
                 rxn.reverse = reactions[0]
                 return True
 
-    def calculate_degeneracy(self, reaction):
+    def calculate_degeneracy(self, reaction, resonance=True):
         """
         For a `reaction`  with `Molecule` or `Species` objects given in the direction in which
-        the kinetics are defined, compute the reaction-path degeneracy.
+        the kinetics are defined, compute the reaction-path degeneracy. Can specify whether to consider resonance.
 
         This method by default adjusts for double counting of identical reactants. 
-        This should only be adjusted once per reaction. To not adjust for 
-        identical reactants (since you will be reducing them later in the algorithm), add
-        `ignoreSameReactants= True` to this method.
+        This should only be adjusted once per reaction. 
         """
-        # Check if the reactants are the same
-        # If they refer to the same memory address, then make a deep copy so
-        # they can be manipulated independently
         reactants = reaction.reactants
-        same_reactants = 0
-        if len(reactants) == 2:
-            if reactants[0] is reactants[1]:
-                reactants[1] = reactants[1].copy(deep=True)
-                same_reactants = 2
-            elif reactants[0].is_isomorphic(reactants[1]):
-                same_reactants = 2
-        elif len(reactants) == 3:
-            same_01 = reactants[0] is reactants[1]
-            same_02 = reactants[0] is reactants[2]
-            if same_01 and same_02:
-                same_reactants = 3
-                reactants[1] = reactants[1].copy(deep=True)
-                reactants[2] = reactants[2].copy(deep=True)
-            elif same_01:
-                same_reactants = 2
-                reactants[1] = reactants[1].copy(deep=True)
-            elif same_02:
-                same_reactants = 2
-                reactants[2] = reactants[2].copy(deep=True)
-            elif reactants[1] is reactants[2]:
-                same_reactants = 2
-                reactants[2] = reactants[2].copy(deep=True)
-            else:
-                same_01 = reactants[0].is_isomorphic(reactants[1])
-                same_02 = reactants[0].is_isomorphic(reactants[2])
-                if same_01 and same_02:
-                    same_reactants = 3
-                elif same_01 or same_02:
-                    same_reactants = 2
-                elif reactants[1].is_isomorphic(reactants[2]):
-                    same_reactants = 2
+        reactants, same_reactants = check_for_same_reactants(reactants)
 
         # Label reactant atoms for proper degeneracy calculation
-        ensure_independent_atom_ids(reactants, resonance=True)
+        ensure_independent_atom_ids(reactants, resonance=resonance)
         molecule_combos = generate_molecule_combos(reactants)
 
         reactions = []
         for combo in molecule_combos:
             reactions.extend(self._generate_reactions(combo, products=reaction.products, forward=True,
-                                                      react_non_reactive=True))
+                                                      prod_resonance=resonance, react_non_reactive=True))
 
         # remove degenerate reactions
         reactions = find_degenerate_reactions(reactions, same_reactants, template=reaction.template,
-                                              kinetics_family=self)
+                                              kinetics_family=self, resonance=resonance)
 
         # log issues
         if len(reactions) != 1: