Example: LT, Frequencies, TS, and IRC: HCN¶
- For a sequence of intermediates, each defined by a fixed angle H-C-N between the linear extremes HCN and CNH, the remaining geometrical parameters are optimized, giving a Linear Transit point-by-point scan of the energy curve of the Hydrogen atom traveling from one end of the CN fragment to the other. This is a useful way to get a reasonable first guess of the Transition State.
- At the approximate TS a Frequencies calculation is performed to obtain a fairly accurate Hessian for the next calculation.
- A TS search is carried out, using the computed Hessian. As variation, the TS search is repeated, first with the automatic (internal) Hessian (based on force fields) and then also with a constraint applied.
- A full IRC scan of the full path, starting from the TS, down to the two minima.
The first calculation is a Linear Transit where the Hydrogen atom moves from one side of CN to the other by a parametrized step-by-step change of the angle H-C-N. The other coordinates of the system are optimized along the path.
In the atoms block, one coordinate value is represented by an identifier (th). In the geovar block this is assigned two values, implying that it is a Linear Transit parameter. The initial and final values for the parameter are given.
Since the geometry block does not have OPTIM SELECTED, all other coordinates are optimized for each of the 10 Linear Transit points.
The subkey iterations in the geometry block carries two arguments: the first is the maximum number of optimization steps (per LT point). The second is the number of LT points to compute in this run: 4. This implies that only a part of the 10-point path defined by the LT parameter(s) will be scanned. The remainder will be done in a follow-up run to illustrate usage of the restart facility.
$ADFBIN/adf <<eor Title HCN Linear Transit, first part NoPrint SFO, Frag, Functions, Computation Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.3 0 0 3 H 1 2 0 1.0 th 0 End Basis Type DZP End Symmetry NOSYM BeckeGrid Quality good End Geometry Branch Old LinearTransit 10 Iterations 30 4 Converge Grad=3e-2, Rad=3e-2, Angle=2 END Geovar th 180 0 End End Input eor mv TAPE21 t21.LT rm logfile
The NoPRINT key turns off a lot of default output. There are several PRINT and NOPRINT options; see the User’s Guides for details.
Since the geometry changes from linear to planar (and finally back to linear again), the symmetry must be given explicitly in the input file. Otherwise the program would find a C(lin) symmetry for the initial geometry and assume that this symmetry is preserved throughout. This would of course result in an error abort when the first LT step is carried out, breaking the linear symmetry.
The here specified symmetry (NOSYM: no symmetry at all) is not the true symmetry of the complete path C(s) but a subgroup. It is always allowed to specify a lower symmetry than the actually present symmetry. Such may be necessary (for instance when the true symmetry cannot be handled by adf) or in special cases required for reasons of analysis. Generally speaking, however, we recommend to use the highest symmetry possible (given the case at hand and taking into account the symmetries recognizable by ADF) to boost performance.
Convergence thresholds in the geometry block are set less tight than the defaults: we need only a reasonable estimate of the path, but no highly converged geometries.
At the end of the run the tape21 result file is saved and renamed t21.LT to serve as restart file for the follow-up calculation.
$ADFBIN/adf <<eor Title HCN Linear Transit NoPrint SFO,Frag,Functions,Computation Restart t21.LT Fragments N t21.N C t21.C H t21.H End Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.3 0 0 3 H 1 2 0 1.0 th 0 End symmetry NOSYM BeckeGrid Quality good End Geometry Branch Old LinearTransit 10 Converge Grad=3e-2, Rad=3e-2, Angle=2 END Geovar th 180 0 End End Input eor rm TAPE21 logfile
From the restart file, supplied with the key restart, the program reads off that the first 4 points of the LT path have been done already and the scan is continued with LT point #5. The same path definition is supplied again, including the original starting values for the coordinates. The actual starting coordinates (for LT point #5) are read from the restart file. The input values, however, serve to define and verify consistency of the defined LT path and must therefore be supplied correctly.
The key noprint is used to suppress major parts of standard output: all information pertaining to the sfo analysis, all build-from-fragments information, and the lists of elementary functions in the basis sets and fit sets.
Frequencies at the estimated Transition State
From the results of the Linear Transit run we can sketch the energy barrier that H passes over when going from one side of the molecule to the other. This yields a reasonable guess for the Transition State.
To check that the so-obtained estimate is adequate we compute the frequencies in that geometry: one of them should be imaginary.
Apart from serving as a check that the TS estimate is not too bad, the computed Hessian will also serve in the follow-up calculation to obtain the true TS.
$ADFBIN/adf <<eor Title HCN Frequencies in LT max (approx), moderate precision NoPrint SFO,Frag,Functions,Computation BeckeGrid Quality good End Fragments N t21.N C t21.C H t21.H End Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.186 0 0 3 H 1 2 0 1.223 70 0 End Geometry Frequencies End End Input eor mv TAPE21 t21.Freq
Inspection of the output file shows that one of the frequencies is imaginary, as expected (printed as negative), signaling the proximity of the Transition State.
The TAPE21 result file of the calculation is renamed and saved. Later we will use it as a ‘restart’ file for a TS search, namely to supply the computed Hessian as the initial ‘guess’ of the Hessian in the (TS) optimization run.
Now carry out the Transition State search, starting from the lt-derived guess.
In this first attempt to find the TS, no use is made of the tape21 result file from the Frequencies run. That will be done in the next calculation.
$ADFBIN/adf <<eor Title HCN Transition State, automatic initial Hessian NoPrint SFO,Frag,Functions,Computation BeckeGrid Quality good End Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.186 0 0 3 H 1 2 0 1.223 70 0 End Fragments N t21.N C t21.C H t21.H End Geometry TransitionState End End Input eor rm TAPE21 logfile
The TS-search run type is specified in the geometryblock.
No symmetry is specified; the program determines the symmetry to be C(s) and consequently carries out the ts search in that symmetry.
TS search, using the Hessian from the Frequencies run
$ADFBIN/adf <<eor Title HCN Transition State, initial Hessian from Freq run NoPrint SFO,Frag,Functions,Computation Restart t21.Freq Save TAPE13 BeckeGrid Quality good End Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.186 0 0 3 H 1 2 0 1.223 70 0 End Fragments N t21.N C t21.C H t21.H End Geometry TransitionState End End Input eor mv TAPE13 t13.TS rm TAPE21 logfile
The CheckPoint file TAPE13, at normal termination automatically deleted by the program, is here saved, using the SAVE key. TAPE13 is as good a restart file as TAPE21 is, but it is a lot smaller. TAPE21 contains a large amount of information for analysis purposes, while TAPE13 contains essentially only restart-type data.
The input is identical to the previous one, except for the restart file. This is used here to provide the Hessian computed in the Frequencies run as the start-up Hessian for the ts optimization. At the same time the atomic coordinates are read off from the restart file and override the values in the input file. This latter aspect could have been suppressed; see the User’s Guide for using the restart key.
Constrained TS search
Finally the ts search where one coordinate is kept frozen, to illustrate a constrained optimization.
$ADFBIN/adf <<eor Title HCN constrained TS search NoPrint SFO,Frag,Functions,Computation Restart t21.Freq BeckeGrid Quality good End Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 rNC 0 0 3 H 1 2 0 1.223 70 0 End GeoVar rNC=1.186 F End Fragments N t21.N C t21.C H t21.H End Geometry TransitionState End End Input eor rm TAPE21 logfile rm t21.Freq
The geovar key specifies that the nc distance, rNC has the initial value 1.15 and remains frozen (‘F’).
The fact that the optimization is now carried out in a different subspace of atomic coordinates does not prevent us from using the t21.Freq restart file to supply the initial Hessian.
IRC scan of the reaction path
The IRC calculation is split in three steps, to illustrate the Restart facility applied to the IRC functionality.
In the first only a few points are computed, along one of the two paths leading from the TS to the adjacent minima. Since no explicit directives are given in the input to specify the direction of the first path, the so-called ‘forward’ path is taken. The definition of which is ‘forward’ and which is ‘backward’ is in fact quite arbitrary and is determined by the program. See the User’s Guide for details.
The saved TAPE13 file from one of the TS calculations is used as restart file. This provides (a) the optimized coordinates of the TS as starting point, (b) the initial Hessian to guide the point-by-point optimizations along the IRC path, and (c) the eigenvector of the lowest Hessian eigenvalue to define the initial direction of the IRC path.
The TAPE13 file from this partial IRC scan is saved to serve as start-up file for the next calculations, which will continue the IRC scan.
In the Geometry key block, the run type is set to IRC and the ‘Points’ option is used to limit the number of IRC points to compute.
$ADFBIN/adf << eor Title HCN IRC partial path (forward) NoPrint SFO,Frag,Functions, Computation BeckeGrid Quality good End Restart t13.TS Save TAPE13 Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.186 0 0 3 H 1 2 0 1.223 70 0 End Fragments N t21.N C t21.C H t21.H End Geometry IRC Points=5 End End Input eor mv TAPE13 t13.IRC_1 rm TAPE21 logfile
The IRC is continued in the next calculation, using the TAPE13 file from the previous one as restart file. From this file, the program reads the IRC path information computed sofar. By default, it would continue on the ‘forward’ path, since that was not yet finished. However, in the Geometry key block, we now specify not only that a limited number of points is to be computed in this run (5 again), but we instruct the program also to compute only points on the ‘backward’ path.
$ADFBIN/adf << eor Title HCN IRC partial part (backward) NoPrint SFO,Frag,Functions, Computation Restart t13.IRC_1 Save TAPE13 BeckeGrid Quality good End Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.186 0 0 3 H 1 2 0 1.223 70 0 END Fragments N t21.N C t21.C H t21.H End Geometry IRC Points=5 Backward End End Input eor mv TAPE13 t13.IRC_2 rm TAPE21 logfile
In the third IRC run, the IRC scan is finished. We start with the TAPE13 file from the previous run and set a maximum of 70 IRC points to compute (which turns out to be sufficient for the complete IRC scan). The program starts on the forward path, continuing where the first (not the previous) had stopped after 5 points, completes the forward path, and then continues on the backward path, starting where the second IRC run had stopped. Both paths are finished and a summary of the path characteristics is printed in the final part of the output.
$ADFBIN/adf << eor Title HCN IRC completion NoPrint SFO,Frag,Functions, Computation Restart t13.IRC_2 BeckeGrid Quality good End Atoms Internal 1 C 0 0 0 0 0 0 2 N 1 0 0 1.186 0 0 3 H 1 2 0 1.223 70 0 End Fragments N t21.N C t21.C H t21.H End Geometry IRC Points=70 End End Input eor