Different TSs for symmetry equivalent reactions

[adragthest]

Hello,

I have been running two autodE reaction calculations, which should lead to the same products and have the same TS due to symmetry; a symmetry which I hadn't realised initially. Thus the SMILES strings and the conditions of the reaction are exactly the same between the two runs. I am doing the calculations using:

rxn.calculate_reaction_profile(with_complexes=False)

Interestingly although the reactants and products are the same in both cases (in terms of conformations), the conformations of the TS are different. It seems that in the second run autodE failed to locate what I assume is the global minimum in the conformational space of the TS.

I understand that certain parts of autodE are of stochastic nature, but this didn't seem to affect the results for the reactants and products. Obviously something has been done differently in the TS location/TS conformer generation phase of the calculations.

Do you know why this behaviour might occur and how to troubleshoot it? How are the initial reactant association complexes created? In other words I am trying to understand how the two reactants are brought together to create the initial TS structure; for non symmetric systems there are many ways one of the reactant can approach the other, and all can lead to a viable TS. Is it through NCI? If not, how is it done? This might help me understand the reason for this behaviour.

Thank you!

P.S. I am using the latest version of autodE (1.4.3)

[t-young31]

Hi @adragthest – all the conformer searchers are indeed stochastic so some variance can be expected with reactions with significant conformation freedom.

Do you know why this behaviour might occur and how to troubleshoot it?

The approach to conformer searches for minima and TSs are different for organic systems. The former uses RDKit's ETKDGv3 algorithm while TSs use a less informed randomise+minimise approach. The maximum number of generated conformers is the same (set with Config.num_conformers) but the space searched in TSs is probably bigger. The default maximum number of conformers attempts to balance computational cost with a reasonable search of the space, so increasing it (Config.num_conformers) should reduce the variance between runs.

How are the initial reactant association complexes created?

It depends on the reaction. However, TS conformer searches don't proceed by searching the association complex conformational space.

Hope this helps!

[adragthest]

@t-young31 Thank you very much for your feedback. I think I have probably confused you because I didn't pose my question correctly; although your feedback on the conformers is certainly useful.

What I was trying to figure out is how autodE creates the very first 3D model of the TS after the conformers of the reactants and products have been created and ranked and before any attempt to precisely find the TS has been performed (e.g. via adaptive path). How does autodE move from a 2D graph representation of the initial TS guess to an initial 3D atom arrangement? It might be that there are different "configurations" of reactants and products that might lead to equally valid TSs; I was wondering how autodE creates these, how it assesses their likelihood of occurring and chooses which to use for the downstream workflows (e.g. finding a TS using adaptive path). Put differently, there are many ways the two reactants can come together and these might be equally valid TS guesses for giving the desired products. Is there any way the user might have a control on how these are created and also in which files can one find these alternative "association" complexes that are used for the TS locations?

It seems that although these "association" complexes are not conformers in the normal sense but rather different initial starts for the location of different TSs autodE can indeed locate them and promote them as TS candidates but in a rather stochastic fashion: in different runs different arrangements of the two reactants can occur, but what autodE progresses downstream seems to be stochastic. Is there any way the user can guide autodE to consider all the alternatives (or at least more alternatives) and consistently promote the same TS for the subsequent workflow steps in a reproducible among the runs fashion?

Thank you!

[t-young31]

What I was trying to figure out is how autodE creates the very first 3D model of the TS after the conformers of the reactants and products have been created

So this depends on the reaction type – it's probably better explained in the paper, but briefly alignment is only performed if 1 bond breaks+1 forms (code) by minimising an forcefield that includes the alignment and repulsion, for the lowest energy conformer of all the reactants.

It might be that there are different "configurations" of reactants and products that might lead to equally valid TSs; I was wondering how autodE creates these, how it assesses their likelihood of occurring and chooses which to use for the downstream workflows

Certainly! However, there isn't any enumeration of reactant/product complexes as the space is just too big (e.g n_reactant1_confs x n_reactant2_confs x n_arrangements for an SN2). I'm fairly sure in the paper we show the lack of correlation between minima and TS conformer rankings, so it's not like you can choose the lowest N and proceed from there.

Is there any way the user can guide autodE to consider all the alternatives (or at least more alternatives) and consistently promote the same TS for the subsequent workflow steps in a reproducible among the runs fashion?

The alignment is reasonably non-stochastic from my experience, with a limited sample set of course! You could generate some reactant complexes and call get_ts_adaptive_path on them, but as I said, the space is huge and the rate-governing energy difference is to the lowest energy TS.

[adragthest]

Many thanks for your feedback!

[t-young31]

No worries! Feel free to reach out if you've got more questions