1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78import numpy as np
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
def f2(x,y):
return np.exp(x*x+(y-2)*(y-2))
def hx1(x,y):
return 2*x*np.exp(x*x+(y-2)*(y-2))
def hx2(x,y):
return 2*(y-2)*np.exp(x*x+(y-2)*(y-2))
X1 = np.arange(0,1,0.1)
X2 = np.arange(0,1,0.1)
X1, X2 = np.meshgrid(X1, X2) # 生成xv、yv,将X1、X2变成n*m的矩阵,方便后面绘图
#pdb.set_trace()
Y = np.array(list(map(lambda t : f2(t[0],t[1]),zip(X1.flatten(),X2.flatten()))))
Y.shape = X1.shape # 1600的Y图还原成原来的(40,40)
x1 = 1
x2 = 1
alpha = 0.01
GD_X1 = [x1]
GD_X2 = [x2]
GD_Y = [f2(x1, x2)]
y_change = f2(x1, x2)
iter_num = 0
while (y_change > 1e-20):
tmp_x1 = x1 - alpha * hx1(x1, x2)
tmp_x2 = x2 - alpha * hx2(x1, x2)
tmp_y = f2(tmp_x1, tmp_x2)
f_change = np.absolute(tmp_y - f2(x1, x2))
y_change=f_change
x1 = tmp_x1
x2 = tmp_x2
GD_X1.append(x1)
GD_X2.append(x2)
GD_Y.append(tmp_y)
iter_num += 1
print(u"最终结果为:(%.5f, %.5f, %.5f)" % (x1, x2, f2(x1, x2)))
print(u"迭代过程中X的取值,迭代次数:%d" % iter_num)
print(GD_X1)
fig = plt.figure(facecolor='w', figsize=(20, 18))
ax = Axes3D(fig)
ax.plot_surface(X1, X2, Y, rstride=1, cstride=1, cmap=plt.cm.jet)
ax.plot(GD_X1, GD_X2, GD_Y, 'ko-')
ax.set_xlabel('x')
ax.set_ylabel('y')
ax.set_zlabel('z')
ax.set_title(u'函数;\n学习率:%.3f; 最终解:(%.3f, %.3f, %.3f);迭代次数:%d' % (alpha, x1, x2, f2(x1, x2), iter_num))
plt.show()
```js