用python实现二次函数的求导、求梯度、Hesse矩阵、求模

优化算法经常要用到导数、梯度、Hesse矩阵等,因此编写了一个类用于实现这些功能

 

建立一个Function类,构造函数的参数是一个函数

其中part的功能是求偏导,var_index表示是第几个变量,val表示这些变量的值

diff的功能是方便一元函数求导

私有函数__diff_是为了hesse编写,传入要求导的变量,返回一个求导后的Function类

hesse函数利用__diff_函数计算Hesse矩阵,返回一个matrix

为方便梯度进行线性代数运算,在定义grad函数时,返回值转化成了numpy.array类型

 1 from numpy import *
 2 
 3 
 4 class Function:
 5     def __init__(self, _f):
 6         self.fun = _f
 7 
 8     def value(self, val):
 9         return self.fun(val)
10 
11     def part(self, var_index, val):
12         a = self.fun(val)
13         b = a + 1
14         i = 0
15         e = 2 ** 10 - 1
16         e1 = 2 ** 10
17         while 10 ** (-6) < e < e1 or i > -6:
18             e1 = e
19             a = b
20             val_ = list(val)
21             val_[var_index] += 10 ** i
22             m = self.fun(val_)
23             n = self.fun(val)
24             b = (m - n) / 10 ** i
25             i -= 2
26             e = abs(b - a)
27         return a
28 
29     def part_2(self, x_index, y_index, val):
30         return self.__diff_(x_index).__diff_(y_index).value(val)
31 
32     def diff(self, val):
33         a = self.fun(val)
34         b = a + 1
35         i = 0
36         e = 2 ** 10 - 1
37         e1 = 2 ** 10
38         while 10 ** (-6) < e < e1 or i > -6:
39             e1 = e
40             a = b
41             val_ = val + 10 ** i
42             m = self.fun(val_)
43             n = self.fun(val)
44             b = (m - n) / 10 ** i
45             i -= 2
46             e = abs(b - a)
47         return a
48 
49     def grad(self, val):
50         g = array(val).astype('float')
51         for i in range(0, g.size):
52             g[i] = self.part(i, val)
53         return array(g)
54 
55     def __diff_(self, index):
56         def diff_f(vals):
57             vals_ = list(vals)
58             vals_[index] = vals_[index] + 10 ** (-6)
59             m = self.fun(vals_)
60             n = self.fun(vals)
61             return (m - n) / 10 ** (-6)
62         return Function(diff_f)
63 
64     def hesse(self, val):
65         v = mat(val)
66         G = mat(dot(v.T, v)).astype('float')
67         for i in range(0, v.size):
68             for j in range(0, v.size):
69                 p = self.part_2(i, j, val)
70                 G[i, j] = p
71         return G
72 
73     def norm(self, val):
74         s = 0
75         for x in self.grad(val):
76             s += x ** 2
77         return sqrt(s)

 

posted @ 2018-06-08 08:54  "kisetsu  阅读(9785)  评论(0编辑  收藏  举报