Posts Tagged ‘Theoretical chemistry’

An Ambimodal Trispericyclic Transition State: the effect of solvation?

Thursday, May 2nd, 2019

Ken Houk’s group has recently published this study of cycloaddition reactions, using a combination of classical transition state location followed by molecular dynamics trajectory calculations,[1] and to which Steve Bachrach’s blog alerted me. The reaction struck me as being quite polar (with cyano groups) and so I took a look at the article to see what both the original[2] experimental conditions were and how the new simulations compared. The reaction itself is shown below.


Turns out that chloroform was used as solvent (also benzene), whilst the transition state calculations and the subsequent molecular dynamics trajectories were modelled for the gas phase. The key observation is that if TS1 is used as the starting point for trajectory calculations, only 87% lead to the product predicted by classical transition state theory (3), as revealed by a classical intrinsic reaction coordinate (IRC) calculation. The remaining 13% lead to 4 and 5, the ambimodal effect. So here, I want to explore what effect including a continuum solvent on the computation of TS1 and its IRC might have on the classical (non-dynamic) model.

Firstly, the model for TS1 as reported[3], ωB97XD/6-31G(d) (FAIR data at DOI: 10.14469/hpc/5590). The basis set is modest by today’s standards, but is largely imposed by the need to use it for very large numbers of trajectory calculations. I was able to copy/paste the coordinates from the reported supporting information and then to replicate the IRC at this level (gas phase), also shown in the SI.

There is a feature in this IRC I want to expand upon (red arrow above). It occurs at an IRC value of ~-0.9 on the axis below.

It can be seen more clearly if the RMS gradient norm is plotted, the value of which drops to almost zero at IRC -0.9. Had it reached exactly 0.0, we would have had an intermediate formed. As it is we have what is called a hidden intermediate. The origins of this intermediate can be more readily inferred from this dipole moment plot along the IRC. At TS1, the dipole moment is > 10D. A rule of thumb I have often used is that if a TS has a DM > 10, then one cannot ignore solvation any more and a gas phase model must be augmented.

Here are the same plots, but now with an added solvent field for water, an extreme polarity.

The IRC now stops at -2, being a high energy true (ionic) intermediate (rather than a hidden one). The IRC stops because gradient norm has now reached 0.0 at IRC -2, again an indication of a true intermediate.

The dipole moment has been increased from 10 in the gas phase to around 15 at TS1 and it continues to increase until the ionic intermediate is reached. These differences from the gas phase plots are induced entirely by applying a continuum solvent model.

Next an intermediate solvent, chloroform, being one of the solvents used for the actual reaction. This time the gradient norm almost reaches a value of 0.0, avoiding it only by a whisker! A barely hidden intermediate.

The dipole moment totters around 13.5D, before finally collapsing as the ionic intermediate itself collapses to a neutral molecule again. Benzene as solvent (not shown here) reaches an intermediate dipole moment of about 12D. It too can stabilize an ionic intermediate noticeably even though it is not ionic itself.

I want to also briefly explore what effect if any the use of a relatively small basis set (6-31G(d)) has on the shape of the IRC. Below is a repeat of the gas phase IRC using the Def2-TZVPP basis, which is about twice the size of the smaller one (and hence is around 16 times slower to compute). The gradient norm shows that the “hidden intermediate” region around IRC -1 is a little more prominent (flatter).

So we see that both a solvent model (as a continuum field) and a larger basis set can increase the degree of “hidden intermediate” character in the classical reaction coordinate for this cycloaddition reaction, to the extent that if water is used as model solvent an actual discrete albeit shallow ionic intermediate forms. As Houk puts it, solvation induces a conversion from an entropic intermediate to an enthalpic one.[4] Molecular dynamics trajectories however have a propensity for not settling into quite shallow intermediates (those with escape barriers of < 3 kcal/mol, as would be the case here).

It will indeed be interesting to see the extent, if any, that either of the augmented models shown above affect the calculated distribution of molecular dynamics trajectories compared to those obtained using a gas phase model.

References

  1. X. Xue, C.S. Jamieson, M. Garcia-Borràs, X. Dong, Z. Yang, and K.N. Houk, "Ambimodal Trispericyclic Transition State and Dynamic Control of Periselectivity", Journal of the American Chemical Society, vol. 141, pp. 1217-1221, 2019. https://doi.org/10.1021/jacs.8b12674
  2. C.Y. Liu, and S.T. Ding, "Cycloadditions of electron-deficient 8,8-disubstituted heptafulvenes to electron-rich 6,6-disubstituted fulvenes", The Journal of Organic Chemistry, vol. 57, pp. 4539-4544, 1992. https://doi.org/10.1021/jo00042a039
  3. O.M. Gonzalez-James, E.E. Kwan, and D.A. Singleton, "Entropic Intermediates and Hidden Rate-Limiting Steps in Seemingly Concerted Cycloadditions. Observation, Prediction, and Origin of an Isotope Effect on Recrossing", Journal of the American Chemical Society, vol. 134, pp. 1914-1917, 2012. https://doi.org/10.1021/ja208779k

Tags:, , , , , , , , , , , , ,
Posted in reaction mechanism | 2 Comments »

Reaction coordinates vs Dynamic trajectories as illustrated by an example reaction mechanism.

Monday, March 20th, 2017

The example a few posts back of how methane might invert its configuration by transposing two hydrogen atoms illustrated the reaction mechanism by locating a transition state and following it down in energy using an intrinsic reaction coordinate (IRC). Here I explore an alternative method based instead on computing a molecular dynamics trajectory (MD).

I have used ethane instead of methane, since it is now possible to envisage two outcomes:

An animation of the IRC starting from the located transition state is shown below (DOI: 10.14469/hpc/2331). This is based purely on the computed potential energy surface of the molecule. The IRC is computed from the forces experienced on the atoms as they are displaced from an initial set of coordinates corresponding to the located transition state and then following the direction indicated by the eigenvectors of the negative force constant required of a transition state. Importantly, there is no time component; the path is based entirely on energies and forces.

Next, a molecular dynamics simulation (ωB97XD/6-31G(d,p), DOI: 10.14469/hpc/2330).  This uses the ADMP method, which requests a classical trajectory calculation using the “atom-centered density matrix propagation molecular dynamics model”. This integrates kinetic energy contributions from the molecular vibrations and so explicitly now includes a time component. In this example the evolution of the system from the transition state is charted over a period of 100 femtoseconds (1000 integrated steps). As it happens this is a relatively short period of evolution; sometimes periods of picoseconds may be required to get a realistic model, especially if one is also dealing with explicit solvent molecules (of which perhaps 500 might be required).

Such explicit inclusion of the kinetic energy from molecular vibrations in effect allows the molecule to “jump” over shallow barriers. In this case, the barrier for a [1,2] hydrogen shift from the methyl group to the adjacent carbene (watch atom 8). Simultaneously, the path taken by two hydrogens no longer corresponds to their transposition but to their elimination as a hydrogen molecule. So this quite different outcome from the IRC is very probably also a much more realistic one.

If the MD method is so much more realistic than the IRC, then why is it not always used? The simple answer is computational time! For this very small molecule and using quite a modest basis set (6-31G(d,p)), the relatively short 1000 time steps took about three times as long to compute as the IRC. The factor gets worse as the size of the molecule increases and the number of time steps for a realistic result increases. Perhaps, as technology gets better and new computer architectures emerge, MD simulations of ever increasingly complex reactions will become far more common. In ten years time, I expect most of the examples on this blog will use this method!

Tags:, , , , , , , , , , , , , , , ,
Posted in reaction mechanism | 2 Comments »

VSEPR Theory: A closer look at trifluorothionitrile, NSF3.

Saturday, January 16th, 2016

The post on applying VSEPR ("valence shell electron pair repulsion") theory to the geometry of ClF3 has proved perennially popular. So here is a follow-up on another little molecue, F3SN. As the name implies, it is often represented with an S≡N bond. Here I take a look at the conventional analysis.

trifluoorothionitrile

This is as follows:

  1. Six valence electrons on the central S atom.
  2. Three F atoms contribute one electron each.
  3. One electron from the N σ-bond.
  4. Donate two electrons from S to the two π-bonds.
  5. Eight electrons left around central S, ≡ four valence shell electron pairs.
  6. Hence a tetrahedral geometry.
  7. The bond-bond repulsions however are not all equal. The SN bond repels the three SF bonds more than the S-F bonds repel each-other.
  8. Hence the N-S-F angle is greater than the F-S-F angle, a distorted tetrahedron.

Now for a calculation[1];  ωB97XD/Def2-TZVP, where the wavefunction is analysed using ELF (electron localisation function), which is a useful way of locating the centroids of bonds and lone pairs (click on diagram below to see 3D model).

Trifluorosulfonitrile

  • At the outset one notes that there are six ELF disynaptic basins surrounding the central S, integrating to a total of 7.05e. The sulfur is NOT hypervalent; it does not exceed the octet rule.
  • These six "electron sub-pair" basins are arranged octahedrally around the sulfur. The coordination is NOT tetrahedral, as implied above.
  • The three S-N basins have slightly more electrons (1.25e) than the three S-F basins (1.10e), resulting in …
  • the angle subtended at the S for the SN basins being 96° (a bit larger than octahedral) whilst the angle subtended at the S for the SF basins being smaller (89.9°). This matches point 7 above, but is achieved in an entirely different manner.
  • As a result, the N-S-F angle (122.5°) is larger than the ideal tetrahedral angle and the F-S-F angle (93.9°) is smaller, an alternative way of expressing point 7 above.
  • The S≡N triple bond as shown above does have some reality;  it is a "banana bond" with three connectors rather than two. Each banana bond however has only 1.25e, so the bond order of this motif is ~four (not six) but nevertheless resulting in a short S-N distance (1.406Å) with multiple character.

So we have achieved the same result as classical VSEPR, but using partial rather than full electron pairs to do so. We got the same result with ClF3 before. So perhaps this variation could be called "valence shell partial electron pair repulsions" or VSPEPR.

References

  1. H.S. Rzepa, "F 3 N 1 S 1", 2016. https://doi.org/10.14469/ch/191808

Tags:, , , , , , , , , , ,
Posted in Hypervalency | 110 Comments »