最优化方法及其matlab程序设计 马昌凤 课后答案

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

最优化方法-习题解答

张彦斌计算机学院2014年10月20日

Contents

1第一章最优化理论基础-P13习题1(1)、2(3)(4)、3、42第二章线搜索算法-P27习题2、4、63第三章最速下降法和牛顿法P41习题1，2，34第四章共轭梯度法P51习题1，3，6(1)5第五章拟牛顿法P73-26第六章信赖域方法P86-8

7第七章非线性最小二乘问题P98-1，2，68第八章最优性条件P112-1，2,5,6

9第九章罚函数法P132，1-(1)、2-(1)、3-(3),610第十一章二次规划习题11P178-1（1），5

14710121418232629

1第一章最优化理论基础-P13习题1(1)、2(3)(4)、3、4

1.验证下列各集合是凸集:

(1)S={(x1,x2)|2x1+x2≥1,x1 2x2≥1};需要验证：

根据凸集的定义，对任意的x(x1,x2),y(y1,y2)∈S及任意的实数λ∈[0,1],都有λx+(1 λ)y∈S.

即,(λx1+(1 λ)y1,λx2+(1 λ)y2)∈S证：由x(x1,x2),y(y1,y2)∈S得到，

{

2x1+x2≥1,x1 2x2≥1

(1)

2y1+y2≥1,y1 2y2≥1

1

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

把(1)中的两个式子对应的左右两部分分别乘以λ和1 λ,然后再相加，即得

λ(2x1+x2)+(1 λ)(2y1+y2)≥1,λ(x1 2x2)+(1 λ)(y1 2y2)≥1

合并同类项，

2(λx1+(1 λ)y1)+(λx2+(1 λ)y2)≥1,(λx1+(1 λ)y1) 2(λx2+(1 λ)y2)≥1

证毕.

2.判断下列函数为凸（凹）函数或严格凸（凹）函数：

2

(3)f(x)=x21 2x1x2+x2+2x1+3x2

(2)

(3)

首先二阶导数连续可微，根据定理1.5，f在凸集上是（I）凸函数的充分必要条件是 2f(x)对一切x为半正定；（II）严格凸函数的充分条件是 2f(x)对一切x为正定。

(

f(x)=

半正定矩阵（4）

2

2 2 22

)

(4)

41

2

f(x)= 12

30

30 4

(5)

正定矩阵

1T

3.证明f(x)=xGx+bTx为严格凸函数当且仅当Hesse矩阵G正定。

证明：根据严格凸函数定义证明。

对任意x=y,及任意实数λ∈(0,1)都有f(λx+(1 λ)y)<λf(x)+(1 λ)f(y).充分性：Hesse矩阵G正定=》严格凸函数.

TTf(λx+(1 λ)y)=1(λx+(1 λ)y)G(λx+(1 λ)y)+b(λx+(1 λ)y)1TTTTλf(x)+(1 λ)f(y)=λ(1xGx+bx)+(1 λ)(yGy+by)

1TTλf(x)+(1 λ)f(y) f(λx+(1 λ)y)=λ(1xGx)+(1 λ)(yGy) 111TTT[(λx)TG(λx)+1(1 λ)yG(1 λ)y+λxG(1 λ)y+(1 λ)yGλx]111TTTT=1λxG(1 λ)x+(1 λ)yGλy λxG(1 λ)y (1 λ)yGλx

2

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

1TT

=1λxG(1 λ)(x y)+(1 λ)yGλ(y x)1=λ(1 λ)(x y)TG(x y)>0G正定保障了严格不等式成立。反之，必要性：严格凸函数=》Hesse矩阵G正定.

类似，当对任意x=y,及任意实数λ∈(0,1)都有f(λx+(1 λ)y)<λf(x)+(1 λ)f(y).

1TT

λf(x)+(1 λ)f(y) f(λx+(1 λ)y)=λ(1xGx)+(1 λ)(yGy) 111TTT[(λx)TG(λx)+1(1 λ)yG(1 λ)y+λxG(1 λ)y+(1 λ)yGλx]111TTTT=1λxG(1 λ)x+(1 λ)yGλy λxG(1 λ)y (1 λ)yGλx1TT=1λxG(1 λ)(x y)+(1 λ)yGλ(y x)1=λ(1 λ)(x y)TG(x y)>0

4.若对任意x∈ n及实数θ>0都有f(θx)=θf(x),证明f(x)在 n上为凸函数的充要条件是 x,y∈ n,f(x+y)≤f(x)+f(y)证明：根据严格凸函数定义证明。

定义：对任意x=y,及任意实数λ∈(0,1)都有f(λx+(1 λ)y)≤λf(x)+(1 λ)f(y).

充分条件： x,y∈ n,有f(x+y)≤f(x)+f(y)

对任意x=y,及任意实数λ∈(0,1)都有f(λx+(1 λ)y)≤f(λx)+f((1 λ)y)利用f(θx)=θf(x),

f(λx+(1 λ)y)≤f(λx)+f((1 λ)y)=λf(x)+(1 λ)f(y).充分性证毕；

必要性：f(x)在 n上为凸函数=》 x,y∈ n,f(x+y)≤f(x)+f(y)根据定义有对任意x=y,及任意实数λ∈(0,1)都有f(λx+(1 λ)y)≤λf(x)+(1 λ)f(y).

不妨取λ=1,则111f(x+(1 1)y)≤f(x)+(1 )f(y).利用f(θx)=θf(x),

11f((x+y))=f(x+y)≤1(f(x)+f(y))

x,y∈ n,f(x+y)≤f(x)+f(y)证毕!

3

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

2第二章线搜索算法-P27习题2、4、6

黄金0.618算法：

function[s,phis,k,G,E]=golds(phi,a,b,delta,epsilon)%输入：phi是目标函数，a，b是搜索区间的两个端点%delta，epsilon分别是自变量和函数值的容许误差

%输出:s,phis分别是近似极小点和极大值，G是n×4矩阵，%其第k行分别是a,p,q,b的第k次迭代值ak,pk,qk,bk,%E=[ds,dphi],分别是s和phis的误差限%%%%%%%%%%t=(sqrt(5)-1)/2;h=b-a;

phia=feval(phi,a);phib=feval(phi,b);

p=a+(1-t)*h;q=a+t*h;phip=feval(phi,p);phiq=feval(phi,q);k=1;G(k,:)=[a,p,q,b];

while(abs(phib-phia)>epsilon)|(h>delta)if(phip<phiq)

b=q;phib=phiq;q=p;phiq=phip;h=b-a;p=a+(1-t)*h;phip=feval(phi,p);else

a=p;phia=phip;p=q;phip=phiq;h=b-a;q=a+t*h;phiq=feval(phi,q);end

第2题：

4

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

k=k+1;G(k,:)=[a,p,q,b];end

ds=abs(b-a);dphi=abs(phib-phia);if(phip<=phiq)s=p;phis=phip;else

s=q;phis=phiq;endE=[ds,dphi];

运行：[s,phis,k,G,E]=golds(inline(′s3 2 s+1′),0,3,0.15,0.01);结果ak,pk,qk,bk

01.14591.85413.000000.70821.14591.854100.43770.70821.14590.43770.70820.87541.14590.70820.87540.97871.14590.70820.81150.87540.97870.70820.77210.81150.87540.77210.81150.8359

0.8754

[s,phis,k,G,E]=golds(inline(′s3 2 s+1′),0,3,0.15,0.001); GG=

5

(6)

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

01.14591.854100.70821.145900.43770.70820.43770.70820.87540.70820.87540.97870.70820.81150.87540.70820.77210.81150.77210.81150.83590.77210.79650.81150.79650.81150.8208

第4题：

3.00001.85411.14591.14591.14590.97870.87540.87540.83590.8359

(7)

clearall;[s,phis,k,ds,dphi,S]=qmin(inline(′s3 2 s+1′),0,3,1e 2,1e 4); ss=0.8165第6题

functionf=fun(x)

f=100 (x(2) x(1)2)2+(1 x(1))2;functiongf=gfun(x)

gf=[ 400 (x(2) x(1)2) x(1) 2 (1 x(1)),200 (x(2) x(1)2)]′;functionmk=armijo(xk,dk)beta=0.5;sigma=0.2;m=0;mmax=20;while(m¡=mmax)

if(fun(xk+betam dk)<=fun(xk)+sigma betam gfun(xk)′ dk)mk=m;break;

6

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

endm=m+1;end

alpha=betamknewxk=xk+alpha*dkfk=fun(xk)newfk=fun(newxk)

clearall;xk=[-1,1]’;dk=[1,1]’;mk=armijo(xk,dk)alpha=0.0020newxk=-0.99801.0020fk=4newfk=3.9956mk=9

3第三章最速下降法和牛顿法P41习题1，第1题：

functionf=funone(x)

7

2，3

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

f=3 x(1)2+2 x(2)2 4 x(1) 6 x(2);functiongf=gfunone(x)gf=[6 x(1) 4,4 x(2) 6]′;

x0=[0,1]’;[xvalk]=grad(’funone’,’gfunone’,x0)x=0.66671.5000val=-5.8333k=10第2题：(1)牛顿法

functionf=funtwo1(x)

f=4 x(1)2+x(2)2 8 x(1) 4 x(2);functiongf=gfuntwo1(x)gf=[8 x(1) 8,2 x(2) 4]′;

x0=[0,1]’;[xvalk]=grad(’funtwo1’,’gfuntwo1’,x0)x=12val=-8

8

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

k=2

（2）阻尼牛顿法

functionHe=Hesstwo(x)n=length(x);He=zeros(n,n);He=[8,0;0,2];

x0=[0,1]’;[xvalk]=dampnm(’funtwo1’,’gfuntwo1’,’Hesstwo’,x0)x=12val=-8k=1第3题.

functionf=fun(x)

f=(x(1) 2)4+(x(1) 2 x(2))2;functiongf=gfun(x)

gf=[4 (x(1) 2)3+2 (x(1) 2 x(2)), 4 (x(1) 2 x(2))]′; clearall;

x0=[03]’;[v,val,k]=grad(’fun’,’gfun’,x0)

9

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

x=2.01391.0070val=3.7685e-008k=2111

4第四章共轭梯度法P51习题1，3，6(1)

1.证明向量α1=(1,0)T和α2=(3, 2)T关于矩阵

()23

A=

35

T

Aα2=0.共轭.验证α1

TT

3.设f(x)=1xHx+bx,其中

()()423

H=,b=

243

(8)

(9)

(1)证明d0=(1,0)T与d1=( 1,2)T关于H共轭;

(2)以x0=(0,0)T为初始点,d0和d1为搜索方向，用精确线搜索求f的极小点. 验证(1)dT0Hd1=0.(2)首先，g(X)= f(X)=HX+b=

()()()42x13

+

24x23

(10)

用定理4.1，也就是算法4.1产生的迭代序列，则每一步迭代点xk+1都是f(x)在x0和

∑k

方向d0,d1,...,dk所张成的线性流形，Sk={x|x=x0+i=0αidi, αi}中的极小点，特别地，xn=x = G 1b是问题的唯一极小点.

T

精确线搜索得到步长因子αk具有如下性质，gk+1dk={

+1=X k+αkdkXk

Tgk+1dk=0

(11)

10

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

+1=X k+αkdk,即X 1=X 0+α0d0;Xk

1)=g1=g(X 1)=HX 1+b,;gTd0=0g(X1 1=( 3/4,0)T,f(X 1)= 9/8;α0= 3/4,X

同理，利用(11)迭代，即

{

2=X 1+α1d1X

Tg2d1=0 2=( 1/2, 1/2)T,f(X 2)= 3/2,α1= 1/4;X 2)<f(X 1),f(X

2=( 1/2, 1/2)T定理4.1保证了极小点为X

2T

6.(1)f(x)=4x21+4x2 4x1x2 12x2,取初始点x0=( 0.5,1);g(x)= f(x)=Gx+b,G(x)= 2f(x)=G;

(12)

共轭方向的构造过程，取初始方向d0= g0,令x1=x0+α0d0,其中 f(x1)Td0=T

d0=0,在x1处，用f在x1的负梯度方向 g1与d0的组合来生成d1,即d1=g1

g1+β0d0,然后选取系数β0,使得d1与d0关于G共轭，即令dT1Gd0=0确定β0.因此，β0=

T

g1Gd0

,g001

g0=G(x1 x0)=α0Gd0,

T

di=0(i=0,1)利用定理4.1可知g2

计算过程：

()()()

8 4x10

g(x)=Gx+b=+

48x2 12

(13)

(G=(

d0= g(x0)= Gx b=

(

x1=x0+α0d0=

T

f(x1)Td0=g1d0=0

8 4

48 48)+α0

)(

)

(14)

) (=

(

)=)

(16)(

)(15)

8 4

0.5182)

0 12

82

0.51

(

8α0 0.52α0+1

[(

Tg1d0=

8

4 48

)(

8α0 0.52α0+1

11

)+

(

0 12

)]T(

82

)=0

(17)

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

α0=17/104,x1=(21/26,69/52)T≈(0.80769,1.32692)Tg1=(15/13, 60/13)Tβ0=

Tg1Gd000

=225/676≈0.33284

d1= g1+β0d0=(3315/2197,23205/4394)T≈(1.5088757,5.281065)T;

T

x2=x1+α1d1; f(x2)Td1=g2d1=0

(x2=

(g2=

(255α1)/169+21/26

(1785α1)/338+69/5215/13 (1530α1)/169(6120α1)/169 60/13

)

(18)

)

(19)

由此可以求出α1=0.127450980392157;极值点为X2=(1,2)T;

5第五章拟牛顿法P73-2

2.DFP程序算法调用

极值点x=( 0.2203×10 6, 0.1599×10 6);极小值val=1.2527×10 13

附程序：

function[x,val,k]=dfp(fun,gfun,x0)

%功能：用DFP算法求解无约束问题：minf(x)

%输入：x0是初始点，fun，gfun分别是目标函数及其梯度%输出：x,val分别是近似最优点和最优值，k是迭代次数。maxk=1e5;%给出最大迭代次数ρ=0.55;σ=0.4; =1e-5;k=0;n=length(x0);%Hk=inv(feval(’Hess’,x0));%Hk=eye(n);

12

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

Hk=[21;11];while(k<maxk)

gk=feval(gfun,x0);%计算梯度

if(norm(gk)< ),break;end%检验终止准则dk=-Hk*gk;%计算搜索方向m=0;mk=0;

while(m<20)%用Armijo搜索求步长

′

if(feval(fun,x0+ρm dk)<feval(fun,x0)+σ ρm gk dk)

mk=m;break;endm=m+1;end%DFP校正x=x0+ρmk dk;

sk=x-x0;yk=feval(gfun,x)-gk;if(sk’*yk>0)

Hk=Hk-(Hk*yk*yk’*Hk)/(yk’*Hk*yk)+(sk*sk’)/(sk’*yk);end

k=k+1;x0=x;end

val=feval(fun,x0);（I）当

H0=

(21

11

)

(20)

13

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

[x,val,k]=dfp(’fun’,’gfun’,[1,-1]’)x=1.0e-006*

-0.220306134442640-0.159928197216675val=

1.252658776679855e-013k=4

（II）当采用%Hk=inv(feval(’Hess’,x0));[x,val,k]=dfp(’fun’,’gfun’,[1,-1]’)x=00val=0k=1

6

8(1)

第六章信赖域方法P86-8

>>gk=[-6-3]’;Bk=[4-4;-48];>>dta=1;

14

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

>>[d,val,lam,k]=trustq(gk,Bk,dta)d=

0.8702817912195740.492554154744547val=

-5.928777686124834lam=

5.158202203432865k=5

>>dta=2;

[d,val,lam,k]=trustq(gk,Bk,dta)d=

1.7265693820449381.009434577568092val=

-10.321239036609670lam=

1.813689513237923k=7>>dta=5;

[d,val,lam,k]=trustq(gk,Bk,dta)d=

3.749999980155628

15

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

2.249999987787719val=

-14.624999999999998lam=

8.078453007598365e-009k=4

(2)

>>gk=[1-3-2]’;Bk=[3-12;-120;204];dta=1;

[d,val,lam,k]=trustq(gk,Bk,dta)d=

-0.2626433660099540.8374331274466090.479295543075525val=

-2.501140183861169lam=

1.268746535391740k=7>>dta=2;

[d,val,lam,k]=trustq(gk,Bk,dta)

16

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

d=

-0.3333333333333821.3333333333290360.666666666665635val=

-2.833333333333333lam=

6.261736529506079e-012k=5>>dta=5;

[d,val,lam,k]=trustq(gk,Bk,dta)d=

-0.3333333333334331.3333333331807220.666666666628600val=

-2.833333333333333lam=

2.286320834416492e-010k=4

17

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

7第七章非线性最小二乘问题P98-1，2，6

2

f1(x)=x31 2x2 1=0f2(x)=2x1+x2 2=0

1.设有非线性方程组

(21)

(1)列出求解这个方程组的非线性最小二乘问题的数学模型；最小二乘问题的数学表达式：minx∈Rnf(x)=

1∥F(x)∥=

1∑m

i=1

fi2(x)

(2)写出求解该问题的高斯-牛顿法迭代公式的具体形式：

(

)

(22)

Jk=F′(x(k))=( F1(x(k)),···, Fm(x(k)))T=

TT

dGN= [JkJk] 1JkF(xk)=k

3x21,k

2 4x2,k

1

[(

3x21,k 4x2,k

21

)(

3x21,k2 4x2,k

1

)] 1(

3x21,k 4x2,k

21

)(

)2

x3 2x 11,k2,k

2x1,k+x2,k 2

(23)

(3)初始点取为x0=(2,2)T,迭代三次：迭代公式：Xk+1=Xk+dGNkX1=X0+dGN0=3.1071428571428593.785714285714287X2=X1+dGN1=5.1574316407151187.685136718569831X3=X2+dGN2=8.766682264589718

18

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

16.466635470820520

2.

解答：（1）测得的t1,t2和y一共5组数据，分别代入关系式

y=

x1x3t1

1+x1t1+x2t2

(24)

F1(x)=x1x3 0.13(1+x1+x2) F2(x)=2x1x3 0.22(1+2x1+x2)F3(x)=x1x3 0.08(1+x1+2x2) F4(x)=2x1x3 0.13(1+2x1+2x2)

F5(x)=0.1x1x3 0.19(1+0.1x1)

(1)最小二乘问题模型表示为minx∈Rnf(x)=(2)高斯牛顿迭代公式的具体公式为：

TT

dGN= [JkJk] 1JkF(xk)k

1

1x3

0.13=x 12 2x1x3 0.22= 12

x0.08=x12

2x1x3 0.13= 12 0.1xx0.19=1

(25)

(26)

∑m

i=1

∥F(x)∥=

1Fi2(x)

∑5

2

6.利用LM方法的matlab程序求解minf(x)=1i=1ri(x)其中 22

r1(x)=x2 1+x2+x3 1 r2(x)=x1+x2+x3 1

22

r3(x)=x21+x2+(x3 2) 1 r4(x)=x1+x2 x3+1 22r5(x)=x31+3x2+(5x3 x1+1) 36t

(27)

t为参数，可取t=0.5,1,5等，注意当t=1时，x =(0,0,1)T是全局极小

点，这时问题为零残量，比较不同参数的计算效果。function[x,val,k]=lmm(Fk,JFk,x0)%功能:用L-M方法求解非线性方程组:F(x)=0

%输入:x0是初始点,Fk,JFk分别是求F(xk)及F’(xk)的函数%输出:x,val分别是近似解及——F(xk)——的值,k是迭代次数.maxk=1000;%给出最大迭代次数

19

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

ρ=0.55;σ=0.4;µk=norm(feval(Fk,x0));k=0;epsilon=1e-6;n=length(x0);while(k<maxk)

fk=feval(Fk,x0);%计算函数值jfk=feval(JFk,x0);%计算Jacobi阵gk=jfk’*fk;

dk= (jfk′ jfk+µk eye(n))\gk;%解方程组Gk*dk=-gk,计算搜索方向if(norm(gk)¡epsilon),break;end%检验终止准则m=0;mk=0;

while(m<20)%用Armijo搜索求步长

newf=0.5 norm(feval(Fk,x0+ρm dk))2;oldf=0.5 norm(feval(Fk,x0))2;if(newf<oldf+sigma ρm gk′ dk)mk=m;break;endm=m+1;end

x0=x0+ρmk dk;muk=norm(feval(Fk,x0));k=k+1;endx=x0;val=0.5*µ2k;

20

最优化方法及其matlab程序设计 马昌凤版 课后答案 杭电课件

%gval=norm(gfun(x));%%目标函数(I)t=0.5

functiony=Fk(x)

y(1)=x(1)2+x(2)2+x(3)2 1;y(2)=x(1)+x(2)+x(3) 1;y(3)=x(1)2+x(2)2+(x(3) 2)2 1;y(4)=x(1)+x(2) x(3)+1;

y(5)=x(1)3+3 x(2)2+(5 x(3) x(1)+1)2 36 0.5;y=y(:);

%%%%Jacobi阵%%%%%%functionJF=JFk(x)JF=[2*x(1),2*x(2),2*x(3);1,1,1;

2*x(1),2*x(2),2*(x(3)-2);1,1-1;

3 x(1)2 2 (5 x(3) x(1)+1),6 x(2),10 (5 x(3) x(1)+1)];>>x0=[1,1,1]′;[x,val,k]=lmm(′Fk′,′JFk′,x0)x=

0.339361063668441-0.2001835788046710.714384339944574val=

0.486062168183995

21

本文来源：https://www.bwwdw.com/article/5u01.html