Hessian Checking from gradient

hessianChecking.py

from __future__ import print_function
import numpy as np

def hessianChecking(grad, hessian, w, *args):
    N = w.shape[0]

    H = hessian(w, *args)
    numericalHessian = np.zeros_like(H)

    # check symmetry of Hessian first
    assert (np.allclose(H, H.T, atol=1e-5))

    eps = 1e-4
    for i in range(N):
        for j in range(N):
            # create a small pertubation
            w1 = np.copy(w)
            w1[j] -= eps

            w2 = np.copy(w)
            w2[j] += eps
            
            loss1 = grad(w1, *args)[i]
            loss2 = grad(w2, *args)[i]

            numericalHessian[i,j] = (loss2 - loss1)/(2*eps)
        
    assert (np.allclose(numericalHessian, H, atol=1e-5))

def quadGrad(w, Q, b):
    return np.matmul(Q, w) + b

def quadHessian(w, Q, b):
    return Q

Q = np.array([[1,0],[0,1]]).astype(float)
b = np.array([[3,4]]).T.astype(float)
w = np.array([[1,20]]).T.astype(float)

hessianChecking(quadGrad, quadHessian, w, Q, b)