  Second derivatives of mathematical models for real-world phenomena are
fundamental ingredients of a wide range of numerical simulation methods
including parameter sensitivity analysis, uncertainty quantification, nonlinear
optimization and model calibration. The evaluation of such Hessians often
dominates the overall computational effort. The combinatorial {\sc Hessian
Accumulation} problem aiming to minimize the number of floating-point
operations required for the computation of a Hessian turns out to be
NP-complete. We propose a dynamic programming formulation for the solution of
{\sc Hessian Accumulation} over a sub-search space. This approach yields
improvements by factors of ten and higher over the state of the art based on
second-order tangent and adjoint algorithmic differentiation.
