Interpret multiplicity from reaciton smiles

[hwpang]

Interpret multiplicity from reaction smiles for consistency. Open for discussion on how to handle when reactant multiplicity != product multiplicity (this occurs for carbene/nitrene)

[xiaoruiDong]

@hwpang , @kspieks Thanks for bringing it up. In my opinion, there needs to be

• an argument for a user to assign multiplicity. this can have a default value, like None. If value is assigned, then check, probably saturate the molecule, and use user's value; otherwise, infer the multiplicity from the smiles.
• an argument for a user to choose if they want to keep either the larger or the smaller multiplicity whenever the infered multiplicity is equal or larger than 3. It is possible there will be different (TSs) depending on the multiplicity.
• Maybe we need to think about what cases results in this inconsistent problem. Usually it happens for disproportionation reactions, e.g, C=C + C=C <=> C[CH2] + [CH]=C, or singlet carbene reactions C1CC1 = [CH2] + C=C ; the left side is 1 and the right side is 3. I guess we will use 1 for this case. Maybe we can always keep the smaller one; but definitely more examples can be helpful to consider this problem more thoroughly.
• More generally, I think user should store and provide the multiplicity info along with reaction smiles when they generate the list of reactions to compute. Infer from smiles is helpful, but since smiles doesn't have the multiplicity info built-in, it is unreliable.

[hwpang]

@xiaoruiDong Thanks for the comments! I will work on this soon!

[alongd]

Hi, if I may comment, I think that generally the multiplicity is conserved across a reactions (it's called an "adiabatic" reaction), other than special cases where we have inter-surface crossing.

For example, the multiplicity of C=C + C=C <=> C[CH2] + [CH]=C is 1 for the reactants, and 1 for the products as well: the two products have opposite radical spins (they are entangled), so overall it is 1.

For C1CC1 = [CH2] + C=C, I think that the [CH2] species is a singlet carbene, so it's multiplicity is also 1.

There are definitely issues with determining the correct multiplicity for reactions from the SMILES, but I don't think that the "reactant multiplicity != product multiplicity" case is important for "conventional" (adiabatic) reactions (only relevant for inter-surface crossing reactions).

[xiaoruiDong]

@alongd Thanks for the comment. I should have made myself clearer. The reactions we are considering are multiplicity-reserved reactions. The major concern here is that if the user doesn't supply the multiplicity information but only provides SMILES, what are we gonna do. Hao-Wei and Kevin are suggesting parsing from SMILES, so the potential issue will be

1. the molecule parsed from SMILES using RDkit may have undesired (at least for some cases) multiplicities. E.g., [CH2], RDKit always makes a molecule with a multiplicity of 3. More generally, it always doesn't try to match possible electron pairs. E.g., when parsing a SMILES complex C[CH2].[CH]=C, the complex generated will have a multiplicity of 3.
2. Due to 1, we can easily run into cases where generated reactant complex and product complex have different multiplicity. E.g., C=C + C=C <=> C[CH2] + [CH]=C. The reactant C=C.C=C will resulted in a complex with a multiplicity of 1, while the product generated from RDKit will have 3. Since they are different, it is not determinisitic about what multiplicity to use when launching TS jobs. The previous proposal is using reactant's multiplicity, however, the user may be unaware of what multiplicity he is actually using, and the number used may not be the one he actually wants to use. E.g, [O][O] + C[CH2] = [O]O + C=C, if only use the reactant multiplicity it will be 4 (generated from RDKit), but he may actually want to use 2.
3. So I propose. It is always best if the user knows what multiplicity to use. But if he doesn't, we can have an option to decide if the larger or the smaller one is used whenever there is a multiplicity-inconsistency. There can also be cases like [CH3] + [O] = [CH2] + [OH], both sides will be parsed as multiplicity 4 molecule complexes, but the actual one the user want to use is 2, we need to decide if we want to fix it. So, in my opinion, we need to detailedly solve the potential issues listed in 3.

[alongd]

Hi, thanks for the detailed explanation, and sorry for barging in... When users define species, can they specify the multiplicity, or can they only give SMILES/InChI? You are right that cases like (3) are problematic. For the [CH2] + [OH] = [CH3] + [O] reaction we get different products for the different multiplicities (singlet or triplet [O]). Since both reactions are feasible, there's no real way to know what was originally meant from just SMILES. So perhaps there should be a way to specify in advance the multiplicity of species, especially ones like [O] or [CH2]?

[hwpang]

@xiaoruiDong @alongd Thanks for the valuable discussion! We have changed the code as suggested by Xiaorui, where the user can specify the multiplicity of the reaction. If not specified, the multiplicity will be interpreted from the reaction smiles. If the interpreted multiplicity is inconsistent between the reactants and the products, it will choose the multiplicity based on the option use_smaller_multiplicity.