The second derivatives of the energy required for the Hessian, eq.~\ref{eq:hessian}, can be calculated either via finite differences, generating what are known as \textit{numerical} derivatives, or using perturbation theory, generating so-called \textit{analytical} derivatives~\cite{Pulay1969,Jorgensen1983,Baroni2001,Kadantsev2005}. 

A property of the calculations of the energy derivatives is that the numerical cost does not depend on the direction they are calculated. This can be readily seen in the case of finite differences, as the cost of calculating \(E(\vec{R}^1,\ldots,\vec{R}^j+\Delta^j,\ldots,\vec{R}^N)\) is essentially the same as computing \(E(\vec{R}^1 + \Delta^1,\ldots,\vec{R}^j+\Delta^j,\ldots,\vec{R}^N + \Delta^N)\). As discussed previously, this ability to compute matrix elements at a comparable numerical cost in any desired basis is an essential requirement of our method.

%If on the contrary, there is a preferred basis where the matrix is cheaper to compute, to %compute it in a different basis might offset the reduction in cost offered by compressed %sensing. 

A second property of both numerical and analytical derivatives that appears in variational  quantum chemistry formalisms like DFT or Hartree-Fock is that each calculation yields a full column of the Hessian, rather than a single matrix element. Again, this is easy to see in finite difference computations. We can write the second derivative of the energy as a first derivative of the force
%
\begin{equation}
\label{eq:force}
\frac{\partial^2E(\vec{R}^1,\ldots,\vec{R}^N)}{\partial R^A_i \partial R^B_j} =
- \frac{\partial F^B_j(\vec{R}^1,\ldots,\vec{R}^N)}{\partial R^A_i} \ .
\end{equation}
%
By the Hellman-Feynman theorem~\cite{Hellmann1937,Feynman1939}, a single energy calculation yields the forces acting over \emph{all} atoms, so the evaluation of eq.~\ref{eq:force} by finite differences for fixed \(A\) and \(i\) yields the derivatives for all values of \(B\) and \(j\), a whole column of the Hessian.  An equivalent result holds for analytic derivatives obtained via perturbation theory~\cite{Baroni2001,Kadantsev2005}. Thus, our compressed sensing procedure for this particular application focuses on measuring random columns of the quantum mechanical Hessian rather than individual random entries.

The full compressed sensing procedure applied to the calculation of a quantum mechanical Hessian is therefore implemented as the following:

\begin{enumerate}
\item Calculate approximate vibrational modes using molecular mechanics.
\item Transform the approximate modes using the DCT matrix.
\item Randomly select a few of the transformed modes.
\item Calculate the energy second derivatives along these random modes to yield a set of random columns of the quantum mechanical Hessian.
\item Apply compressed sensing to rebuild the full quantum mechanical Hessian in the basis of approximate vibrational modes.
\item Transform the full quantum mechanical Hessian back into the atomic coordinate basis.
\item Diagonalize the quantum mechanical Hessian to obtain the vibrational modes and frequencies.
\end{enumerate}

Fig.~\ref{fig:Anthracene} illustrates the results of applying our Hessian recovery procedure to anthracene (\(\textrm{C}_{14}\textrm{H}_{10}\)), a moderately-sized polyacene consisting of three linearly fused benzene rings.  The top panel illustrates the vibrational frequencies obtained by the compressed sensing procedure outlined above for different extents of undersampling of the true quantum mechanical Hessian.  Even sampling only 25\% of the columns yields vibrational frequencies that are close to the true quantum mechanical frequencies, and much closer than the molecular mechanics frequencies.  The middle panel illustrates the error in the vibrational frequencies from the true quantum mechanical frequencies.  Sampling only 30\% of the columns gives rise to a maximum frequency error of less than 3 cm\(^{-1}\), and sampling 35\% of the columns yields nearly exact recovery.  The bottom panel illustrates the error in the normal modes.  Once again, sampling only 30\% of the columns gives accurate recovery of all vibrational normal modes to within 1\%.  In short, our compressed sensing procedure applied to anthracene reduces the number of expensive quantum mechanical computations by a factor of three. The additional cost of the molecular mechanics computation and the compressed sensing procedure, which take a few seconds, is negligible in comparison with this reduction in computational time in the computation of the Hessian that for anthracene takes on the order of hours.