Obtaining the Hessian

For saddle point searches, the updating must allow the Hessian to develop negative eigenvalues, and the Powell or updates based on combining several methods are usually employed.  

The use of approximate Hessians within the NR method is known as pseudo-Newton-Raphson or variable metric methods. It is clear that they do not converge as fast as true NR methods, where the exact Hessian is calculated in each step, but if for example five steps can be taken for the same computational cost as one true NR step, the overall computational effort may be less. True NR methods converge quadratically near a stationary point, while pseudo-NR methods display a linear convergence. Far from a stationary point, however, the true NR method will typically also only display linear convergence. 

Pseudo-NR methods are usually the best choice in geometry optimizations using an energy function calculated by electronic structure methods. The quality of the initial hessian of course affects the convergence when an updating scheme is used. The best -choice is usually an exact Hessian at the first point, however, this may not be the most 

---

We are now able to obtain the Hessian matrix of the objective function S(k) which is denoted by H and is given by the following equation... 

Our next step for achieving a second order optimization procedure of the energy functional is to obtain the Hessian contribution, denoted by af due to the interactions between the quantum and classical subsystems. This is effectively done by performing linear transformations using configuration state function trial vectors and orbital trial vectors. The trial vectors are denoted and we obtain the following expressions... 

Compute the energy and gradients of this structnre. Obtain the Hessian matrix as a guess or by analytical or numerical computation. 

Model calculations were performed on the VAMP [24], DMOL [25, 26], and CASTEP [27] modules of the Materials Studio program package from Accelrys. Full geometry optimizations and vibrational frequency analyses were carried out in all electron approximation using in DMOL the BLYP [28, 29] functional in conjunction with the double-numeric-basis set with polarization functions (DNP) and the IR models were calculated from the Hessians [30], In CASTEP the gradient-corrected (GGA) PBE [31] functional was selected for the density functional theory (DFT) computations with norm conserving and not spin polarized approach [32], In the semi-empirical VAMP method we used the PM3 parameterization [33] from the modified neglect of diatomic differential overlap (NDDO) model to obtain the Hessians for vibrational spectrum models [30],... 

We next solve the secular equation F — I = 0 to obtain the eigenvalues and eigenvectors o the matrix F. This step is usually performed using matrix diagonalisation, as outlined ii Section 1.10.3. If the Hessian is defined in terms of Cartesian coordinates then six of thes( eigenvalues will be zero as they correspond to translational and rotational motion of th( entire system. The frequency of each normal mode is then calculated from the eigenvalue using the relationship ... 

The above formula is obtained by differentiating the quadratic approximation of S(k) with respect to each of the components of k and equating the resulting expression to zero (Edgar and Himmelblau, 1988 Gill et al. 1981 Scales, 1985). It should be noted that in practice there is no need to obtain the inverse of the Hessian matrix because it is better to solve the following linear system of equations (Peressini et al. 1988)... 

The Hessian matrix H(r) is defined as the symmetric matrix of the nine second derivatives 82p/8xt dxj. The eigenvectors of H(r), obtained by diagonalization of the matrix, are the principal axes of the curvature at r. The rank w of the curvature at a critical point is equal to the number of nonzero eigenvalues the signature o is the algebraic sum of the signs of the eigenvalues. The critical point is classified as (w, cr). There are four possible types of critical points in a three-dimensional scalar distribution ... 

INMs are obtained by diagonalizing the Hessian matrix generated by expanding the system potential energy to quadratic order in mass-weighted coordinates. 

There are a few points with respect to this procedure that merit discussion. First, there is the Hessian matrix. With elements, where n is the number of coordinates in the molecular geometry vector, it can grow somewhat expensive to construct this matrix at every step even for functions, like those used in most force fields, that have fairly simple analytical expressions for their second derivatives. Moreover, the matrix must be inverted at every step, and matrix inversion formally scales as where n is the dimensionality of the matrix. Thus, for purposes of efficiency (or in cases where analytic second derivatives are simply not available) approximate Hessian matrices are often used in the optimization process - after aU, the truncation of the Taylor expansion renders the Newton-Raphson method intrinsically approximate. As an optimization progresses, second derivatives can be estimated reasonably well from finite differences in the analytic first derivatives over the last few steps. For the first step, however, this is not an option, and one typically either accepts the cost of computing an initial Hessian analytically for the level of theory in use, or one employs a Hessian obtained at a less expensive level of theory, when such levels are available (which is typically not the case for force fields). To speed up slowly convergent optimizations, it is often helpful to compute an analytic Hessian every few steps and replace the approximate one in use up to that point. For really tricky cases (e.g., where the PES is fairly flat in many directions) one is occasionally forced to compute an analytic Hessian for every step. 

The normal modes of vibration u are obtained as solutions of the Hessian eigenvalue problem,... 

Let us now turn to the second derivatives. The Hessian matrix is divided into three parts the orbital - orbital part (oo) the configuration - configuration part (cc) and the so called Cl coupling part (co). For the (cc) part we obtain ... 

---