Polarizable Continuum Model (PCM) in G03
This is my first post on a series I have in mind regarding frequent questions on the CCL regarding the use of some computational chemistry software, mostly Gaussian. Readers are still encouraged to contact the Gaussian Help Desk for further (and more accurate) help.
Gaussian 03, the popular electronic structure calculation suite of programs, includes the necessary modules for performing calculations in a solvated environment using the continuum models approximations. Among such models, the Polarizable Continuum Model (PCM) is one of the most widely used methods since it meets a good compromise between accuracy and computation time. Nevertheless, Gaussian may not be the best option for performing such calculations (as opposed to other programs as COSMO) but it still can be very useful when used properly. Unfortunately there is a lack of specific info in the literature regarding the usage of the different variables involved in the cavity generation for G03; the newest version, G09, includes some improvements on the corresponding codes making PCM calculations more achievable. While browsing the CCL archives, it is common to find more questions than answers and usually the same questions are posted over and over by different users over time. This post will get updated as needed.
I hope with this post I can summarize most of the common problems found in Gaussian regarding implicit solvation calculations as well as their respective solutions. Some of the solutions come from Gaussian technical support itself, so my best advise is always to address your questions directly to them. Keywords are typed in capital letters, variables in italics.
Implicit solvent calculations imply the generation of a vacuum cavity inside a continuous and homogeneous dielectric field. The simplest model to do this is Onsager’s in which the molecule is treated as a dipole inside a spherical cavity (SCRF=DIPOLE in Gaussian use along the VOLUME keyword to generate the optimum radius for such cavity.) PCM calculations generate a cavity that relates more closely to the molecule’s shape by placing spheres on each atom or groups of atoms. Check the Gaussian link at the bottom of this post for further info; this is a troubleshooting post, not a tutorial on PCM.
Some common errors and their solutions
In order to get a better definition of a cavity it has been recommended to use the option SCRF=(READ,model,SOLVENT=solvent) with the following parameters to be read at the end of the input file:
Additionally we may include a third line indicating the kind of radii to be used on each atom to generate a sphere around it, the default option is Radii=UA0 (Topological United Atoms model) which treats functional groups as a single sphere. Including this line with Radii=UFF; Radii=Pauling or Radii=Bondi will treat each atom independently, which is very useful to use when some H atoms lye outside the UAKS sphere. The error message associated with this problem is: “Error message, treat H atom explicitly” see below
-> BldSpC: Error generating genealogic tree for sphere 309 at level 15
According to Gaussian’s Help Desk, this is a numerical error in the generation of the cavity. The use of fewer spheres (implicit H atoms for instance) is recomended, so if you are using RADII=PAULING or BONDI, delete that line. It is also recommended to use the NOSYMM keyword on the route section. This problem seems to have been addressed in G09.
-> Too many tesserae. Increase the MxTs.
Try using the TSNUM keyword in the route section as SCRF=(TSNUM=num,…) This will modify the number of tiles to describe each sphere that makes the cavity.
-> AdVTs1: ISph= 500 is engulfed by JSph= 501 but Ae( 500) is not yet zero! Error in link301
Generation of cavity fails. Try using a different radii model (RADII=…) and/or the NOSYMMCAV keyword at the end of the file, via the SCRF=(Read,…) option. Also using the NOSYMM keyword in the route section can work. Once again using the OFAC=0.8 and RMIN=0.5 parameters is useful.
-> UA0: Hydrogen 40 is unbound. Keep it explicit at all point on the …
-> UA0: potential energy surface to get meaningful results.
The location of a certain H atom (number 40 in this case) lies outside the cavity placed on a functional group, so it must be treated explicitly by either changing the RADII= model or by placing a sphere on that particular atom alone through the SPHEREONH=40 (40 for this example) option via SCRF=(Read,…)
Additional remarks and suggestions
- The use of spheres on functional groups is suggested for calculating energies, but for geometry optimization the use of a more sophisticated model in generating the cavity is encouraged.
- Always pay attention to the value of the density lying outside the cavity, i.e. inside the dielectric. In G03 this value is labeled as “error on total polarization charges”. As a rule of thumb this value should be less than 0.05 for the calculation to be acceptable.
- Just in case you are using a very old version of Gaussian, be aware that the keyword COSMORS doesn’t launch a COSMO-RS calculation (thermodynamics of solutes and solvents) but a CPCM calculation in a format that can be post-processed by COSMO software.
- PCM calculations are highly parametrized so it’s useful to always have an experimental reference to which you can validate your choices in each calculation.
If you found interesting or helpful information in this post, please leave a comment however short. This will encourage me to keep gathering and posting this kind of information which in turn may be of help for other users, thanks.