How is B3LYP implemented in Gaussin 0*, GAMESS-US, Molpro, ... etc?

Question

Specifically I want to extend work involving B3LYP started with Gaussian 03 but continued with GAMESS-US. The energies provided by the default B3LYP methods are not the same. There is a discussion about this in the GAMESS-US manual (Further Information section):

Note that B3LYP in GAMESS is based in part on the VWN5 electron gas correlation functional. Since there are five formulae with two possible parameterizations mentioned in the VWN paper about local correlation, other programs may use other choices, and therefore generate different B3LYP energies. For example, NWChem's manual says it uses the "VWN 1 functional with RPA parameters as opposed to the prescribed Monte Carlo parameters" as its default. Should you wish to use this VWN1 formula in a B3LYP hybrid, simply choose "DFTTYP=B3LYP1".

It says the defaults are different between GAMESS and NWCHEM and that there is an option to get GAMESS to do the same calculation as NWCHEM is doing by default.

How do I get a G03 and a GAMESS B3LYP calculation to agree?

What are the differences between various software packages' default implementation of B3LYP and their capabilities, i.e. can their B3LYP definitions/implementations be adjusted?

Answer 1

Aesin already answered part of your question. I can provide you with some more information on GAMESS(US).

It is possible to make GAMESS(US) use the same 'type' of B3LYP as Gaussian 03. For this, you need to specify "DFT=B3LYP1" as you already mentioned in your question. This selects B3LYP with VWN formula 1 RPA local correlation which, to the best of my knowledge, is identical to what they call VWN formula III in some other programs (like Gaussian 03).

Of course, choosing the same functional in both programs is not the only requirement to get identical results. Some other things to consider are:

• Basis set. Both programs have to use the exact same basis set. If you're using an internally stored basis set in Gaussian (like, for example, 6-31G(d,p)), you can make Gaussian print the basis set details by adding the keyword GFINPUT to the route section. GAMESS(US) prints basis set details in its main output.

• Grid size. By default, Gaussian 03 uses a (75,302) grid while GAMESS(US) uses a (96,302) grid. In Gaussian, the grid size can be controlled by the INT keyword. In GAMESS(US), you should take a look at the NRAD and NLEB keywords in the $DFT group. The type of grid will also make a difference, but to the best of my knowledge, GAMESS(US) and Gaussian make use of similar grids.

• Integral cutoff. Both programs neglect very small integrals, as this will speed up the calculation without having a significant impact on accuracy. However, the cutoff factors between the two programs can be different, which can lead to slightly different results. In Gaussian 03, you can control the cutoff factor with IOP(3/27). In GAMESS(US), you can use the ICUT keyword in $CONTRL.

• SCF convergence. Gaussian normally uses EDIIS and CDIIS for the SCF procedure while GAMESS(US) uses DIIS or SOSCF. Both methods should converge towards the same solution, provided that your case is not too complex for DFT. However, you should specify very tight convergence criteria if you'd like to compare energies obtained with both programs.

• With regard to geometry optimizations: Gaussian and GAMESS(US) use very different coordinate systems, geometry optimizers, and convergence criteria. Making both programs optimize to the exact same geometry is difficult, maybe even impossible.

There might be other minor differences you've got to take into account. Maybe it is best to start with a Hartree-Fock calculation, and see if the two programs yield the same SCF energy -- this takes differences in functional and DFT grid out of the equation.

Hope this helps.

Answer 2

The Gaussian implementation of B3LYP uses the VWN3 functional, according to the manual.

Making Gaussian use the VWN5 functional instead for it is a bit tricky, but can apparently be done by adding all the following to the route line:

• bv5lyp - to specify which functional components - Becke exchange, and VWN5 local, LYP non-local correlation.
• iop(3/76=1000002000) - 20% HF exchange, plus
• iop(3/77=0720008000) - 72% Becke non-local exchange, plus 80% Slater local exchange, plus
• iop(3/78=0810010000) - 81% LYP non-local correlation, plus 100% of the V5LYP VWN5 local correlation.

(You can see why people try to avoid using the IOP keyword.) More information on using these is on the aforementioned DFT keywords page of the Gaussian manual, under 'User-Defined Models'.

I'm not that familiar with GAMESS, but it doesn't seem to have the option of using the VWN3 version of B3LYP, so it doesn't seem as though you can go the other way.

As for these and adaptability in other packages, I know Turbomole has listed both B3LYP (using VWN5) and B3LYP_Gaussian (using VWN3), and the manual for ADF suggests you can only use VWN5 for B3LYP there, but you can tweak the amount of HF exchange if that's something you want to do.

Answer 3

NWChem's built-in B3LYP is supposed to agree with Gaussian's, modulo the grid and tolerance issues noted in Thom's answer. You can prescribe any functional form for which the constituents are supported using the explicit XC interface: http://www.nwchem-sw.org/index.php/Density_Functional_Theory_for_Molecules#XC_and_DECOMP_--_Exchange-Correlation_Potentials.

I recognize that the question has already been answered but I wanted to add the NWChem details for completeness since the question indicates a desire to make GAMESS=NWChem=Gaussian.

As a side note, Dalton supports B3LYP and B3LYP-G. The latter agrees with Gaussian while the former is the more canonical version that may agree with GAMESS.