ReactionMechanismGenerator / RMG-Py

Python version of the amazing Reaction Mechanism Generator (RMG).
http://reactionmechanismgenerator.github.io/RMG-Py/
Other
388 stars 227 forks source link

reverse reactions to CH2(S) #159

Closed bbuesser closed 1 year ago

bbuesser commented 10 years ago

After extending the H_abstraction family I'm running into the following problem:

For the reaction NO2 + CH2(S) -> HONO + CH RMG-Py finds the reverse reaction HONO + CH -> NO2 + CH2(T)

RMG stops because it does not find a valid reverse reaction as the products of the reverse are not equal to the reactants of the forward reaction (singlet-triplet difference)

I would propose to add code that checks if one-centered biradicals are involved and then guides the products of the reverse reaction to the correct state (singlet or triplet). I would make the changes in __generateProductStructures of family.py.

What is your opinion?

rwest commented 10 years ago

Hi Beat, How did you resolve this?

keceli commented 10 years ago

I am not sure if it is related but I have the following error while running the example in scoop directory.

AssertionError: Expecting one matching reverse reaction, not 0. Forward reaction <Molecule "[C]#C"> + <Molecule "O=[C]O"> <=> <Molecule "C#CO"> + <Molecule "[C]=O"> : Reaction(reactants=[Molecule(SMILES="[C]#C"), Molecule(SMILES="O=[C]O")], products=[Molecule(SMILES="C#CO"), Molecule(SMILES="[C]=O")], pairs=[[Molecule(SMILES="[C]#C"), Molecule(SMILES="[C]=O")], [Molecule(SMILES="O=[C]O"), Molecule(SMILES="C#CO")]])

Interestingly, if I draw these molecules based on SMILES, hydrogen is missing on the products side.

I was not getting this error with RMG-database at commit bb84e1e61f6f3fefc9294714c9b26dd575c14b19. So some change in the database might be the reason.

On Fri, Nov 15, 2013 at 2:02 PM, Richard West notifications@github.comwrote:

Hi Beat, How did you resolve this?

— Reply to this email directly or view it on GitHubhttps://github.com/GreenGroup/RMG-Py/issues/159#issuecomment-28594516 .

bbuesser commented 10 years ago

I'm sorry, I've closed that issue the wrong way. I wanted to cite the preliminary solution in commit bbuesser/RMG-Py@f0d2f235458678e317a3c69b97a02792eb5fa214

It is a quick fix that checks if CH2(S) is a reactant in the forward reaction and then checks if the product of the reverse reaction includes a CH2(T). If that's the case it changes the CH2(T) to CH2(S).

bbuesser commented 10 years ago

I just realized that it could that it will lead to the wrong results for the special case of CH2(S) + CH2(T). Propaply I have to break the second for loop once it has found a CH2(T).

Or we need a more general solution if Murat's problem is related.

shamelmerchant commented 10 years ago

I get a very similar error

AssertionError: Expecting one matching reverse reaction, not 0. Forward reaction

Molecule "[O][O]" + Molecule "HN" <=> Molecule "[O]O" + Molecule "N"

Reaction( reactants=[Molecule(SMILES="[O][O]"), Molecule(SMILES="HN")], products=[Molecule(SMILES="[O]O"), Molecule(SMILES="N")], degeneracy=4, pairs=[[Molecule(SMILES="[O][O]"), Molecule(SMILES="[O]O")], [Molecule(SMILES="HN"), Molecule(SMILES="N")]])

bbuesser commented 10 years ago

Shamel's problem on the H-abstraction from NH(singlet) is the same as discussed above for CH2 (singlet) and calls for a more general solution of the problem for reverse reactions of singlet biradicals (coming soon)

bbuesser commented 10 years ago

I have included NH(S) into the quick fix of bbuesser/RMG-Py@f0d2f23 for the CH2(S) case. We expect the same problem as discussed above for NH(S) reacting with NH(T).

bbuesser commented 10 years ago

I think I have identified another, related problem by introducing 1-centered biradicals to the H-abstraction family. RMG finds of course also the special case of CH2(T) + CH = CH + CH2(T). However for the reverse reaction it recognizes that nothing happens and removes the reverse reaction from the reaction list which leads to an error in the methylformate example. I propose as a quick solution bbuesser/RMG-Py@2528de23. Probably we will soon find the same problem for N + NH = NH + N and implement a more general solution.

rwest commented 10 years ago

See also #153. Seem to be a lot of singlet/triplet problems. Can we solve them all at once?

bbuesser commented 10 years ago

Finding all possible electronic states in product structures

The current assumptions for newly found products is that they attain the highest possible multiplicity as this in most cases the most stable species and resulting often in the fastest kinetics. However if we start with a singlet biradical reactant in the forward direction we find with this assumption only the triplet biradical as product of the reverse reaction. I became aware of this because I have added the nitrogen atom as a triradical to the three important reaction families H_abstraction, Disproportionation and R_Addition_MultipleBond.

I think making it very general is the simplest solution. For that I have added the electronic states of quartet (Q) and quintet (V) a575e4fe1d3c9b981cc28a95078dd64181af0269 and a function to set the multiplicity 5ea2e403299281d4c6b22cb10597d0a482df849a. I have also added all electronic states to the reaction families dd492978410bbac7294a08f6abd28e8840140130. I didn’t want to continue hardcode specific exceptions so I have reverted the changes of 2528de2388cd99a95e1e7eb154b06cfe2f6bbed5. I have modified the function __generateProductStructures() to find all possible electronic states for the newly found product structures and add all possible combinations of these products to a product structures list 59e0a7bebecddf0bf328e16722953ee7ead0e8d5. This list can also contain spin-forbidden reactions! But instead of hardcoding the removal of them I think we should add their often very slow kinetics to the library.

I have extended __generateReactions() 59e0a7bebecddf0bf328e16722953ee7ead0e8d5 to make a list of reactions using the new list of product structures. That way it will find the reverse reaction leading to the correct electronic state of the reactants of the forward reaction.

Because we don’t have thermo or kinetic library values for CH(D), C(T), and C(S) I have added them so far to the list of forbidden structures d908a75ce30e8fa95f420a1649902d44f647d8bb and cd9ac0da948fbdf8386e6d8bf26b4a4addce289e.

If you all agree with this solution I will merge it into my master which then can be merged into the official master. I think it solves this #159 and #153 apart from the augmented' InChI's part.

I have added these changes to my spin_conservation branch. @shamelmerchant could you also make a test run with the spin_conserving branch bbuesser/RMG-Py@59e0a7bebecddf0bf328e16722953ee7ead0e8d5 and bbuesser/RMG-databse@fb757925c8875149e9024650868e0e8e5bac8243 to see if this gives similarly reasonable results as my master?

rwest commented 10 years ago

I like the idea of a general solution, but might this take a lot of time generating lots of "possible" but unlikely, unstable, and possibly spin-forbidden structures, for which we have no way to estimate thermo or kinetics, thus making everything worse? Will these structures count as different species? Will we seek kinetics for each of them?

bbuesser commented 10 years ago

I have just this week been working on this part 819b28a7c0e86e17d63edb9fba971af1d8c59b9a and wanted to merge these improvements in my recent pull request.

I did not notice that RMG was slowed down by this. Up to three radicals there are only 2 possible multiplicities in total.

Yes, these new structures are new species because they have different thermo and kinetics.

I think RMG should find all reactions even the spin-forbidden ones and we should control if we want spin forbidden reactions at all or with their slow rates through the reaction families/libraries. However I have recently noticed that there are reaction families that use the same rate for spin-conserving and the way slower spin-forbidden versions (e.g. 1+2_CycloAddition). Bill asked me to fix that.

However RMG handles multiplicity in its own "special" way by treating multiplicity as an atom property instead of as a species property. That can be quite confusing and I think Aron has an issue where the formation of 2-centered biradicals looses the multiplicity information of the combined reactants by storing multiplicity as an atom property and assuming that the triplet state is more stable even if the reaction can only be fast to a less stable singlet. We have encountered a similar problem running the molecule described in #179 and therefore I have limited the general solution to the (labeled) atoms participating in the reaction 819b28a7c0e86e17d63edb9fba971af1d8c59b9a.

I think we should treat multiplicity as a species properties.

rwest commented 10 years ago

Great! There is much for us to learn. I propose an RMG Study Group session on electronic structures, radicals, free electrons, multiplicity, spin, etc. Would you be able to lead it? Whether you do Nitrogen at the same time (I'm guessing it's highly related) is up to you.

bbuesser commented 10 years ago

Yes, I think I would be able to present an overview of current topics in that direction that would be interesting to discuss with the group. Let me know when there is a free slot in the RMG Study Group schedule you would like to fill.

github-actions[bot] commented 1 year ago

This issue is being automatically marked as stale because it has not received any interaction in the last 90 days. Please leave a comment if this is still a relevant issue, otherwise it will automatically be closed in 30 days.