KineticsError: Unable to calculate degeneracy in reaction family 1,2_Insertion_carbene

[DetlevCM]

This error occurred on a forked build - so I then ran it on (an older) regular RMG-Py install and obtained the same error...

At the end of the mechanism, RMG just fails creating an unexpected biradical as well as phenol and a cresol...

Error: <Molecule "[CH2]">
Error: <Molecule "OC1C=CC=CC=1">
Error: <Molecule "CC1=CC=CC(O)=C1">
Traceback (most recent call last): 
  File "/home/irsrvshare2/R10/L0277/DCM_RMG/VPYTH_uptd/RMG-Py/rmg.py", line 165, in <module>
    rmg.execute(inputFile, output_dir, **kwargs)
  File "/home/irsrvshare2/R10/L0277/DCM_RMG/VPYTH_uptd/RMG-Py/rmgpy/rmg/main.py", line 541, in execute
    bimolecularReact=self.bimolecularReact)
  File "/home/irsrvshare2/R10/L0277/DCM_RMG/VPYTH_uptd/RMG-Py/rmgpy/rmg/model.py", line 748, in enlarge
    self.processNewReactions(self.react(database, self.core.species[i]), self.core.species[i], None)
  File "/home/irsrvshare2/R10/L0277/DCM_RMG/VPYTH_uptd/RMG-Py/rmgpy/rmg/model.py", line 649, in react
    reactionList.extend(database.kinetics.generateReactionsFromFamilies([moleculeA], products=None, only_families=only_families))
  File "/home/irsrvshare2/R10/L0277/DCM_RMG/VPYTH_uptd/RMG-Py/rmgpy/data/kinetics/database.py", line 419, in generateReactionsFromFamilies
    reactionList.extend(family.generateReactions(reactants))
  File "/home/irsrvshare2/R10/L0277/DCM_RMG/VPYTH_uptd/RMG-Py/rmgpy/data/kinetics/family.py", line 1244, in generateReactions
    reactionList.extend(self.__generateReactions(reactants, forward=False))
  File "/home/irsrvshare2/R10/L0277/DCM_RMG/VPYTH_uptd/RMG-Py/rmgpy/data/kinetics/family.py", line 1474, in __generateReactions
    reaction.degeneracy = self.calculateDegeneracy(reaction)
  File "/home/irsrvshare2/R10/L0277/DCM_RMG/VPYTH_uptd/RMG-Py/rmgpy/data/kinetics/family.py", line 1272, in calculateDegeneracy
    'but generated {2}').format(reaction, self.label, len(reactions)))
rmgpy.data.kinetics.common.KineticsError: Unable to calculate degeneracy for reaction <Molecule "[CH2]"> + <Molecule "OC1C=CC=CC=1"> <=> <Molecule "CC1=CC=CC(O)=C1"> in reaction family 1,2_Insertion_carbene. Expected 1 reaction but generated 2

Output from the RMG job including error at the end
out.log.txt

(This could potentially be fixed in a newer build of RMG - and it may also be related to a similarly funny issue in RMG Java where two identical biradicals lead to an error termination.)

[rwest]

I think related to (or perhaps even same as) #617 and #545. It's another "unable to calculate degeneracy" error. Maybe #142 and others too

[DetlevCM]

Comparing the issues you linked it seems like that is spot on - just triggered by a different group.
(Error messages look similar enough to my eyes)

[connie]

These errors are usually related to problems within the kinetic groups where the reaction can get formed in one direction but not the other. You'd probably have to debug while in the rmgpy.data.kinetics module and the specific reaction family that the reaction fails in.

[DetlevCM]

Thanks for the extra information. Any more explicit suggestions as to how I could go about doing that?

What do you think about adding forbidden species for specific groups? (Karl mentioned they can be defined per group.)
Maybe have forbidden species tied to a temperature range?

[connie]

I added this line of code: https://github.com/ReactionMechanismGenerator/RMG-Py/blob/master/rmgpy/data/kinetics/family.py#L1330 which should actually debug the case where a forbidden structure is found in one direction but not the other. So RMG will progress if such a case is found and simply skip the reaction because if it's forbidden in one direction, it should probably be forbidden in the other direction.

Regarding forbidden species, there are 2 types: global and family specific. Globally forbidden structures are in RMG-database/input/forbiddenStructures.py and are forbidden to be formed in all families. So RMG will never generate them unless they are used in a library. The family specific groups are appended to the normal RMG-database/input/kinetics/family/family_name/groups.py file. This will prevent that family from having reactants of the forbidden group type, and usually this is done because we know the kinetics from that family does not apply well in that specific case (big hindrance, or cases like Diels_alder on aromatic benzene). We currently do not have forbidden species tied to a specific temp range. Forbidden structures are actually a poor shortcut to truly improving the database, because the better solution would be to have actual kinetics (likely very slow), built into the database. The prevention of finding these paths can lead to trouble down the road in the case of low temperature or reaction systems where these reactions, though having slow kinetics, do matter.

[DetlevCM]

@connie so I guess "hacking out" the troublemaking species via forbidden species isn't a good idea.

I have turned off the carbene library for now to see where it goes...
As to your code, I'm not 100% sure I follow how it applies to this case - or are you saying that RMG can only calculate the rate in one direction in this case?

@rwest Thank you for changing the title - this is much more informative.

[connie]

In the case of the 1,2_Insertion_carbene family, it is not reversible. So I think it's not a problem of finding rates in both directions.

From looking at your case more specifically, it is because RMG is detecting 2 transition states, producing 2 reactions rather than one. The reaction families recognizes the kekulized benzene double vs single bonds differently so it's treating the cyclic transition state as 2 distinct transition states, even though there should only be one since the bonds in benzene are equivalent. Since we consider the resonance isomers (both kekulized and aromatized) a single species, RMG thinks there should only be one reaction. So therein lies the conflict. There should only be 1 reaction, but the family says there's two and can't collapse them together.

I think the long term fix would have to be to create separate aromatic families, or to create new aromatic groups within the current families. Otherwise you could hack a fix by taking the kinetics of only the first reaction.

[DetlevCM]

@connie Thanks for the detailed reply (coming back to it after the weekend) which means I now at least understand where the error is coming from.

I guess I will need to see which course of action I will take with my colleagues.

[DetlevCM]

OK, so I switched off the carbene library as there is no literature evidence to suggest it would apply to my case in the liquid phase. (not dealing with ionic liquids or the relevant initiators)

And I got an error again - which to me sounds very much like the issue you described. In fact, RMG reports that it expected 1 reaction and found two.

This is the error output (don't think the rest of the log is needed)

 Adding species C=C1[CH]C([O])=CC=C1(49357) to model core
 Generating thermodynamics for new species...
 Generating kinetics for new reactions...

 Summary of Model Enlargement
 ---------------------------------
 Added 1 new core species
      C=C1[CH]C([O])=CC=C1(49357)
 Created 0 new edge species
 Moved 1 reactions from edge to core
      C=C1[CH]C([O])=CC=C1(49357) <=> C=C1C=CC=C2OC12(40915)
 Added 0 new core reactions
 Created 0 new edge reactions

 After model enlargement:
      The model core has 68 species and 183 reactions
      The model edge has 53179 species and 105764 reactions


 Error: Expecting one matching reverse reaction, not 2 in reaction family H_Abstraction for forward reaction <Molecule "Cc1cccc(O)c1"> + <Molecule "[CH2]C1=CC=CC([O])=C1"> <=> <Molecule "CC1=CC=CC([O])=C1"> + <Molecule "[CH2]c1cccc(O)c1">.

 Error: There is likely a bug in the RMG-database kinetics reaction family involving a missing group, missing atomlabels, forbidden groups, etc.
 Reactant
 1  C u0 p0 c0 {2,S} {9,S} {10,S} {11,S}
 2  C u0 p0 c0 {1,S} {4,B} {5,B}
 3  C u0 p0 c0 {4,B} {6,B} {8,S}
 4  C u0 p0 c0 {2,B} {3,B} {12,S}
 5  C u0 p0 c0 {2,B} {7,B} {15,S}
 6  C u0 p0 c0 {3,B} {7,B} {13,S}
 7  C u0 p0 c0 {5,B} {6,B} {14,S}
 8  O u0 p2 c0 {3,S} {16,S}
 9  H u0 p0 c0 {1,S}
 10 H u0 p0 c0 {1,S}
 11 H u0 p0 c0 {1,S}
 12 H u0 p0 c0 {4,S}
 13 H u0 p0 c0 {6,S}
 14 H u0 p0 c0 {7,S}
 15 H u0 p0 c0 {5,S}
 16 H u0 p0 c0 {8,S}

 Reactant
 multiplicity 3
 1  C u0 p0 c0 {2,S} {3,D} {7,S}
 2  C u0 p0 c0 {1,S} {4,D} {9,S}
 3  C u0 p0 c0 {1,D} {6,S} {10,S}
 4  C u0 p0 c0 {2,D} {5,S} {8,S}
 5  C u0 p0 c0 {4,S} {6,D} {12,S}
 6  C u0 p0 c0 {3,S} {5,D} {11,S}
 7  C u1 p0 c0 {1,S} {13,S} {14,S}
 8  O u1 p2 c0 {4,S}
 9  H u0 p0 c0 {2,S}
 10 H u0 p0 c0 {3,S}
 11 H u0 p0 c0 {6,S}
 12 H u0 p0 c0 {5,S}
 13 H u0 p0 c0 {7,S}
 14 H u0 p0 c0 {7,S}

 Product
 multiplicity 2
 1  C u0 p0 c0 {2,S} {8,S} {9,S} {10,S}
 2  C u0 p0 c0 {1,S} {3,S} {4,D}
 3  C u0 p0 c0 {2,S} {5,D} {11,S}
 4  C u0 p0 c0 {2,D} {6,S} {12,S}
 5  C u0 p0 c0 {3,D} {7,S} {13,S}
 6  C u0 p0 c0 {4,S} {7,D} {14,S}
 7  C u0 p0 c0 {5,S} {6,D} {15,S}
 8  H u0 p0 c0 {1,S}
 9  H u0 p0 c0 {1,S}
 10 H u0 p0 c0 {1,S}
 11 H u0 p0 c0 {3,S}
 12 H u0 p0 c0 {4,S}
 13 O u1 p2 c0 {5,S}
 14 H u0 p0 c0 {6,S}
 15 H u0 p0 c0 {7,S}

 Product
 multiplicity 2
 1  C u0 p0 c0 {3,B} {4,B} {7,S}
 2  C u0 p0 c0 {3,B} {5,B} {8,S}
 3  C u0 p0 c0 {1,B} {2,B} {9,S}
 4  C u0 p0 c0 {1,B} {6,B} {12,S}
 5  C u0 p0 c0 {2,B} {6,B} {10,S}
 6  C u0 p0 c0 {4,B} {5,B} {11,S}
 7  C u1 p0 c0 {1,S} {13,S} {14,S}
 8  O u0 p2 c0 {2,S} {15,S}
 9  H u0 p0 c0 {3,S}
 10 H u0 p0 c0 {5,S}
 11 H u0 p0 c0 {6,S}
 12 H u0 p0 c0 {4,S}
 13 H u0 p0 c0 {7,S}
 14 H u0 p0 c0 {7,S}
 15 H u0 p0 c0 {8,S}

 Error: Debugging why no reaction was found...
 Error: Checking whether the family's forbidden species have affected reaction generation...
 Error: Still experiencing error: Expecting one matching reverse reaction, not 2 in reaction family H_Abstraction for forward reaction <Molecule "Cc1cccc(O)c1"> + <Molecule "[CH2]C1=CC=CC([O])=C1"> <=> <Molecule "CC1=CC=CC([O])=C1"> + <Molecule "[CH2]c1cccc(O)c1">.

expected-1-reac-got-2.txt

Incidentally, RMG Java has an issue with the species called Reactant2 as well and fails, however only with a less informative Species = Species message.

[DetlevCM]

Followup query:

Would it be possible to disable a biradical in a reaction family as to avoid the problematic species altogether? Obviously I cannot disable global biradicals as this clashes with oxygen - however that is the only biradical that is reasonable for the liquid phase systems I am looking at.

[connie]

I would not recommend using RMG-Java anymore, since it is no longer developed. Also it has a tendency to silently squash chemistry errors from my experience.

Yes, you can disable the birad in that specific reaction family if you just want to run it again. Simply convert the entire molecule into a group and paste it into that family's groups.py file as an additional forbidden group.

[DetlevCM]

Disabling a single biradical wouldn't help as RMG is very creative at finding a different one. One can at
times anticipate one of them, but then it just finds another...
I.e. coming back to the previous query, it isn't so much the issue of removing one species but the whole class of biradicals except for oxygen.

The other issue that remains is the point you originally raised about two reactions being found where one is expected.
Do you know where in the code the error is generated? How good or bad would it be to just average the predicted values from the two reactions to obtain one set of parameters?

[connie]

Banning birads entirely sounds dangerous.

The error is raised here: https://github.com/ReactionMechanismGenerator/RMG-Py/blob/master/rmgpy/data/kinetics/family.py#L1327

If len(reactions) > 1 you could try just returning the avg or even the first of the list as a hack.

[DetlevCM]

@connie Thanks.

I'm looking at the code and am trying to figure out how RMG-Py is actually written....

Is reactions an array of all found reactions applicable to a combination of species?
And rxn is then an individual entry that contains the reactants and products for one found reaction?

So where are the Arrhenius parameters estimated? Does this take place elsewhere?

[goldmanm]

Yes and yes. Reaction objects also contain kinetics and other specific information, like degeneracy. Arrhenius parameters are estimated in the rmgpy/data/kinetics/family.py file.

We just pulled in code that allows multiple transition states. Feel free to reopen if this doesn't solve the issue