共轭梯度法

共轭梯度法

对于任意形式的目标函数f(X)，在极值点X附近展开成泰勒级数，且取前三项，有

T1***T2**

f(X)fXfX.XXXX.fX.XX2

*

2**

因在极值点X处fX*0，而fXH(X)为f(X)在X的二阶偏导数矩阵，

*

*



*

即Hessian矩阵，故

1*T**

f(X)fX*XX.H(X).XX 2

对于二次函数来说，若令

2fX*x1

则

2

a,

2fX*x1x2

b,

2fX*x2

2

c

ab

H(X*),而fX*d1—常数 

bc

则，得到

1

f(X)d1+x1-x1*

2



2

1

=d1+x1-x1*

2



x1-x*ab1x2-x*

2*

bcx2-x2

ax1-x*bx2-x*

12

x2-x*2bx1-x*cx2-x*

12













1

d1ax1-x1*

2



2bx1-x1*



2

*x2-x*cx-x222

x1*x1*

由上式可知，当XX时，得到目标函数的极小值

*xx22

f(X)fX*d1，当f(X)d2,d2,...时，则有等值线族。

令f(X)d2，代入上式，则有

f(X)d2d1

*

1

ax1-x1*2



2

2bx1-x1*



2

* x2-x*cx-x222



所以目标函数f(X)在X点附近的等值线方程为

*

ax1-x1



2

*

2bx1-x1



*

x2-x*cx-x222



2

d0

式中，d2(d1d2)常数。由极小值点存在的充分必要条件为二阶导数行列式小于零，

即

b2ac0

在解析几何中已经证明，当b2ac0时为椭圆，所以二元函数的等值线在极值点附近是近似的同心椭圆族。

1. 算法原理

共轭梯度法是利用目标函数梯度逐步产生共轭方向作为线搜索方向的方法，每次搜索方向都是在目标函数梯度的共轭方向，搜索步长通过一维极值算法确定。

T

设二次函数为f(X)CbX

1T

XAX，其中C为常数，b,X为n维列向量，A2

为对称正定矩阵，用共轭梯度法求f(X)的极小点：

共轭梯度法探索的第一步是沿负梯度方向。即X

(k)

点按S

(k)

fX(k)方向找到

X(k1)，然后沿着与上一次探索方向S(k)相共轭的方向S(k1)进行探索直达到最小点X*。

令S

(k1)

fX(k1)kS(k)。

(k)

S(k)的一部分即kS(k)，加上新的负梯上式的意义就是以原来的负梯度fX



(k1)

度fX，构造S(k1)。



在上式中的选择，应使n维欧氏空间E中的两个非零向量S(k)与S(k1)关于矩阵A共轭。即

(k1)(k)

SAS0

T

kn

(k0,1,2,...n1)

因

f(X)CbTX

若令

1T

XAX，故有f(X)bAX 2

g(k)fX(k)bAX(k) g(k1)fX(k1)bAX(k1)

二式相减，得

g(k1)g(k)A(X(k1)X(k))

设在第k1次迭代中

X(k1)X(k)(k)S(k)

代入上式，得

g(k1)g(k)A(k)S(k),(k0,1,2,...n1)

式中(k)为第k1次迭代的最优步长。

由式S

(k1)T

(k)

AS0

(k0,1,2,...n1)和式

g(k1)g(k)A(k)S(k),(k0,1,2,...n1)，得

(k1)(k1)

g(k))0 S(g

T

将式S

(k1)

fX(k1)kS(k)和式g(k1)fX(k1)bAX(k1)代入上式，得

[g(k1)kS(k)]T(g(k1)g(k))0

因为g(k1),g(k),...,g(0)是一正交系，故有[g(k1)]Tg(k)0或[g(k1)]TS(k)0

故上式可简化为

[g(k1)]Tg(k1)(k)[g(k)]Tg(k)0

得

(k)

g[g]g



[g(k)]Tg(k)g(k)

(k1)T

(k1)

(k1)2

2



fX(k1)fX

(k)

22



(k)用一维探索最优化方法确定，即求

minf(X(k)S(k))f(X(k)(k)S(k))

a

或用解析法，使

df(X(k)(k)S(k))

0

d

求得

(k)

。

(k1)

或由式g

g(k)A(k)S(k),(k0,1,2,...n1)，得

g(k1)g(k)A(k)S(k)

即

fX(k1)fX(k)(k)AS(k)

又由于进行一维最优化搜索，故有

fX(k1)S(k)0 

T

由上两式可求得

(k)

fX(k)S(k)

(k)T(k)SAS

T

如此，即可得到共轭梯度法的一组计算公式为

X(k1)X(k)(k)S(k)

(k)

fX(k)S(k)fX(k)S(k)

 (k)T(k)(k)T(k)(k)

SASSHXS

S(k1)fX(k1)kS(k)

TT



(k)

g[g]g

(k)T(k)

[g]gg(k)

(k1)T

(k1)

(k1)2

2



fX

(k1)(k)

fX





22

2. 算法步骤

用共轭梯度法求无约束多维极值问题minf(x),xRn的算法步骤如下： （1） 给定初始点x(0)，及精度0；

(0)(0)

（2） 若f(x)，停止，极小值点为x，否则转步骤（3）；

（3） 取p

(0)

f(x(0))，且置k0；

(k)

(k)

t0

（4） 用一维搜索法求tk，使得f(xktp)minf

k())

xtkp令，(，

x(k1)x(k)tkp(k)，转步骤5；

(k1)

)，停止，极小值点为x(k1)，否则转步骤（6）（5） 若f(x；

（6） 若k1n，令x

(0)

x(n)，转步骤（3），否则转步骤（7）；

k（7） 令p(k1)f(x(k1))kp(k)，

f(x(k1))f(x)

(k)

22

，置kk1，转步骤（4）。

3. 算法的MATLAB实现

在MATLAB中编程实现的共轭梯度法函数为：minGETD

功能：用共轭梯度法求解多维函数的极值。

调用格式：[x,minf]minGETD(f,x0,var,eps) 其中，f：目标函数；

x0：初始点； var：自变量向量；

x：目标函数取最小值时的自变量值； minf：目标函数的最小值。 共轭梯度法的MATLAB程序代码如下： function [x,minf]=minGETD(f,x0,var,eps) %目标函数:f; %初始点:x0;

%自变量向量:var;

%目标函数取最小值时的自变量值:x; %目标函数的最小值:minf; format long; if nargin==3 eps=1.0e-6; end

x0=transpose(x0); n=length(var); syms l;

gradf=jacobian(f,var); %梯度方向 v0=Funval(gradf,var,x0); p=-transpose(v0); k=0; while 1

v=Funval(gradf,var,x0); tol=norm(v); if tol

y=x0+l*p;

yf=Funval(f,var,y); [a,b]=minJT(yf,0,0.1); xm=minHJ(yf,a,b); x1=x0+xm*p;

vk=Funval(gradf,var,x1); tol=norm(vk); if tol

end

if k+1==n x0=x1; continue; else

lamda=dot(vk,vk)/dot(v,v);

p=-transpose(vk)+lamda*p; %共轭方向 k=k+1; x0=x1; end end

minf=Funval(f,var,x); format short; 例：

用共轭梯度法求函数f(t,s)(t3)2s2的极小值，其中初始值取x0(2,6)。 解：在MATLAB命令窗口中输入：

syms t s;

f=(t-3)^2+s^2;

[x,mf]=minGETD(f,[-2 6],[t s])

所得结果为： x =

3.0000 0.0000 mf =

2.0116e-037

22

例：试用共轭梯度法求解目标函数f(X)x1x2x1x210x14x260的极小值。设初

始点X

(0)

[x1(0)

(0)T

x2]00。

T

解：原函数的梯度为

f(X)[

f(X)

x1f(X)TT

]2x1x2102x2x14 x2f(X(0))104

T

第一次探索的方向：

S(0)g(0)104

T

2f(X)x2

1

2f(X)H(X)2

f(X)x2x1

T

2f(X)

x1x221



2f(X)12

2

x2

T

由式

(k)

fX(k)S(k)fX(k)S(k)

，得 (k)T(k)(k)T(k)(k)

SASSHXS

(0)

fX(0)S(0)fX(0)S(0)

(0)T(0)(0)T

HX(0)S(0)SASS

10

104

11640.763157894

2110152

104

124

TT

由此得

X

求X

(1)

X

(0)

S

(0)(0)

010T0.7631578947.631578943.052631576 04

(1)

点处的梯度为

g(1)2.210526305.52631579

T

(k)

由式

g[g]g



[g(k)]Tg(k)g(k)

(k1)T

(k1)

(k1)2

2



fX(k1)fX

(k)

22



求

(0)

：

(0)

第二次探索的方向：

g

(1)2

2

g(0)



35.42659278

0.305401661

116

S(1)g(1)(0)S(0)0.84349036.74792243

T

由式

(k)

fXS



(k)T(k)SAS

(k)

T

(k)

fXS(k)

T求探索步长为：

(k)(k)(k)

SHXS

(k)

T

0.8434903

2.210526305.52631579

6.74792243(1)

a

210.8434903

0.84349036.747922436.74792243

12

35.426592830.[1**********].10825193

由此得

X

(2)

X

(2)

aS

(1)(1)

7.99999995.9999999

T

T

x1(2)(2) x2

g(2)0.00000020.0000002

g(2)

则，极小值点为

2

x1(2)T

(2)7.99999995.9999999 x2