Lyapunov-certified Hessian-damped optimizer
Source paper: A Refined Parameter Condition in the Lyapunov Analysis of IGAHD arXiv:2608.28088 ⓘ · analyzed Sep 2, 2026
AI-generated research hypothesis, automatically tested. Not peer-reviewed.
Idea description
Replace the momentum update in a gradient optimizer by inertial motion plus a gradient-difference term, which discretely approximates Hessian-driven damping. Choose the damping coefficient and step size using the paper's refined stability inequality instead of the older restrictive bound, and adapt them whenever the estimated smoothness changes.
Formulas
Mathematical statement
The paper considers a convex function f: R^d -> R with L-Lipschitz gradient, meaning ||∇f(x)-∇f(y)|| <= L||x-y||, a step size s satisfying 0 < s <= 1/L, and a Hessian-driven damping parameter β > 0. Retaining two exact coefficients in the Lyapunov dissipation makes the relevant normalized symmetric 2x2 matrices converge to a positive-definite limit when (β L sqrt(s)-1)^2 < 1+sL(1-sL). Equivalently, the sufficient admissible region is 0 < β L sqrt(s) < 1 + sqrt(1+sL(1-sL)); the simpler conservative condition β < 2/(L sqrt(s)) is also sufficient. Use the gradient-difference approximation ∇f(x_k)-∇f(x_{k-1}) ≈ ∇²f(x_k)(x_k-x_{k-1}) as a discrete Hessian-damping force. The transfer hypothesis is that staying inside this region prevents oscillatory inertial divergence and permits a larger useful damping range than the old condition β < 2 sqrt(s).
Implementation notes
Integrate this at the optimizer update for each parameter tensor, not inside the network forward pass. Maintain x_prev, x, g_prev, and g; use α=3 initially, s as the gradient step size, and β as the gradient-difference damping coefficient. At step k, compute g=grad(loss,x), estimate local smoothness from the secant ratio r=||g-g_prev||/(||x-x_prev||+ε), update L_est=max(decayed_L_est,r), and clip L_est to a robust percentile or exponential moving maximum to avoid one-batch outliers. Set s=min(s, 0.99/L_est). Compute bmax=[1+sqrt(1+s*L_est*(1-s*L_est))]/(L_est*sqrt(s)); choose β=min(β_target, 0.9*bmax), or use β=0.5*0.9*bmax for a less aggressive first run. Then execute x_next=x+(1-α/k)*(x-x_prev)-s*g-β*(g-g_prev), and shift x_prev=x and g_prev=g. At k=1, initialize x_prev=x and g_prev=g, so the damping term is zero. Estimate the mathematical quantities L_est and the admissible boundary empirically; do not claim a global L unless the model has a certified smoothness bound. First test on a 2-layer MLP and a small Transformer on MNIST or CIFAR-10, comparing AdamW, Nesterov momentum, fixed-β IGAHD, and the adaptive controller at equal optimizer steps, parameters, and measured FLOPs. Pre-register three mechanism tests: (1) on convex quadratics with known largest Hessian eigenvalue L, divergence should begin near βL√s=1+sqrt(1+sL(1-sL)), with boundary error below 10%; (2) at equal stable step size, the proposed method should reduce loss faster than β=0 in the intermediate damping range, measured by area under the training-loss curve; (3) reducing s while holding β fixed should remain stable once the refined condition holds, whereas the old β<2√s rule would reject those settings. Remove only the β(g-g_prev) term as the ablation. The transfer is falsified if the adaptive method has no stability boundary near the predicted inequality on quadratics, or if gradient-difference damping consistently worsens equal-budget convergence after accounting for its extra gradient-memory cost.
Verification
This idea has not been verified yet.
Verification happens in two stages: Stage 1 — a mechanism check on a toy system confirms the claimed mathematical phenomenon reproduces; Stage 2 — a benchmark implements the idea on a real (small) neural network task and compares it against a tuned baseline over 8 paired seeds with a permutation test.
Artifacts
Artifacts unavailable.