In order to check whether a stationary point, i.e., a point on the potential energy surface with a vanishing gradient, is a local minimum or not, the second derivatives of the energy with respect to nuclear displacement need to be calculated. This matrix of all possible second derivatives is usually referred to as the Hessian matrix. The Hessian can successively be diagonalized to obtain harmonic force constants, vibrational frequencies, and normal modes. Second derivatives of the energy can, of course, be calculated numerically (which has the same advantages and disadvantages as in the case of the gradient itself), or analytically. In the following, we will describe how the analytical Hessian is calculated at the level of Hartree--Fock (HF) theory, which leads to the coupled-perturbed HF equations. Successively, it is explained how a vibrational analysis is carried out based on the Hessian matrix calculated at an arbitrary level of theory.