G09
go back to Main Page, Computational Resources, Chemistry & More, Computational Codes, GAUSSIAN
Tricks
Dispersion
One of the new features in Gaussian G09, is the possibility to include dispersion. This can be done in two different ways.
- By using the B97D functional, which includes dispersion. With this you can optimize or do single point calculations.
- With the iop: iop(3/124=3),which allows to compute the Dispersion on a given structure optimized at your favourite level of theory (to add Dispersion as single point). The good point of this is that can be done using any set of basis and method (using SDDall you can save lot's of computer time).
The issue is that the first and second revision of Gaussian A1 and A2 do NOT perform the Dispersion calculation properly. You need the revision B1. This revision is not available in kimik2 and you will need to use CESCA to work on it.
DO NOT RUN GAUSSIAN B97D or iop(3/124=3) with REVISIONS A1 and A2 = KIMIK2
IF YOU RUN THEM ON CESCA CHANGE TO REVISION B1
ONIOM in Gaussian09
For convergence problems (of the optimization) when using ONIOM in gaussian09 you can try the QuadMac option with the opt keyword (or Opt=QuadMacro). This option takes into account the coupling between atoms in the model system and the atoms only in the MM layer in order to produce more accurate steps than the regular microiterations algorithm. Also, sometimes for low frequency transition structures, it might be useful, but sometimes it is not useful.
From gaussian manual:
- For problem convergence cases, the main alternative is Opt=QuadMac, which does a quadratic step in the coordinates of all the atoms. Such an optimization can use either an updated approximate Hessian for the internal coordinates or an analytically computed Hessian. It computes products of the low-level (MM) Hessian with vectors as needed. If there are still convergence problems, then try Opt=(QuadMac,CalcFC) or Opt=(QuadMac,CalcAll).
- For transition structures, the QuadMac option helps to ensure that you move in the direction of a proper transition structure, so Opt=(QuadMac,TS,CalcFC) is usually a good choice.
Transition state searches with ONIOM
Gaussian09 uses microiterations in the optimisation of transition states with ONIOM. This is different from Gaussian03. The use of microiterations makes the calculation faster, but in many cases the exact transition state can be difficult to find (the structure "jumps" around in an area very close to the true transition state, but the optimisation procedure does not succeed in localising the transition state). A solution is to turn off the microiterations with the "noMicro" option.
ONIOM tips
- You can use iop(4/33=2) to print out all MM contributions.
- To get out the model system from an oniom calculation include iop(1/33=2). In the output file, search for: "ONIOM: generating point 2 high level on model system". All atoms without the label -1 belong to the model system (this list includes also the link hydrogens).
- There are two options If you want to check other force fields which are not implemented in Gaussian09
Option1: You can do one-dimensional potential energy surface scans over a reaction coordinate, where the constrained structure optimisations at the potential energy surface can be performed with ONIOM(DFT:UFF), ONIOM(DFT:DREIDING) or ONIOM (DFT:AMBER), where the iop(1/33=2) option helps to obtain QM and MM parts of the molecule. Then you can do MM optimisations on the QM and MM parts of the molecule separately with MACROMODEL or TINKER, where they have various force filed parameters. Finally, total energy can be calculated using;
Etot(QM,MM) = EQM(QM) + EMM(QM,MM) − EMM(QM) (the total energy of the ONIOM approach)
Option2: ´External´ keyword.
There are two examples in the G09 test suite which make use of the ´External´ keyword.
(1) Test726 uses a shell script to run g09 as an external process with UFF. The script ´extuffex´ is stored in the g09/bsd directory.
(2) Test769 uses a perl script, mygau, which is stored in the g09/bsd subdirectory. Note that you will likely need to edit mygau to change the path to perl, rather than /usr/local/bin, because most machines store it in /usr/bin.
In both cases you could use ´External´ keyword in the input file.
Eg:-
- ONIOM(B3LYP/6-31G(d):AM1:External="Script_name") Opt
How to use external MM program:
This is not fully implemented in G09. The Gau_External is an example where we use Gaussian as an external program, specifically 'Gau_External' script helps to run MOPAC as external procedure for ONIOM. You would need to make a copy of this script and modify it to create an input file for other MM programes (Eg. TINKER), and run the MM program with this input, and then use the output from the MM programe to generate result files suitable for processing by Gaussian.
- In Gaussview: When selecting MM regions in QM/MM methods you can change the status (H/L) of several atoms at same time, by:i) select the upper atom in the list ii) go to the last one and press control+shift 3) change the last atom status to low level.
- GaussView can help assigning atom types and it will also add point charges according to the AMBER convention for atoms in standard protein residues. For other molecules that are not standard residues, there is no information about appropriate charges (GaussView will try to assign reasonable atom types), so one would need to derive charges for the atoms in such molecule, or look in the literature to see if other researchers have derived charges and parameters for such molecule.
- Only the "amber96" version is internally hardwired in the G09 code and that is the one selected when using the AMBER keyword. The "amber.prm" file in the "g09" directory is simply a text file that shows this same parameter file that is internally hardwired in the code, but the program does not read that file. It is included for illustrative purposes so one see examples of how to specify MM parameters manually. The files that you mentioned (amber98.prm and oplsaa.prm) unfortunately are not ready to be used.
Remember that you can make modifications to the internally hardwired parameter file by adding the modifications to the input file, for instance:
- ONIOM(<Method>/<Basis Set>:AMBER=SoftFirst) geom=connectivity
<Title>
<Charge and Multiplicity> <Coordinates>
<Connectivity table>
<MM Parameter modifications>
With the "SoftFirst" option, the program will use the internally hardwired parameter file for AMBER and add the modifications given in the input file, with the latter taking preference over the internally hardwired parameters in case of duplicates.
If you do not want for G09 to use the internally hardwired parameter file but to read the full parameter table from an external file, you can do so with the "SoftOnly" option, which takes only the MM parameter specifications given in the input file and nothing from the internally hardwired parameter file. For example:
- ONIOM(<Method>/<Basis Set>:AMBER=SoftOnly) geom=connectivity
<Title>
<Charge and Multiplicity> <Coordinates>
<Connectivity table>
@/home/fernando/MM/file_with_parameters.prm /N
The "@" character at the beginning of the line means to include the content of a text file in the input, that way one does not need to type or paste the content of another file in the input, the program will automatically look for that other file and include it in this section of the input file. By default, the full content of the "included file" will be printed at the beginning of the output file of this job. If you do not want this behavior, you can use the "/N" at the end of the line as shown above, which means that the file will be included in the input but the content would not be dumped in the beginning of the output file.
Note that there are some differences in the format of the MM parameter specification between AMBER and Gaussian, so AMBER ".dat" files with the whole parameter table specification cannot be used directly in Gaussian, they need to be converted to the correct format. The documentation about MM methods in the Gaussian manual contains detailed descriptions of the input formats for each type of MM parameter. Below you can see a few examples that compare how they appear in the AMBER parameter files and how they should be input in Gaussian.
a) Bond Stretches:
-AMBER
C -CM 410.0 1.444
-GAUSSIAN
HrmStr1 C CM 410.0 1.444
b) Angle Bends:
-AMBER
CB-C -NA 70.0 111.30
-GAUSSIAN
HrmBnd1 CB C NA 70.0 111.30
c) Torsions. Wildcards in AMBER format are represented by "X" while the asterisk "*" is used in Gaussian instead. The format used in Gaussian is a compact way to write a series.
-AMBER
X -CK-NB-X 2 20.00 180.0 2.
-GAUSSIAN
AmbTrs * CK NB * 0 180 0 0 0.000 20.000 0.000 0.000 2.0
-AMBER
H -N -C -O 1 2.50 180.0 -2.
H -N -C -O 1 2.00 0.0 1.
-GAUSSIAN
AmbTrs H N C O 0 180 0 0 2.000 2.500 0.000 0.000 1.0
d) Improper Torsions:
-AMBER
X -X -C -O 10.5 180. 2.
-GAUSSIAN
ImpTrs * * C O 10.5 180.0 2.0
e) van der Waals.
AMBER uses a line in the parameter file that lists atoms with equivalent van der Waals parameters. In Gaussian, there is no such an entry and all the van der Waals parameters must be listed explicitly.
-AMBER
HS 0.6000 0.0157
-GAUSSIAN
VDW HS 0.6000 0.0157
You can do more comparisons by looking at the "parm96.dat" parameter file in the AMBER distribution in relation to the "amber.prm" file in the G09 distribution. Note that Gaussian does not read from this "amber.prm" file, so you should not make any changes to it. The parameters are internally stored in the program. The file "amber.prm" is only included in the G09 distribution for informative purposes.
Issues concerning the coordinate system
Some extra checks of the completeness of the coordinate system were added in Gaussian09. Thus in some cases you might get the error message "Error in internal coordinate system" where it did not appear in Gaussian03.
Typical cases where this error message appear, are for linear system. E.g. when three covalently bond atoms are arranged on a line. For organometallic systems a typical case is for carbonyl ligands coordinated to a metal centre, M-CO. If your system contains any linearly arranged atoms and you get this error message, a solution might be to add linear bend coordinates to these atoms. In your modredundant section you write:
- L atom1 atom2 atom3 -1 A
Solvent Effects
As with Gaussian03, the SCRF keyword requests that a calculation be performed in the presence of a solvent by placing the solute in a cavity within the solvent reaction field.
- The integral equation formalism variant IEFPCM, is the default SCRF method. It has not changed from Gaussian03, BUT the formalism used and its implementation has changed. That is: you will get different results in G03 and G09, using the same method.
- The default RADII used in G03 was UAO , whilst now is UFF with spheres placed by default on all Hydrogen atoms. No need to use the SPHEREONH= keyword, except if you use UAO or another radii model that does not have them explicitly.
- It seems to give no "convergence failure" problems, but maybe its too early to say so.
- It is able to perform frequency calculations in solvent, giving enthalpies, free energies, ZPE corrections...
- The G09 output for PCM gives only one energy which corresponds, by default, to the electrostic energy. The non-electroestatic components of the energy of solvation are not included. If you wish to include them you will need to include "scrf=(read,PCM,...)" in the command line, and Cav Dis Rep (in separate lines) at the end of the input.
- The SMD method includes automatically the non-electrostatic terms.
The RI-DFT method
From the Gaussian manual:
Gaussian 09 provides the density fitting approximation for pure DFT calculations. This approach expands the density in a set of atom-centered functions when computing the Coulomb interaction instead of computing all of the two-electron integrals. It provides significant performance gains for pure DFT calculations on medium sized systems too small to take advantage of the linear scaling algorithms without a significant degradation in the accuracy of predicted structures, relative energies and molecular properties. Gaussian 09 can generate an appropriate fitting basis automatically from the AO basis, or you may select one of the built-in fitting sets.
What all this means is that for all calculations that have more than about 1000 basis functions, it will speed up the calculation greatly if you use the RI-DFT method. Generally, the greater the number of basis functions, the better it compares to normal DFT. The accuracy of this approximation is very good, and generally you can afford to use a much higher basis set than you would with a normal DFT calcualtion (so you get rid of BSSE errors etc.).
How to set up an RI-DFT calculation.
You have to use a pure DFT functional... no hybrids allowed...
You also have to specify TWO basis sets... The first one is the normal basis set specification (GENECP is supported for complicated setups with transition metals and ECPs), then the second one is the "density fitting basis", which you can specify yourself or Gaussian can generate automatically.
For example:
# PBEPBE/TZVP/TZVPFit
would specify a calculation with the PBEPBE functional, TZVP "normal" basis set, and TZVPFit "density fitting" basis set. The / sign between each section is essential.
For an example of an input file, and more information on how to set up the calculation, see the following page... RI-DFT input