![\fbox{
\parbox[h][\height][l]{12cm}{
\small
\noindent
{\bf Reference literature:...
...ewblock {\em J.Phys.Chem.}, {\bf
100},\hspace{0.25em}19771, (1996).
\end{list}}}](img18.png)
A short input file for the location of a molecular geometry corresponding to minimum molecular energy (using a second-order optimization method) and performing a vibrational analysis at the stationary point, as well as determining the nuclear shielding constants at the optimized geometry, will look like (assuming a Hartree-Fock wave function):
**DALTON INPUT .OPTIMIZE *OPTIMIZE .2NDORD **WAVE FUNCTIONS .HF **PROPERTIES .VIBANA .SHIELD **END OF DALTON INPUT
The keyword .OPTIMIZE by default signals a search for a minimum, that is, a stationary point on the molecular surface with Hessian index 0. With the *OPTIMIZE module the second-order method has to be explicitly requested ( .2NDORD) [22], as first-order methods are the default. At the final, optimized geometry the input in the **FINAL section requests a vibrational analysis and the evaluation of nuclear shieldings.
In second-order methods the Hessian is calculated at every
geometry. The analytical Hessian naturally gives a better description
of the potential energy surface than the approximate Hessians of
first-order methods, and fewer steps will usually be needed to reach
the minimum [30]. However, the price one has to pay for this, is that each
iteration uses significantly more CPU time. Another advantage of
second-order methods, is that the program will automatically break the symmetry
of the molecule, if needed, in order to reach a minimum
(unless the user specifies that symmetry should not be broken). If
one wants to minimize a water molecule starting from a linear
geometry and using molecular symmetry, a first-order method will
``happily'' yield a linear geometry with optimized bond lengths, whereas
the second-order method will correctly decide that symmetry should be
reduced to C
and minimize the energy of the molecule.
The bottom-line is, however, that while second-order methods
definitely are the more robust, they will generally be outperformed by
the much more time-efficient first-order methods.
If you are convinced that the symmetry of your starting geometry is the correct, you may take advantage of the fact there will only be forces in the totally symmetric modes of the molecules, and thus only evaluate the Hessian in the totally symmetric Cartesian symmetry coordinates. This can be achieved with the input:
**DALTON INPUT .WALK .TOTSYM **WAVE FUNCTIONS .HF **PROPERTIES .VIBANA .SHIELD **END OF DALTON INPUT
In case a vibrational analysis has been requested at the optimized geometry, the program will recalculate the complete molecular Hessian at the optimized geometry in order to be able to evaluate all frequencies. Note that only the *WALK module, indicated by the keyword .WALK, supports the .TOTSYM keyword. We also note that any requests for static polarizabilities, Cioslowski population analysis, dipole gradients or vibrational circular dichroism (VCD) will be ignored as all the necessary property vectors will not be available. However, at the final optimized geometry, the full Hessian will be evaluated, as may then also these properties.
For direct calculations--where the two-electron integrals are not stored on disc--the cost of a second-order geometry optimization algorithm may become prohibitively expensive, despite the rather few iterations normally needed to converge the geometry. In these cases, it may be advantageous to use first-order geometry optimization routines as illustrated in the following input:
**DALTON INPUT .OPTIMIZE **WAVE FUNCTIONS .HF **PROPERTIES .VIBANA .SHIELD **END OF DALTON INPUT
This will use the default method of .OPTIMIZE: a first-order optimization [23] in redundant internal coordinates [31,32,30] using the BFGS Hessian update. This means that only the energy and the gradient is calculated in each iteration, and an approximate Hessian is obtained using the gradients and the displacements. The initial Hessian is diagonal in the redundant internal coordinates [33]. To obtain the properties at the optimized geometry, the Hessian has to be calculated. By looking at the vibrational frequencies one can then verify that a true minimum has been reached (as all frequencies should be real, corresponding to a positive definite Hessian).
Several first-order methods have been implemented in DALTON, the recommended method being the Broyden-Fletcher-Goldfarb-Shanno (BFGS) Hessian updating scheme. Though other updates may perform better in certain cases, the BFGS update generally seems to be the most reliable and thus is the default method for minimization [30].
Without the calculation of the Hessian at the optimized geometry, one can never be sure that the geometry is indeed a minimum, and this is the main problem with first-order methods.
Both first- and second-order methods are based on trust-region optimization. The trust-region is the region of the potential energy surface where our quadratic model (based on an analytical or approximate Hessian) gives a good representation of the true surface. This region is given the shape of a hypersphere. In the trust-region based optimization algorithm, the radius of the trust-region is automatically updated during the optimization by a feedback mechanism. Occasionally, however, the potential surface may show locally large deviations from quadratic form. This will result in a large disagreement between the predicted energy change and the energy calculated at the new geometry. If this deviation is larger than a given threshold, the step will be reduced through a decrease in the trust radius and a simple line search will be performed. A quadratic function is fitted to the rejected energy and the previous energy and gradient. The minimum of this function is then used as the new step, and this process may be repeated.
Another consideration concerning optimizations, is the choice of coordinate system. Second-order methods does not seem to be very sensitive to this, and the *WALK module thus only uses Cartesian coordinates. However, the coordinate system may be crucial to the performance of first-order methods, especially when starting from approximate Hessians. *OPTIMIZE provides the normal Cartesian coordinates and redundant internal coordinates, the latter being highly recommended when using first-order methods. Redundant internal coordinates consists of all bond lengths, angles and dihedral angles in a molecule, and the redundancies are removed through projection. The default for second-order methods is Cartesian coordinates, while redundant internal coordinates is the default for first-order methods, but it is possible to override this with the keywords .REDINT and .CARTES respectively.
The
*OPTIMIZE module decides when a minimum is reached through a
set of convergence criteria. These will be
set automatically depending on the accuracy of the wave function
using three different criteria, namely the change in energy
between two iterations, the norm of the molecular
gradient and the norm of the
step. The convergence threshold for the change
in energy is the maximum of
Hartree and two times the
convergence threshold for the wave function, whereas for the norm of
the gradient and step it is the maximum of
a.u. and
two times the
convergence threshold for the wave function. Two out of these
three criteria have to decrease below their threshold before the program
declares convergence. These criteria ensure a good geometry
suitable for vibrational analysis, but may be unnecessary tight if
one only wants the energy.
The convergence criteria can be controlled through the input file. One
should keep in mind that while second-order methods are rather
insensitive to changes in these thresholds due to their quadratic
convergence, first-order methods may run into severe trouble if the
criteria are too close to the accuracy of the wave function and the
gradient. The quickest way to get somewhat looser convergence
criteria, is to use the keyword
.BAKER enforcing the convergence
criteria of Baker[34] (see the Reference Manual for
details). This will normally give energies converged to within
Hartree.
In case of non-variational
wave functions, where a numerical
gradient is used, the threshold for
convergence is reset to
a.u., due to possible numerical
instabilities occurring in the
finite difference methods employed to
determine the molecular gradient.
We also note that for the optimization of excited states where the state of interest is not the lowest of its symmetry, most of the parameters controlling the trust-region optimization method will be reset to much more conservative values, as one in such calculations will loose the possibility of back-stepping during the wave function optimization, as well as to ensure that the start wave function at the new geometry is within the local region of the optimized wave function.
Calculations involving large basis sets may be very time-consuming, and it becomes very important to start the optimization from a decent geometry. One way to solve this is through preoptimization, that is, one or more smaller basis sets are first used to optimize the molecule. *OPTIMIZE supports the use of up to ten different basis sets for preoptimization. This may also be used as a very effective way to get the optimized energy of a molecule for a series of basis sets. One may also want to calculate the energy of a molecule using a large basis set at the optimized geometry of a smaller basis set, and this is also supported in DALTON, but limited to only one single-point basis set. Both these features are illustrated in the input file below:
**DALTON INPUT .OPTIMIZE *OPTIMIZE .PREOPT 2 STO-3G 4-31G .SP BAS d-aug-cc-pVTZ **WAVE FUNCTIONS .HF **END OF DALTON INPUT
As one can see, two small basis sets have been chosen for preoptimization. Note that they will be used in the order they appear in the input file, so one should sort the sets and put the smaller one at the top for maximum efficiency. The ``main'' basis set, that is, the one that will be used for the final optimized geometry, is specified in the file MOLECULE.INP as usual. After the last optimization, the energy of the molecule will be calculated at the final geometry with the d-aug-cc-pVTZ basis. Please note that these features will work only when using basis sets from the basis set library.
It is possible to control the optimization procedure more closely by giving the program instructions on how to do the different parts of the calculation. Below we have listed all the modules that may affect the geometry optimization, as well as a short description of which part of the calculation that module controls. The reader is referred to the section describing each module for a closer description of the different options. These sections may, with the exception of the *OPTIMIZE and *WALK section, be specified in any of the modules **START, **EACH STEP and **PROPERTIES. The *OPTIMIZE and *WALK is to be placed in the **DALTON module, and will apply to the entire calculation.