Thursday, December 30, 2010

In my group we regularly use semi-empirical QM methods when we need to minimize large structures. Lately we've had some good results with Jimmy Stewart's PM6 method (he also go by the, IMO, cool name of Mr. MOPAC). Recently, the 25 year old AM1 model of Dewar et al. has been reparametrized, by the group of Mr. MOPAC, into the RM1 method. Supposedly RM1 beat the [euphemism] out of older parameterizations, such as AM1, PM3 and PM5.

In short: Right now the big players in the field of semi-empirical methods seem to be the well-known PM6 and newcomer RM1. The approximations in the two are very alike with the most notable difference between the two, probably being the introduction of d-orbitals on main group elements in PM6, which are not present in AM1 and it's derivative reparametrizations.

I'm not aware of any program which does RM1 in parallel out-of-the-box. If you do know of any open source or free program capable of this, do leave a comment here! I know that Gaussian 09 does PM6 and AM1 and will let you run on multiple cores too. The down side is that I don't have a clue of how much speed up you get from running semi-empirical calculations in parallel in Gaussian09. From my experience, I've seen good scaling with parallel DFT calculations in Gaussian, but I've never really tested the semi-empirical methods in parallel this way. Step one is of course to get RM1 up and running in Gaussian. This proved to be a near impossible task.

Here I present how to make Gaussian09 do RM1 calculations. This is something I've been trying to figure out for about a week now, and solely using the information found on the RM1 website combined with the sub-parGaussian09 manual, I've had little luck, until I contacted Mahmoud A. A. Ibrahim, who kindly supplied me with a g09 input-file. These things should be really simple, but when the documentation is erroneous and incomplete, there is not much to do other than old fashioned trial-and-error and guesswork. However, since I felt that my time was more valuable than this, I more or less gave up until I five minutes later got a reply. So big thanks to Ibrahim for having sorted this out! Below is an outline of how it is done. I hope the Gaussian-staff will include this in the manual, because it is by no means straightforward to see how and where to get the parameters and the correct names and units. For instance, the GCore parameters are not mentioned in the manual, so it can be hard to specify these correctly, and I believe that some of the entries in the MOPAC to Gaussian parameter table are erroneously paired up.

Do remember, that Gaussian is ultra picky about blank lines, to be sure to have exactly one blank line between sections and your file must be terminated with two or more blank lines.

A short explanation is in order. The method is specified as AM1, because AM1 and RM1 use identical formulas. The Input keyword tells AM1 to read the parameters from the input stream (i.e. the stuff in the bottom part of the file). I'm including the Print and NoGenerate keywords as a safety check. Print tells Gaussian to print the used parameters in the output file. This way we're sure, that our optimization used the right keywords. NoGenerate tells Gaussian to not use the native AM1 parameters if a parameter is missing in the input stream. This way we're likely to get some sort of failure in the calculation if our parameters are not read in correctly.

For those who worry about too much clutter in their input files, it is of course safe to leave out the parameters for elements which are not present in the molecule.

EDIT: To honor the fact, that this website promotes free software, let me point out, that the RM1 method is public domain and the parameters are fully disclosed. It is also included in the free quantum chemistry programs GAMESS and MOPAC. GAMESS is free, while MOPAC is free for academic or non-profit purposes.