在一般問題的優化中，最速下降法和共軛梯度法都是非常有用的經典方法，但最速下降法往往以」之」字形下降，速度較慢，不能很快的達到最優值，共軛梯度法則優於最速下降法，在前面的某個文章中，我們給出了牛頓法和最速下降法的比較，牛頓法需要初值點在最優點附近，條件較為苛刻。
我們選用了64維的二次函式來作為驗證函式，具體參見上書111頁。
採用的三種方法為：共軛梯度方法（fr格式）、共軛梯度法（prp格式）、最速下降法
# -*- coding: utf-8 -*-
"""created on sat oct 01 15:01:54 2016
@author: zhangweiguo
"""import sympy,numpy
import math
import matplotlib.pyplot as pl
from mpl_toolkits.mplot3d import axes3d as ax3
import sd#這個檔案裡有最速下降法sd的方法，參見前面的部落格
#共軛梯度法fr、prp兩種格式
def cg_fr(x0,n,e,f,f_d):
x=x0;y=;y_d=;
n = 1
ee = f_d(x0)
e=(ee[0]**2+ee[1]**2)**0.5
d=-f_d(x0)
a=sympy.symbol('a',real=true)
print '第%2s次迭代：e=%f' % (n, e)
while ne:
n=n+1
g1=f_d(x0)
f1=f(x0+a*f_d(x0))
a0=sympy.solve(sympy.diff(f1[0,0],a,1))
x0=x0-d*a0
ee = f_d(x0)
e = math.pow(math.pow(ee[0,0],2)+math.pow(ee[1,0],2),0.5)
g2=f_d(x0)
beta=(numpy.dot(g2.t,g2))/numpy.dot(g1.t,g1)
d=-f_d(x0)+beta*d
print '第%2s次迭代：e=%f'%(n,e)
return x,y,y_d
def cg_prp(x0,n,e,f,f_d):
x=x0;y=;y_d=;
n = 1
ee = f_d(x0)
e=(ee[0]**2+ee[1]**2)**0.5
d=-f_d(x0)
a=sympy.symbol('a',real=true)
print '第%2s次迭代：e=%f' % (n, e)
while ne:
n=n+1
g1=f_d(x0)
f1=f(x0+a*f_d(x0))
a0=sympy.solve(sympy.diff(f1[0,0],a,1))
x0=x0-d*a0
ee = f_d(x0)
e = math.pow(math.pow(ee[0,0],2)+math.pow(ee[1,0],2),0.5)
g2=f_d(x0)
beta=(numpy.dot(g2.t,g2-g1))/numpy.dot(g1.t,g1)
d=-f_d(x0)+beta*d
print '第%2s次迭代：e=%f'%(n,e)
return x,y,y_d
if __name__=='__main__':
'''g=numpy.array([[21.0,4.0],[4.0,15.0]])
#g=numpy.array([[21.0,4.0],[4.0,1.0]])
b=numpy.array([[2.0],[3.0]])
c=10.0
x0=numpy.array([[-10.0],[100.0]])
'''m=4
t=6*numpy.eye(m)
t[0,1]=-1;t[m-1,m-2]=-1
for i in xrange(1,m-1):
t[i,i+1]=-1
t[i,i-1]=-1
w=numpy.zeros((m**2,m**2))
w[0:m,0:m]=t
w[m**2-m:m**2,m**2-m:m**2]=t
w[0:m,m:2*m]=-numpy.eye(m)
w[m**2-m:m**2,m**2-2*m:m**2-m]=-numpy.eye(m)
for i in xrange(1,m-1):
w[i*m:(i+1)*m,i*m:(i+1)*m]=t
w[i*m:(i+1)*m,i*m+m:(i+1)*m+m]=-numpy.eye(m)
w[i*m:(i+1)*m,i*m-m:(i+1)*m-m]=-numpy.eye(m)
mm=m**2
mmm=m**3
g=numpy.zeros((mmm,mmm))
g[0:mm,0:mm]=w;g[mmm-mm:mmm,mmm-mm:mmm]=w;
g[0:mm,mm:2*mm]=-numpy.eye(mm)
g[mmm-mm:mmm,mmm-2*mm:mmm-mm]=-numpy.eye(mm)
for i in xrange(1,m-1):
g[i*mm:(i+1)*mm,i*mm:(i+1)*mm]=w
g[i*mm:(i+1)*mm,i*mm-mm:(i+1)*mm-mm]=-numpy.eye(mm)
g[i*mm:(i+1)*mm,i*mm+mm:(i+1)*mm+mm]=-numpy.eye(mm)
x_goal=numpy.ones((mmm,1))
b=-numpy.dot(g,x_goal)
c=0f = lambda x: 0.5 * (numpy.dot(numpy.dot(x.t, g), x)) + numpy.dot(b.t, x) + c
f_d = lambda x: numpy.dot(g, x) + b
x0=x_goal+numpy.random.rand(mmm,1)*100
n=100
e=10**(-6)
print '共軛梯度pr'
x1, y1, y_d1=cg_fr(x0,n,e,f,f_d)
print '共軛梯度pbr'
x2, y2, y_d2=cg_prp(x0,n,e,f,f_d)
figure1=pl.figure('trend')
n1=len(y1)
n2=len(y2)
x1=numpy.arange(1,n1+1)
x2=numpy.arange(1,n2+1)
x3, y3, y_d3=sd.sd(x0,n,e,f,f_d)
n3=len(y3)
x3=range(1,n3+1)
pl.semilogy(x3,y3,'g*',markersize=10,label='sd:'+str(n3))
pl.semilogy(x1,y1,'r*',markersize=10,label='cg-fr:'+str(n1))
pl.semilogy(x2,y2,'b*',markersize=10,label='cg-prp:'+str(n2))
pl.legend()
#影象顯示了三種不同的方法各自迭代的次數與最優值變化情況，共軛梯度方法是明顯優於最速下降法的
pl.xlabel('n')
pl.ylabel('f(x)')
pl.show()
從圖中可以看出，最速下降法sd的迭代次數是最多的，在與共軛梯度（fr與prp兩種方法）的比較中，明顯較差。
最速下降法 and 共軛梯度法
註明 程式中呼叫的函式jintuifa.m golddiv.m我在之前的筆記中已貼出 最速下降法 最速下降法求解f 1 2 x1 x1 9 2 x2 x2的最小值，起始點為x0 9 1 演算法根據最優化方法 天津大學出版社 97頁演算法3.2.1編寫 v1.0 author liuxi bit fo...
python實現梯度法 python最速下降法
更多程式設計教程請到 菜鳥教程 高州陽光論壇 人人影視 假設我們已經知道梯度法 最速下降法的原理。現給出乙個算例 如果人工直接求解 現給出python求解過程 import numpy as np from sympy import import math import matplotlib.pyp...
python實現梯度法 python最速下降法
假設我們已經知道梯度法 最速下降法的原理。現給出乙個算例 如果人工直接求解 現給出python求解過程 imwww.cppcns.comport numpy as np from sympy import import math import matpllvoifjipiotlib.pyplot a...