- Search
- Download PDF

(May 16, 2021)

The concept of a reaction path is chemically intuitive (a pathway from
reactants to products) yet somewhat theoretically ambiguous because most
mathematical definitions depend upon the chosen coordinate system. Stationary
points on a potential energy surface are independent of this choice, but the
path connecting them is not, and there exist various mathematical definitions
of a “reaction path”. Q-Chem uses the intrinsic reaction coordinate (IRC)
definition, as originally defined by Fukui,
^{
333
}
J. Phys. Chem.

(1970),
74,
pp. 4161.
Link
which has come to
be a fairly standard choice in quantum chemistry. The IRC is essentially
sequence of small, steepest-descent paths going downhill from the transition state.

The reaction path is most unlikely to be a straight line and so by taking a
finite step length along the direction of the gradient you will leave the
“true” reaction path. A series of small steepest descent steps will zig-zag
along the actual reaction path (a behavior known as “stitching”). Ishida
*et al.*
^{
488
}
J. Chem. Phys.

(1977),
66,
pp. 215.
Link
developed a predictor-corrector algorithm, involving
a second gradient calculation after the initial steepest-descent step, followed
by a line search along the gradient bisector to get back on the path, and this
algorithm was subsequently improved by Schmidt *et al.*.
^{
972
}
J. Am. Chem. Soc.

(1985),
107,
pp. 2585.
Link
This is the method that Q-Chem adopts. It cannot be used for the first
downhill step from the transition state, since the gradient is zero, so instead
a step is taken along the Hessian mode whose frequency is imaginary.

The reaction path can be defined and followed in $Z$-matrix coordinates,
Cartesian coordinates or mass-weighted Cartesian coordinates. The latter
represents the “true” IRC as defined by Fukui.
^{
333
}
J. Phys. Chem.

(1970),
74,
pp. 4161.
Link
If the
rationale for following the reaction path is simply to determine which local
minima are connected by a given transition state, which, is arguably the major
use of IRC algorithms, then the choice of coordinates is irrelevant. In order
to use the IRC code, the transition state geometry and the exact Hessian must
be available. These must be computed via two prior calculations, with
JOBTYPE = TS (transition structure search) and
JOBTYPE = FREQ (Hessian calculation), respectively. Job
control variables and examples appear below.

An IRC calculation is invoked by setting JOBTYPE = RPATH in
the *$rem* section, and additional *$rem* variables are described below. IRC
calculations may benefit from the methods discussed in
Section 9.2 for obtaining good initial guesses for
transition-state structures.

RPATH_COORDS

Determines which coordinate system to use in the IRC search.

TYPE:

INTEGER

DEFAULT:

0

OPTIONS:

0
Use mass-weighted coordinates.
1
Use Cartesian coordinates.
2
Use *Z*-matrix coordinates.

RECOMMENDATION:

Use the default. Note that use of $Z$-matrix coordinates requires that geometries be input in $Z$-matrix format.

RPATH_DIRECTION

Determines the first direction of the eigenmode to follow. This will not usually be
known prior to the Hessian diagonalization.

TYPE:

INTEGER

DEFAULT:

1

OPTIONS:

1
Descend in the positive direction of the eigenmode, then restart in the negative direction.
-1
Descend in the negative direction of the eigenmode, then restart in the positive direction.

RECOMMENDATION:

It is usually not possible to determine in which direction to
go *a priori*, so both directions are automatically
considered. A job that reads in the final geometry from
the reaction path job will use the final step from
the second direction.

RPATH_MAX_CYCLES

Specifies the maximum number of points to find on the reaction path.

TYPE:

INTEGER

DEFAULT:

20

OPTIONS:

$n$
User-defined number of cycles.

RECOMMENDATION:

Use more points if the minimum is desired, but not reached using the default.

RPATH_MAX_STEPSIZE

Specifies the maximum step size to be taken (in 0.001 a.u.).

TYPE:

INTEGER

DEFAULT:

150
corresponding to a step size of 0.15 a.u..

OPTIONS:

$n$
Step size = $n$/1000 a.u.

RECOMMENDATION:

None.

RPATH_TOL_DISPLACEMENT

Specifies the convergence threshold for the step. If a step size is chosen by
the algorithm that is smaller than this, the path is deemed to have reached the
minimum.

TYPE:

INTEGER

DEFAULT:

5000
Corresponding to 0.005 a.u.

OPTIONS:

$n$
User-defined. Tolerance = $n$/1000000 a.u.

RECOMMENDATION:

Use the default. Note that this option *only* controls the
threshold for ending the RPATH job and does nothing to
the intermediate steps of the calculation. A smaller value will
provide reaction paths that end closer to the true minimum.
Use of smaller values without adjusting RPATH_MAX_STEPSIZE,
however, can lead to oscillations about the minimum.

RPATH_PRINT

Specifies the print output level.

TYPE:

INTEGER

DEFAULT:

2

OPTIONS:

$n$

RECOMMENDATION:

Use the default, as little additional information is printed at higher levels.
Most of the output arises from the multiple single point calculations that are
performed along the reaction pathway.

$molecule 0 1 C H 1 1.20191 N 1 1.22178 2 72.76337 $end $rem JOBTYPE ts BASIS sto-3g METHOD hf $end @@@ $molecule read $end $rem JOBTYPE freq METHOD hf BASIS sto-3g SCF_GUESS read $end @@@ $molecule read $end $rem JOBTYPE rpath BASIS sto-3g METHOD hf SCF_GUESS read RPATH_MAX_CYCLES 50 $end