《数值最优化方法》课程教学课件（讲稿，打印版）拟牛顿法

最优化方法及其Matlab程序设计 第五章拟牛顿法 Back Close

1/52 JJ II J I Back Close Å`zê{9Ÿ Matlab ßSO 1 Ÿ [⁄Ó{

第3章所介绍的牛顿法的优点是具有二阶收敛速度，但当Hesse 阵G(x)=Vf(x)不正定时，不能保证所产生的方向是目标函数在 xk处的下降方向.特别地，当G(x)奇异时，算法就无法继续进行下 去.尽管修正牛顿法可以克服这一缺陷，但其中的修正参数的选取 很难把握，过大或过小都会影响到收敛速度.此外，牛顿法的每一迭代 步都需要目标函数的二阶导数，即Hesee阵，对于大规模问题其计算 量是惊人的 本章即将介绍的拟牛顿法克服了这些缺点，并且在一定条件下这 类算法仍然具有较快的收敛速度一超线性收敛速度， S5.1拟牛顿法及其性质 拟牛顿法的基本思想是在基本牛顿法的步2中用Hesee阵G= Vf(xk)的某个近似矩阵Bk取代Gk.通常，Bk应具有下面的三个特 Back Close

2/52 JJ II J I Back Close 1 3 Ÿ§0⁄Ó{`:¥‰k¬ÒÑ›,  Hesse G(xk) = ∇2 f(xk) ÿ½û, ÿUy§)êï¥8IºÍ3 xk ?e¸êï. AO/,  G(xk) ¤.û, é{“Ã{UY?1e . ¶+?⁄Ó{å±é—˘ò"Ä, Ÿ•?ÎÍ µk ¿ ÈJrº, Lå½L—¨Kè¬ÒÑ›. d , ⁄Ó{zòSì ⁄—Iá8IºÍÍ, = Hesee , Èuå5ØKŸOé ˛¥Ø<. Ÿ=Ú0[⁄Ó{é— ˘ ":, øÖ3ò½^áe˘ aé{E,‰kجÒÑ›—áÇ5¬ÒÑ›. §5.1 [⁄Ó{9Ÿ5ü [⁄Ó{ƒgé¥3ƒ⁄Ó{⁄ 2 •^ Hesee Gk = ∇2 f(xk) ,áCq› Bk ì Gk. œ~, Bk A‰ke°náA

点 (1)在某种意义下有B≈Gk,使相应的算法产生的方向近似于 牛顿方向，以确保算法具有较快的收敛速度 (2)对所有的k,B是对称正定的，从而使得算法所产生的方向 是函数∫在xk处下降方向 (3)矩阵B更新规则相对比较简单，即通常采用一个秩1或秩 2矩阵进行校正 下面介绍满足这三个特点的矩阵Bk的构造.设∫：R”→R在开 集DCR”上二次连续可微.那么，∫在xk+1处的二次近似模型为 fa)≈f(+1)+9眼+1(c-xk+H)+)-k+1)PGk+1(-k+1). 对上式求导数得 g(c)≈9k+1+Gk+1(x-ck+1): Back Close

3/52 JJ II J I Back Close :: (1) 3,´ø¬ek Bk ≈ Gk, ¶ÉAé{)êïCqu ⁄Óêï, ±(é{‰kجÒÑ›. (2) ȧk k, Bk ¥È°½, l ¶é{§)êï ¥ºÍ f 3 xk ?e¸êï. (3) › Bk ç#5KÉÈ'{¸, =œ~Ê^òáù 1 ½ù 2 › ?1. e°0˜v˘náA:› Bk E.  f : R n → R 3m 8 D ⊂ R n ˛gÎYåá. @o, f 3 xk+1 ?gCq.è f(x) ≈ f(xk+1) + g T k+1(x − xk+1) + 1 2 (x − xk+1) TGk+1(x − xk+1). È˛™¶Í g(x) ≈ gk+1 + Gk+1(x − xk+1).

令x=xk,位移Sk=xk+1一xk,梯度差y=g张+1一9k,则有 G+15k≈h. 注意到，对于二次函数∫，上式是精确成立的.现在，我们要求在拟牛 顿法中构造出Hesse阵的近似矩阵Bk满足这种关系式，即 Bk+1sk=Yk (5.1) 上式通常称作拟牛顿方程或拟牛顿条件.令Hk+1=B1,则得到拟 牛顿方程的另一个形式： Hk+1yk=Sk (5.2) 其中Hk+1是Hesse阵逆的近似.搜索方向由dk=一Hkgk或Bkdk= 一9k确定.根据Bk(或H)的第三个特点，可令 Bk+1=Bk+Ek:Hk+1=Hk+Dk; (5.3) Back Close

4/52 JJ II J I Back Close - x = xk, †£ sk = xk+1 − xk, F› yk = gk+1 − gk, Kk Gk+1sk ≈ yk. 5ø, ÈugºÍ f, ˛™¥°(§·. y3, ·Çá¶3[⁄ Ó{•E— Hesse Cq› Bk ˜v˘´'X™, = Bk+1sk = yk. (5.1) ˛™œ~°ä[⁄Óêß½[⁄Ó^á. - Hk+1 = B −1 k+1, K[ ⁄Óêß,òá/™: Hk+1yk = sk, (5.2) Ÿ• Hk+1 ¥ Hesse _Cq. |¢êïd dk = −Hkgk ½ Bkdk = −gk (½. ä‚ Bk(½ Hk) 1náA:, å- Bk+1 = Bk + Ek, Hk+1 = Hk + Dk, (5.3)

其中E,Dk是秩1或秩2矩阵.通常将由拟牛顿方程(5.1)（或(5.2）) 和校正规则(5.3)所确立的方法称为拟牛顿法， 下面我们介绍一个对称秩1校正公式.在(5.3)中取E=au以（秩 1矩阵)，其中a∈R,uk∈Rm.由拟牛顿方程(5.1)得 (Bi+aukur)5k=Vk, 即有 a(uk sk)uk =yk -BiSk. (5.4) 上式表明向量平行于(-Bksk),即存在常数B使得=B(y一 Bks).因此有 Ek:aB2(Uk:-BkSk)(Uk:-Bk.Sk)T. 于是，由(5.4)得 aB2[(V1:-Bi.Sk)Tsk](Uk:-BiSk)=(Uk:-BkSk). Back Close

5/52 JJ II J I Back Close Ÿ• Ek, Dk ¥ù 1 ½ù 2 › . œ~Úd[⁄Óêß (5.1) (½ (5.2)) ⁄5K (5.3) §(·ê{°è[⁄Ó{. e°·Ç0òáÈ°ù 1 ˙™. 3 (5.3) • Ek = αuku T k (ù 1 › ), Ÿ• α ∈ R, uk ∈ R n . d[⁄Óêß (5.1)  (Bk + αuku T k )sk = yk, =k α(u T k sk)uk = yk − Bksk. (5.4) ˛™L²ï˛ uk ²1u (yk − Bksk), =3~Í β ¶ uk = β(yk − Bksk). œdk Ek = αβ2 (yk − Bksk)(yk − Bksk) T . u¥, d (5.4)  αβ2 [(yk − Bksk) T sk](yk − Bksk) = (yk − Bksk)

由此，若(-Bxs)Tsk≠0，可取a2[(-Bs)Ts=1,即 a32= 1 E=-BS跳-BS)Y (Uk -BkSk)TSk (Uk -BkSk)TSk 故得对称秩1校正公式如下： B+1=B6+-BAs跳-BS)T (yk-BiSk)TSk (5.5) 类似地，利用拟牛顿方程(5.2)，对H进行对称秩1修正可得 He1=Hg+Sk-HEyk)(Sk-HKVR)T (5.6) (Sk:-HkVK)TVk 有了对称秩1校正公式后，利用它可以构造求解无约束优化问题 的一个拟牛顿算法，步骤如下： 算法5.1（对称秩1算法） 步0给定初始点x0∈R”,终止误差0≤e←1.初始对称正定 阵Ho(通常取单位阵In).令k:=0. Back Close

6/52 JJ II J I Back Close dd, e (yk − Bksk) T sk 6= 0, å αβ2 [(yk − Bksk) T sk] = 1, = αβ2 = 1 (yk − Bksk) T sk , Ek = (yk − Bksk)(yk − Bksk) T (yk − Bksk) T sk . È°ù 1 ˙™Xe: Bk+1 = Bk + (yk − Bksk)(yk − Bksk) T (yk − Bksk) T sk . (5.5) aq/, |^[⁄Óêß (5.2), È Hk ?1È°ù 1 ?å Hk+1 = Hk + (sk − Hkyk)(sk − Hkyk) T (sk − Hkyk) T yk . (5.6) k È°ù 1 ˙™￾, |^ßå±E¶)ÃÂ`zØK òá[⁄Óé{, ⁄½Xe: é{ 5.1 (È°ù 1 é{) ⁄ 0 â½–©: x0 ∈ R n , ™éÿ 0 ≤ ε  1. –©È°½ H0 (œ~¸† In). - k := 0

步1若‖g‖≤e,停算，输出xk作为近似极小点. 步2计算搜索方向dk=-Hgk 步3用线搜索技术求步长Q. 步4令xk+1=xk+adk,由对称秩1校正公式(5.6)确定H+1 步5令k=k十1，转步1. 下面给出基于Armijo搜索的对称秩1算法的Matlab程序， 程序5.1（对称秩1算法程序) function [x,val,k]=sr1(fun,gfun,x0) %功能：用对称秩1算法求解无约束问题：minf(x) %输入：x0是初始点，fun,gfun分别是目标函数及其梯度 %输出：x,val分别是近似最优点和最优值，k是迭代次数. maxk=500;%给出最大迭代次数 rho=0.55;sigma=0.4;epsilon=1e-5; k=0;n=length(x0);Hk=eye(n); while(k<maxk) Back Close

7/52 JJ II J I Back Close ⁄ 1 e kgkk ≤ ε, é, —— xk äèCq4:. ⁄ 2 Oé|¢êï dk = −Hkgk. ⁄ 3 ^Ç|¢E‚¶⁄ αk. ⁄ 4 - xk+1 = xk +αkdk, dÈ°ù 1 ˙™ (5.6) (½ Hk+1. ⁄ 5 - k := k + 1, =⁄ 1. e°â—ƒu Armijo |¢È°ù 1 é{ Matlab ßS. ßS 5.1 ( È°ù 1 é{ßS) function [x,val,k]=sr1(fun,gfun, x0) %ıU: ^È°ù1é{¶)ÃÂØK: min f(x) %—\: x0¥–©:, fun, gfun©O¥8IºÍ9ŸF› %——: x, val©O¥CqÅ`:⁄Å`ä, k¥SìgÍ. maxk=500; %â—ÅåSìgÍ rho=0.55;sigma=0.4; epsilon=1e-5; k=0; n=length(x0); Hk=eye(n); while(k<maxk)

gk=feval(gfun,x0);%计算梯度 dk=-k*gk;%计算搜索方向 if(norm(gk)<epsilon),break;end%检验终止准则 m=0;mk=0; while(m<20)%用Armijo搜索求步长 if(feval(fun,x0+rho m*dk)<feval(fun,x0)+sigma*rho m*gk'*dk) mk=m;break; end m=m+1; end x=x0+rho mk*dk; sk=x-x0;yk=feval(gfun,x)-gk; Hk=Hk+(sk-Hk*yk)*(sk-Hk*yk)'/(sk-Hk*yk)'*yk);%秩1校正 k=k+1; x0=X; end val=feval(fun,x0); 例5.1利用程序5.1求解无约束优化问题 minf(x)=100(x-x2)2+(1-1)2,x=(c1,2)T∈R2. Back Close

8/52 JJ II J I Back Close gk=feval(gfun,x0); %OéF› dk=-Hk*gk; %Oé|¢êï if(norm(gk)<epsilon), break; end %u™éOK m=0; mk=0; while(m<20) % ^Armijo|¢¶⁄ if(feval(fun,x0+rho^m*dk)<feval(fun,x0)+sigma*rho^m*gk’*dk) mk=m; break; end m=m+1; end x=x0+rho^mk*dk; sk=x-x0; yk=feval(gfun,x)-gk; Hk=Hk+(sk-Hk*yk)*(sk-Hk*yk)’/((sk-Hk*yk)’*yk); %ù1 k=k+1; x0=x; end val=feval(fun,x0); ~ 5.1 |^ßS 5.1 ¶)ÃÂ`zØK min f(x) = 100(x 2 1 − x2) 2 + (x1 − 1)2 , x = (x1, x2) T ∈ R 2

该问题有精确解x*=(1,1)T,fx*)=0. 解取终止准则值为‖Vf(xk)川≤10-5，利用程序5.1，取不同的初 始点，数值结果如下表 表5.1对称秩1校正算法的数值结果 初始点(0)迭代次数（目标函数值(f(c) (0,0)T 22 7.0304×10-19 0.5,0.5)7 19 3.8208×10-16 (2,2)7 38 3.3992×10-20 (-1,-1)7 45 8.2927×10-16 (1,10)7 98 1.9321×10-16 (10,10)7 142 2.1578×10-15 说明上述程序的调用方式是： x0=[-1.21]; [x,val,k]=sr1('fun','gfun',x0) 其中fun,gfun分别是求目标函数值及其梯度的M函数文件 Back Close

9/52 JJ II J I Back Close TØKk°() x ∗ = (1, 1)T , f(x ∗ ) = 0. ) ™éOKäè k∇f(xk)k ≤ 10−5 , |^ßS 5.1, ÿ”– ©:, Íä(JXeL. L 5.1 È°ù 1 é{Íä(J. –©: (x0) SìgÍ (k) 8IºÍä (f(xk)) (0, 0)T 22 7.0304 × 10−19 (0.5, 0.5)T 19 3.8208 × 10−16 (2, 2)T 38 3.3992 × 10−20 (−1, −1)T 45 8.2927 × 10−16 (1, 10)T 98 1.9321 × 10−16 (10, 10)T 142 2.1578 × 10−15 `² ˛„ßSN^ꙥ: x0=[-1.2 1]’; [x,val,k]=sr1(’fun’,’gfun’,x0) Ÿ• fun, gfun ©O¥¶8IºÍä9ŸF› M ºÍ©á

§5.2BFGS算法及其Matlab实现 BFGS校正是目前最流行也是最有效的拟牛顿校正，它是由Broy- den,Fletcher,Goldfarb和Shanno在1970年各自独立提出的拟牛顿法， 故称为BFGS算法.其基本思想是：在(5.3)中取修正矩阵E为秩2 矩阵： Er =aukug Bvkvt, 其中u,v∈R”是待定向量，Q,B∈R是待定实数.于是由拟牛顿方 程(5.1)可得 (Bk+au4f+Bku)sk=张， 或等价地 a(u以sJu+B(fs)v%=张-Bksk (5.7） Back Close

10/52 JJ II J I Back Close §5.2 BFGS é{9Ÿ Matlab ¢y BFGS ¥8cÅ61è¥Åk[⁄Ó, ߥd Broy￾den, Fletcher, Goldfarb ⁄ Shanno 3 1970 càg’·J—[⁄Ó{, °è BFGS é{. Ÿƒgé¥: 3 (5.3) •?› Ek èù 2 › : Ek = αuku T k + βvkv T k , Ÿ• uk, vk ∈ R n ¥ñ½ï˛, α, β ∈ R ¥ñ½¢Í. u¥d[⁄Óê ß (5.1) å (Bk + αuku T k + βvkv T k )sk = yk, ½d/ α(u T k sk)uk + β(v T k sk)vk = yk − Bksk. (5.7)