10.7.1 Overview

In order to carry out a frequency analysis users must at a minimum provide a molecule within the $molecule keyword and define an appropriate level of theory within the $rem keyword using the $rem variables EXCHANGE, CORRELATION (if required) (Chapter 4) and BASIS (Chapter 8). Since the default type of job (JOBTYPE) is a single point energy (SP) calculation, the JOBTYPE $rem variable must be set to FREQ.

It is very important to note that a vibrational frequency analysis must be performed at a stationary point on the potential surface that has been optimized at the same level of theory. Therefore a vibrational frequency analysis most naturally follows a geometry optimization in the same input deck, where the molecular geometry is obtained (see examples).

Users should also be aware that the quality of the quadrature grid used in DFT calculations is more important when calculating second derivatives. The default grid for some atoms has changed in Q-Chem 3.0 (see Section 5.5) and for this reason vibrational frequencies may vary slightly form previous versions. It is recommended that a grid larger than the default grid is used when performing frequency calculations.

The standard output from a frequency analysis includes the following.

• •

  Vibrational frequencies.

• •

  Raman and IR activities and intensities (requires $rem DORAMAN).

• •

  Atomic masses.

• •

  Zero-point vibrational energy.

• •

  Translational, rotational, and vibrational, entropies and enthalpies.

Several other $rem variables are available that control the vibrational frequency analysis. In detail, they are:

$molecule
   0  1
   O
   C  1  co
   F  2  fc  1  fco
   H  2  hc  1  hco  3  180.0

   co  =   1.2
   fc  =   1.4
   hc  =   1.0
   fco = 120.0
   hco = 120.0
$end

$rem
   JOBTYPE    opt
   METHOD     edf1
   BASIS      6-31+G*
$end

@@@

$molecule
   read
$end

$rem
   JOBTYPE    freq
   METHOD     edf1
   BASIS      6-31+G*
$end

Numerical values in the following discussion correspond to Example 10.7.1.1 and the referenced partition functions come from the textbook by McQuarrie & Simon. Note that Q-Chem assumes $T=298.15$ K and $P=1.00$ atm by default; see Section 10.7.2 for instructions on how to modify these choices.

The quantity listed as zero-point vibrational energy (ZPVE),

Zero point vibrational energy: 12.695 kcal/mol

corresponds to $\text{ZPVE}=\frac{1}{2}{\sum }_{K}hc{\nu }_K$. Within Q-Chem the vibrational entropy ($S_{\text{vib}}$) is computed from the bottom of the well rather than from the zero-point level, and thus includes the ZPE. This corresponds to a vibrational partition function

where $K$ indexes vibrational modes and ${\mathrm{\Theta }}_K=h{\nu }_K/k_B$ is the vibrational temperature of the $K$th mode. This expression can be altered by a factor of ${\mathrm{\Theta }}_K/2$ to start at the first vibrational level rather than the bottom of the well. The vibrational entropy is then

and for Example 10.7.1.1 it is reported as

Vibrational Entropy: 0.620 cal/mol.K

The vibrational enthalpy includes the ZPVE along with an additional temperature-dependent term:

It is reported as

Vibrational Enthalpy: 12.839 kcal/mol

The rotational and translational enthalpy are multiples of $RT$,

where the final contribution of $RT$ comes from the definition $H=U+PV$. The quantity $H_{\text{total}}$ reported as

Total Enthalpy: 15.209 kcal/mol

Translational and rotational entropies can be computed from the corresponding ideal-gas partition functions.

For thermochemical calculations performed in implicit solvent (using models described in Section 11.2), there is a subtle point with some controversy attached. ⁵⁶⁴ Ho J., Klamt A., Coote M. L.
J. Phys. Chem. A
(2010), 114, pp. 13442. Link ^{, 1116} Ribeiro R. F. et al.
J. Phys. Chem. B
(2011), 115, pp. 14556. Link ^{, 205} Casasnovas R. et al.
Int. J. Quantum Chem.
(2014), 114, pp. 1350. Link ^{, 548} Herbert J. M.
Wiley Interdiscip. Rev.: Comput. Mol. Sci.
(2021), 11, pp. e1519. Link This regards whether vibrational contributions to the free energy should be included or not. In nearly all cases, continuum solvation models are parameterized in order to reproduce experimental free energies of solvation (${\mathrm{\Delta }}_{\text{solv}}{G}^{\circ }$) using rigid gas-phase geometries for the solute molecules. There is a potential double-counting problem if $-T{S}_{\text{vib}}$ from a vibrational frequency calculation in implicit solvent is added, which might be avoided by instead using gas-phase harmonic frequencies for the solute. ⁵⁶⁴ Ho J., Klamt A., Coote M. L.
J. Phys. Chem. A
(2010), 114, pp. 13442. Link However, the difference between these procedures is negligible ($\sim 0.2$ kcal/mol). for the small-molecule data sets that are used to train implicit solvent models, ¹¹¹⁶ Ribeiro R. F. et al.
J. Phys. Chem. B
(2011), 115, pp. 14556. Link and only by using solution-phase vibrational frequencies can one obtain corrections to $S_{\text{vib}}$ arising from solvation-induced changes in geometry, which might be significant for larger molecules. Note that Q-Chem’s harmonic frequency engine assumes a gas-phase molecule (even in the presence of continuum solvent), so that rotational and translational contributions to the free energy are printed out in every case. These should not be included in solution-phase free energy differences.

For electronic structure methods whose analytic 2nd nuclear derivatives are not implemented in Q-Chem, finite-difference hessian calculations based on nuclear gradients will be performed. The Hessian matrix elements are calculated as

where $g_i=\partial E/\partial x_i$. In total, $6N_{\mathrm{atom}}$ gradient calculations at displaced structures are needed to construct the Hessian.

By default, the $6N_{\mathrm{atom}}$ gradient calculations are performed in serial, which is extremely inefficient for molecules of a large number of atoms. Taking advantage of the fact that these finite-difference calculations are independent of each other, Q-Chem provides the option to perform these calculations in molecular segments (i.e., atom blocks). The user will need to specify the number of atoms in each segment using $rem variable FD2ND_BLOCK_SIZE. The total number of segments will then be $\lfloor N_{\mathrm{atom}}/\mathrm{blocksize}\rfloor +1$ (note: the number of atoms in the last segment is allowed to be smaller than the user-specified block size). The user will then need to set up and execute a series of input with FD2ND_BLOCK_INDEX set to differnt values (starting from 0). Note that all these jobs should share the same value of FD2ND_BLOCK_SIZE. On a large computer cluster, these jobs can easily be executed simultaneously. These jobs will produce text files grdntDisp.0, grdntDisp.1, … and dipoleDisp.0, dipoleDisp.1… under the working directory. Note that IDERIV = 1 needs to be added if it is not set automatically.

After all the calculations for molecule segments finish, the user will need to perform an additional job with FD2ND_BLOCK_INDEX set to $-1$. This job will collect the results for all the segments (from the text files in the working directory), construct the Hessian, and complete the vibrational frequency analysis. This approach will greatly speed up the finite-difference frequency calculations for large molecules.

$molecule
0 1
O    -1.48239   -0.12728    0.00000
H    -1.90743    0.74236    0.00000
H    -0.52048    0.04679    0.00000
O     1.37121    0.12361    0.00000
H     1.71353   -0.35106    0.77313
H     1.71353   -0.35106   -0.77313
$end

$rem
JOBTYPE  FREQ
METHOD  B3LYP
BASIS  6-31+G(D)
SYM_IGNORE  TRUE
THRESH      14
SCF_CONVERGENCE  8
IDERIV      1
FD2ND_BLOCK_SIZE  3
FD2ND_BLOCK_INDEX 0
$end

@@@

$molecule
read
$end

$rem
JOBTYPE  FREQ
METHOD  B3LYP
BASIS  6-31+G(D)
SYM_IGNORE  TRUE
THRESH      14
SCF_CONVERGENCE  8
IDERIV      1
FD2ND_BLOCK_SIZE  3
FD2ND_BLOCK_INDEX 1
$end

@@@

$molecule
read
$end

$rem
JOBTYPE  FREQ
METHOD  B3LYP
BASIS  6-31+G(D)
SYM_IGNORE  TRUE
THRESH      14
SCF_CONVERGENCE  8
IDERIV      1
FD2ND_BLOCK_SIZE  3
FD2ND_BLOCK_INDEX -1
$end