Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

Sign up for GitHub

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Jump to bottom

Adapting kinetics for with multiple transition states #323

Closed
wants to merge 5 commits into from
Closed

Adapting kinetics for with multiple transition states #323

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
5 commits
Select commit Hold shift + click to select a range
9ea67ff
Find multiple transition states for reactions that match >1 reaction …
connie Aug 16, 2013
61017c7
Increase degeneracy if reaction template is the same but out of order.
connie Aug 20, 2013
ad7b8fd
Do not add reverse reaction if it already exists.
connie Aug 20, 2013
5d1afe8
Detect degeneracy using the reaction's template, not reaction pairs. …
connie Feb 3, 2015
741b549
Use one of the matching reverse reactions. It is currently difficult …
connie Feb 5, 2015
File filter

Filter by extension

Filter by extension
Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
88 changes: 50 additions & 38 deletions rmgpy/data/kinetics/family.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -903,8 +903,8 @@ def addKineticsRulesFromTrainingSet(self, thermoDatabase=None):
item.kinetics = data
data = item.generateReverseRateCoefficient()

item = Reaction(reactants=[m.molecule[0].copy(deep=True) for m in entry.item.products], products=[m.molecule[0].copy(deep=True) for m in entry.item.reactants])
template = self.getReactionTemplate(item)
item = TemplateReaction(reactants=[m.molecule[0].copy(deep=True) for m in entry.item.products], products=[m.molecule[0].copy(deep=True) for m in entry.item.reactants])
item.template = self.getReactionTemplate(item)
item.degeneracy = self.calculateDegeneracy(item)

new_entry = Entry(
Expand Down Expand Up @@ -1221,16 +1221,13 @@ def generateReactions(self, reactants, failsSpeciesConstraints=None):
# for each reaction, make its reverse reaction and store in a 'reverse' attribute
for rxn in reactionList:
reactions = self.__generateReactions(rxn.products, products=rxn.reactants, forward=True, failsSpeciesConstraints=failsSpeciesConstraints)
if len(reactions) != 1:
logging.error("Expecting one matching reverse reaction, not {0} in reaction family {1} for forward reaction {2}.\n".format(len(reactions), self.label, str(rxn)))
for reactant in rxn.reactants:
logging.info("Reactant")
logging.info(reactant.toAdjacencyList())
for product in rxn.products:
logging.info("Product")
logging.info(product.toAdjacencyList())
raise KineticsError("Did not find reverse reaction in reaction family {0} for reaction {1}.".format(self.label, str(rxn)))
rxn.reverse = reactions[0]
#assert len(reactions) == 1, "Expecting one matching reverse reaction, not {0}. Forward reaction {1!s} : {1!r}".format(len(reactions), rxn)
if len(reactions) == 0:
raise Exception("Could not match to reverse reaction due to no matching pairs. Forward reaction {1!s} : {1!r}".format(len(reactions), rxn))

else:
# Use one of the matching reverse reactions
rxn.reverse = reactions[0]

else: # family is not ownReverse
# Reverse direction (the direction in which kinetics is not defined)
Expand All @@ -1256,12 +1253,16 @@ def calculateDegeneracy(self, reaction):
"""
reactions = self.__generateReactions(reaction.reactants, products=reaction.products, forward=True)
if len(reactions) != 1:
for reactant in reaction.reactants:
logging.error(reactant)
for product in reaction.products:
logging.error(product)
if not reaction.pairs:
reaction.pairs = self.getReactionPairs(reaction)
for possibleReaction in reactions:
# Match the correct reaction if there are multiple reaction pathways
if set(possibleReaction.template) == set(reaction.template): return possibleReaction.degeneracy

raise Exception('Unable to calculate degeneracy for reaction {0} in reaction family {1}.'.format(reaction, self.label))
return reactions[0].degeneracy

else:
return reactions[0].degeneracy

def __generateReactions(self, reactants, products=None, forward=True, failsSpeciesConstraints=None):
"""
Expand Down Expand Up @@ -1396,6 +1397,29 @@ def __generateReactions(self, reactants, products=None, forward=True, failsSpeci
if match:
rxnList.append(reaction)



# Determine the reactant-product pairs to use for flux analysis
# Also store the reaction template (useful so we can easily get the kinetics later)
for reaction in rxnList:

# Restore the labeled atoms long enough to generate some metadata
for reactant in reaction.reactants:
reactant.clearLabeledAtoms()
for label, atom in reaction.labeledAtoms:
atom.label = label

# Generate metadata about the reaction that we will need later
reaction.pairs = self.getReactionPairs(reaction)
reaction.template = self.getReactionTemplate(reaction)

# Unlabel the atoms
for label, atom in reaction.labeledAtoms:
atom.label = ''

# We're done with the labeled atoms, so delete the attribute
del reaction.labeledAtoms

# The reaction list may contain duplicates of the same reaction
# These duplicates should be combined (by increasing the degeneracy of
# one of the copies and removing the others)
Expand Down Expand Up @@ -1433,13 +1457,19 @@ def __generateReactions(self, reactants, products=None, forward=True, failsSpeci
# If we found a match, remove it from the list
# Also increment the reaction path degeneracy of the remaining reaction
if match:
rxnList.remove(reaction)
reaction0.degeneracy += 1
if set(reaction0.template) == set(reaction.template):
# template must match, otherwise we may have found an alternate transition state
reaction0.degeneracy += 1
rxnList.remove(reaction)
else:
index += 1
else:
index += 1

index0 += 1



# For R_Recombination reactions, the degeneracy is twice what it should
# be, so divide those by two
# This is hardcoding of reaction families!
Expand All @@ -1450,28 +1480,10 @@ def __generateReactions(self, reactants, products=None, forward=True, failsSpeci
assert(rxn.degeneracy % 2 == 0)
rxn.degeneracy /= 2

# Determine the reactant-product pairs to use for flux analysis
# Also store the reaction template (useful so we can easily get the kinetics later)
for reaction in rxnList:

# Restore the labeled atoms long enough to generate some metadata
for reactant in reaction.reactants:
reactant.clearLabeledAtoms()
for label, atom in reaction.labeledAtoms:
atom.label = label

# Generate metadata about the reaction that we will need later
reaction.pairs = self.getReactionPairs(reaction)
reaction.template = self.getReactionTemplate(reaction)
if not forward:
reaction.degeneracy = self.calculateDegeneracy(reaction)

# Unlabel the atoms
for label, atom in reaction.labeledAtoms:
atom.label = ''

# We're done with the labeled atoms, so delete the attribute
del reaction.labeledAtoms


# This reaction list has only checked for duplicates within itself, not
# with the global list of reactions
Expand Down
6 changes: 4 additions & 2 deletions rmgpy/rmg/model.py
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -483,8 +483,10 @@ def checkForExistingReaction(self, rxn):
if not rxn.duplicate:
return True, rxn0
else:
return True, rxn0

if set(rxn0.template) == set(rxn.template):
return True, rxn0

# Do not add reverse reaction if it already exists.
if isinstance(family,KineticsFamily) and family.ownReverse:
if (rxn0.reactants == rxn.products and rxn0.products == rxn.reactants):
return True, rxn0
Expand Down