本文已参与「新人创作礼」活动， 一起开启掘金创作之路

题目

符号微分法（失败）

算不出来。Hessian阵算出来是个常数，不知道怎么做：

import numpy as np

def f(x):
    prod1 = 0.2* np.matmul(x,x.transpose()) # 1/2x^tx
    A = np.array([[5,1,0,0.5],[1,4,0.5,0],[0,0.5,3,0],[0.5,0,0,2]])
    prod2 = np.matmul(x.transpose(),A) # xA
    prod3 = 0.25*(np.matmul(prod2,x))**2 # 0.25xAx^t
    print("f(x)=",prod3+prod1) # f(x)

def Df(x): #x+(xAx)Ax
    A = np.array([[5, 1, 0, 0.5], [1, 4, 0.5, 0], [0, 0.5, 3, 0], [0.5, 0, 0, 2]])
    num = np.matmul(x,np.matmul(A,x))
    return x+ np.multiply(num,np.matmul(A,x))

def D2f(x): #1+x(3AA)x
    A = np.array([[5, 1, 0, 0.5], [1, 4, 0.5, 0], [0, 0.5, 3, 0], [0.5, 0, 0, 2]])
    mat= np.matmul(A,A)
    prod1 = np.matmul(x.transpose(), np.matmul(mat,x))
    prod2 = np.array(1)
    return prod1+prod2

Pytorch自动微分

if __name__ == '__main__':

    import torch
    from torch.autograd import Variable

    def f(x):
        A = torch.Tensor([[5, 1, 0, 0.5], [1, 4, 0.5, 0], [0, 0.5, 3, 0], [0.5, 0, 0, 2]])
        z = torch.matmul(x.t(),torch.matmul(A,x))
        y = torch.matmul(x.t(),x)/2 + torch.pow(z,2)/4
        return y

    # 起始迭代点
    x = Variable(torch.Tensor([[1],[1],[1],[1]]), requires_grad=True)

    # 容忍误差
    eps = 1e-5

    # 初始梯度
    g1 = 1000 # 设置为最大
    i = 0

    while(g1 > eps): # 梯度的范数小于eps时收敛
        i = i+1
        # 计算梯度
        grad_x = torch.autograd.grad(f(x), x, create_graph=True,grad_outputs=torch.ones_like(f(x)))
        # 计算Hessian阵
        grad_grad_x = torch.autograd.grad(grad_x[0], x,grad_outputs=torch.ones_like(grad_x[0]))
        # 牛顿法更新
        x = Variable(x.data - grad_x[0].data / grad_grad_x[0].data, requires_grad=True)
        print("第",i,"次迭代：",end='')
        print("x=",x)
        print("y=",f(x).item())
        g1 = torch.as_tensor(torch.clone(grad_x[0]).detach()).norm()
        print("梯度的模|Dy/dx|2=",g1.item())

    print("牛顿法收敛于第",i,"次迭代;","函数最小值y_min=", f(x).item(),";此时x=",x)

小结

微分方法主要分为：（1）符号微分（复杂，需要手动求导）（2）数值微分 （3）自动微分（两者的综合） Pytorch中使用的是自动微分，没想到居然可以用来写数学作业。同学好像主要用的MATLAB中的dff函数实现的。要是手写数值微分的实现的话，还真有一些难度，日后再完善吧。

参考文献： blog.csdn.net/qq_23144435…