I was reading this paper on the mechanism of a set of organometallic hydrosilylations when I realized that I know next to nothing about how to determine reaction pathways and transition states using Gaussian (or any other computational chemistry software, for that matter). It seems like a pretty important thing for any self-respecting p-chemist to know, so I thought I’d chronicle my findings a little.
It’s possible to glean a few mathematical conditions for transition states just from their definition: because they’re defined as peaks on reaction coordinate diagrams, their first derivatives of energy with respect to all normal modes of vibration are zero. Their second derivatives of energy must be negative W.R.T. all normal modes of vibration except one, in order to establish their energy as a maximum. Difficulties can arise because near a maximum the energy curve tends to flatten out, leading to states with very different geometries but similar energy.
Gaussian approximates the peak of an energy curve using a quadratic function via the STQN (synchronous transit-guided quasi-Newton) method. This process is basically the opposite of optimizations for conventional molecules, which rely on the variational principle and energy minimization. A so-called “synchronous transit” function guides the reaction coordinate and molecular energy close to the transition state, then an optimization (i.e. maximization) algorithm takes over. In fact, if you specify Opt=QST2, Gaussian doesn’t even need an initial guess of the transition state structure to approximate it–the reactant and product structures can serve as input. In this case the guess is just “halfway between” the reactants and products.
Once you’ve determined the transition state, it’s possible to follow the reaction pathway using the IRC (intrinsic reaction coordinate) option. IRC=CalcFC tells Gaussian to calculate force constants of the system, which describe the extent of the forces at work in the molecule, for the first point and use these for all points in the reaction path. You can also grab force constants from a frequency calculation of the transition state using IRC=RCFC. Based on this data Gaussian steps through the reaction path (step size and number of steps can be specified) and optimizes the geometry at each point. Using Results->IRC in GaussView, you can even watch a video of the reaction proceeding! I might try this sometime and post the results.