Friday, August 14, 2009

Finding a transition state: amide hydrolysis 2

In a previous post I found two TSs, called 3TSa and 3TSb, at the PM3 level of theory.

In this post I show how the PM3 structure and frequency information can be used as a starting point for finding 3TSa at the RHF/3-21G level of theory.

I use RHF/3-21G as an example of an expensive ab initio method where we want to keep the number of frequency calculation to a minimum, since it takes a long time to compute.

So I extract the PM3 optimized geometry from the output file of the TS search, and add to this the corresponding frequency information in the form of the $HESS group from the .dat file.

The keywords I paste in look like this (many of these keywords have been described in a previous post):

$contrl runtyp=sadpoint $end
$scf dirscf=.t. $end
$statpt opttol=0.0005 nstep=50 $end
$statpt hess=read $end
$contrl nzvar=1 $end
$zmat dlc=.t. auto=.t. $end

Note the $basis group is missing, but that I specify it using Avogadro. It's important to pick the basis set before pasting the remaining keywords in because Avogadro deletes some of the lines in the file when you select something in the menu.

Based on the PM3 guess, GAMESS finds the RHF/3-21G TS in 29 steps.

Actually, I can't be 100% sure the structure is a TS because I haven't verified that the structure has an imaginary frequency. Unlike with PM3 I don't automatically compute the frequencies at the end of the TS search (using hssend=.t.).

The reason is that it is a good idea to verify that the structure looks like a TS before committing CPU time to an expensive frequency calculations. As I showed in a previous post it can happen that TS searches end up find the reactants and products.

Luckily the TS structure has one, and only one imaginary frequency as you can see in this screencast

As you saw the default amount of memory that GAMESS requests (1,000,000 words) is insufficient for the RHF/3-21G frequency calculation, but it tells you how much it needs (~1,100,000 words). A "word" is 8 bytes, so we are talking about 10 MB here. I will discuss this issue in more details in a future post.

Finally, in principle it is possible to adjust the memory keyword in Avogadro (in the Advanced setup), but that option is currently not working.

Now we need to repeat this procedure for 3TSb. The procedure are exactly the same, so I haven't made a screencast. The TS search takes 36 steps.


NUchem said...


Thanks for this tutorial. I was wondering how to restart a job using higher level of theory. Your tutorials are more helpful than you think. Thanks again.


Jan Jensen said...

NUchem, you're very welcome. Glad to hear it's useful

pranaldi said...

Hi Jan,

Thanks for your tutorials. I'm working on molecular chemistry. I have a problem in finding transition state of my reaction. Here, Ph-C(=NH)O-NHCH3 + H-(+)N(Et) --> Ph-C(=NH2(+))O-NHCH3 + N(Et). Do you have any idea to find its transition state?

Thank you,

Jan Jensen said...

First, are you sure there is a TS. Many charge transfer rxns in the gas phase occur without a barrier.

To check, perform optimizations of both dimers, where the proton in question is involved in the N-H..N hydrogen bond. My guess is that both optimization will yield the same structure, i.e. one of the dimers is not a minimum on the PES.

pranaldi said...

I've calculated its transition state, and there was negative frequency. But its vibration doesn't match with my expectation.
Indeed, there's hydrogen bond as you guessed. FYI, the reaction I proposed isn't in gas phase, there's solvent (THF) involvement and acid catalyzed. Do you have any link to the reference about any kind of this reaction?

Here is my reaction :