用python實現二次函數的求導、求梯度、Hesse矩陣、求模

優化算法經常要用到導數、梯度、Hesse矩陣等，因此編寫了一個類用於實現這些功能

 

建立一個Function類，構造函數的參數是一個函數

其中part的功能是求偏導，var_index表示是第幾個變量，val表示這些變量的值

diff的功能是方便一元函數求導

私有函數__diff_是為了hesse編寫，傳入要求導的變量，返回一個求導后的Function類

hesse函數利用__diff_函數計算Hesse矩陣，返回一個matrix

為方便梯度進行線性代數運算，在定義grad函數時，返回值轉化成了numpy.array類型

from numpy import *


class Function:
    def __init__(self, _f):
        self.fun = _f

    def value(self, val):
        return self.fun(val)

    def part(self, var_index, val):
        a = self.fun(val)
        b = a + 1
        i = 0
        e = 2 ** 10 - 1
        e1 = 2 ** 10
        while 10 ** (-6) < e < e1 or i > -6:
            e1 = e
            a = b
            val_ = list(val)
            val_[var_index] += 10 ** i
            m = self.fun(val_)
            n = self.fun(val)
            b = (m - n) / 10 ** i
            i -= 2
            e = abs(b - a)
        return a

    def part_2(self, x_index, y_index, val):
        return self.__diff_(x_index).__diff_(y_index).value(val)

    def diff(self, val):
        a = self.fun(val)
        b = a + 1
        i = 0
        e = 2 ** 10 - 1
        e1 = 2 ** 10
        while 10 ** (-6) < e < e1 or i > -6:
            e1 = e
            a = b
            val_ = val + 10 ** i
            m = self.fun(val_)
            n = self.fun(val)
            b = (m - n) / 10 ** i
            i -= 2
            e = abs(b - a)
        return a

    def grad(self, val):
        g = array(val).astype('float')
        for i in range(0, g.size):
            g[i] = self.part(i, val)
        return array(g)

    def __diff_(self, index):
        def diff_f(vals):
            vals_ = list(vals)
            vals_[index] = vals_[index] + 10 ** (-6)
            m = self.fun(vals_)
            n = self.fun(vals)
            return (m - n) / 10 ** (-6)
        return Function(diff_f)

    def hesse(self, val):
        v = mat(val)
        G = mat(dot(v.T, v)).astype('float')
        for i in range(0, v.size):
            for j in range(0, v.size):
                p = self.part_2(i, j, val)
                G[i, j] = p
        return G

    def norm(self, val):
        s = 0
        for x in self.grad(val):
            s += x ** 2
        return sqrt(s)