reverse reactions to CH2(S)

[bbuesser]

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]

Hi Beat,
How did you resolve this?

[keceli]

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 bb84e1e. 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//issues/159#issuecomment-28594516
.

[bbuesser]

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

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]

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]

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]

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]

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]

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]

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

[bbuesser]

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) a575e4f and a function to set the multiplicity 5ea2e40. I have also added all electronic states to the reaction families dd49297. I didn’t want to continue hardcode specific exceptions so I have reverted the changes of 2528de2. 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 59e0a7b. 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() 59e0a7b 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 d908a75 and cd9ac0d.

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@59e0a7b and bbuesser/RMG-databse@fb757925c8875149e9024650868e0e8e5bac8243 to see if this gives similarly reasonable results as my master?

[rwest]

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?