Explore direct calculation of Hessian in quasi-newton optimisation approach

[daurer]

Starting point should be testing L-BFGS engine (which is similar but approximates the Hessian) on nearfield data and then potentially implement the other approach (https://opg.optica.org/oe/fulltext.cfm?uri=oe-33-14-30543&id=574129) and compare.

[jfowkes]

The main thrust of the approach in the paper is that rather than approximate the Hessian using gradients as we do in L-BFGS, they instead derive equations for the Hessian operator Hx (i.e. the result of applying H to an arbitrary vector x) and the Hessian bilinear form operator x^T H x. When I tried to do this exactly it was a complete mess, but they seem to make some approximations that make the calculations tractable. The main claim they make is that 1) they can speed up nonlinear CG (i.e. ML) by using a more accurate beta rather than the Polak-Ribiere approximation that we currently use, and 2) alternatively they can approximately solve Hx = -g (using e.g. linear CG) to get a Quasi-Newton method that should be more accurate than the L-BFGS approximation we use as it uses the exact Hessian operator in computing the descent direction.