The part that incorporates momentum:

self.g1 *= self.b1
s = -self.step_size * self.inv_hv(self.g1 + self.grad)

Without momentum, the call to self.inv_hv would be just self.inv_hv(self.grad).

self.g1 is my gradient first moment accumulator and self.b1 (from Adam's beta_1) is my momentum parameter (here 0.75 for image synthesis but for neural network training it was 0.9). self.inv_hv() is the two-loop recursion function. This is Nesterov momentum so the search direction is b1 * g1 + self.grad (previous step's gradient) instead of just g1. Nesterov momentum did better in practice on the image synthesis problem than classical momentum. Importantly, I still form y using the actual gradient, not g1. The scaling matrix is another thing altogether from momentum and is applied completely separately, later.

I applied the scaling matrix in a really weird way that I haven't heard of anyone else doing. First I form y as grad - self.grad (the gradient at the current params minus the previous gradient), apply damping, then multiply y by sqrt(self.g2). I noticed, when trying out different damping methods, that reducing the scale of y (leaving s alone) increased the effective learning rate, and vice versa. Further, this effect was per-parameter, meaning I could pick out particular image regions, multiply the corresponding region of y by a constant, and slow the learning rate down for those regions. So I divide y by the Adagrad scaling matrix to apply it. :) I found this to produce more stable behavior than using the scaling matrix as the initial Hessian diagonal in the two-loop recursion! I am not using any kind of step rejection criterion currently so it really has to never take unrecoverably bad steps even if it means going a bit slower.

That is, y is formed as follows:

y = grad - self.grad
y = (1-phi)*y + phi*s
y *= sqrt(self.g2)

phi is a constant, and 0.2 is about right for it. Everyone else I have seen apply damping except oLBFGS uses some sort of rule to dynamically choose phi based on s and y. But that was complicated and, when I tried Powell damping proper, less stable (but faster) than just picking a constant phi. Notably, oLBFGS uses a constant as well (they call it lambda).

self.g2 of course is Adagrad's gradient second moment accumulator. Have you heard at all of this sort of modification of y with a scaling matrix?