a03共轭梯度法编程.docx
《a03共轭梯度法编程.docx》由会员分享，可在线阅读，更多相关《a03共轭梯度法编程.docx（16页珍藏版）》请在冰豆网上搜索。
a03共轭梯度法编程
MK7FCCK282Z3SV
MK6CZBW26Q6LZ9MKHMNXTJJQ4S83
MKMDE2RESG33HS
P90第二章11
（2）
用大M法求解
minw=2X1+X2-X3-X4s.tx1-x2+2x3-x4=2
2x1+x2-3x3+x4=6x1+x2+x3+x4=7
Xi>0,i=1,2,3,4
用matlab求解如下：
f=[2,1,-1,-1,200,200,200];
a=[1-12-1100;21-31010;1111001];b=[267];
lb=zeros（7,1）;[X,fval,eXitflag,output,lambda]=linprog（f,[],[],a,b,lb）Optimizationterminated.
运行结果如下：
3.0000
0.0000
1.0000
3.0000
0.0000
0.0000
0.0000fval=
2.0000exitflag=
1output=
iterations:
7
algorithm:
'large-scale:
interiorpoint'cgiterations:
0
message:
'Optimizationterminated.'lambda=
ineqlin:
[0X1double]eqlin:
[3X1double]upper:
[7X1double]lower:
[7X1double]
从上述运行结果可以得出：
最优解为X=，最小值约为f*=2。
P151
第三章26
用共轭梯度算法求f（x）=（x1-1）A2+5*（x2-x1A2）A2的极小点，取初始点xO=。
26题
用matlab求解如下：
functionmg=MG（）%%共轭梯度法求解习题三第%
clc;clear;
n=2;
x=[2O]';
max_k=1OO;count_k=1;
trace（1,1）=x
（1）;
trace（2,1）=x
（2）;trace（3,1）=f_fun（x）;k=O;
g1=f_dfun（x）;
s=-g1;
whilecount_k<=max_kifk==n
g0=f_dfun（x）;
s=-g0;
k=0;
else
r_min=fminbnd（@（t）f_fun（x+t*s）,-100,100）;
x=x+r_min*s;
g0=g1;
g1=f_dfun（x）;
ifnorm（g1）<10A（-6）
break;
end
m=（norm（g1）A2）/（norm（gO）A2）;
s=-g1+m*s;
count_k=count_k+1;
trace（1,count_k）=x
（1）;
trace（2,count_k）=x
（2）;
trace（3,count_k）=f_fun（x）;
k=k+1;
endendcount_kf=f_fun（x）functiong=f_dfun（x）
g（1,1）=20*x
（1）A3-20*x
（1）*x
（2）+2*x
（1）-2;
g（2,1）=10*x
（2）-10*x
（1）A2;
functionf=f_fun（x）
f=5*x
（1）A4-10*x
（1）A2*x
（2）+x
（1）A2-2*x
（1）+1+5*x
（2）人2;运行结果如下：
k=
59
x0=
1.0000
1.0000
f0=
0
从上述运行结果可以得出：
最优解为x=，最小值约为f*=0。
P151第三章26
用BFGS算法求f（x）=（x1-1）^2+5*（x2-x1^2）^2的极小点，取初始点x0=。
用matlab求解如下：
symsx1x2;
f=（x1-1）A2+5*（x2-x1A2）A2;
v=[x1,x2];df=jacobian（f,v）;
df=df.';x0=[2,0]';g1=subs（df,{x1,x2},{x0（1,1）,x0（2,1）}）;k=0;
H0=[1,0;0,1];while（norm（g1）>0.001&k<300）
ifk==0p=-H0*g1;
elsevk=sk/（sk'*yk）-（H0*yk）/（yk'*H0*yk）;w1=（yk'*H0*yk）*vk*vk';
H1=H0-（H0*yk*yk'*H0）/（yk'*H0*yk）+（sk*sk'）/（sk'*yk）+w1;
p=-H1*g1;
H0=H1;
end
x00=x0;
result=Usearch1（f,x1,x2,df,x0,p）;
arf=result
（1）;
x0=x0+arf*p;
g0=g1;
g1=subs（df,{x1,x2},{x0（1,1）,x0（2,1）}）;p0=p;
yk=g1-g0;sk=x0-x00;
k=k+1;
end;
k
x0f0=subs（f,{x1,x2},{x0（1,1）,x0（2,1）}）
functionresult=Usearch1（f,x1,x2,df,x0,p）
mu=0.001;sgma=0.99;a=0;b=inf;arf=1;
pk=p;
x3=x0;
x4=x3+arf*pk;
f1=subs（f,{x1,x2},{x3（1,1）,x3（2,1）}）;
f2=subs（f,{x1,x2},{x4（1,1）,x4（2,1）}）;
gk1=subs（df,{x1,x2},{x3（1,1）,x3（2,1）}）;
gk2=subs（df,{x1,x2},{x4（1,1）,x4（2,1）}）;
while（f1-f2<=-mu*arf*gk1'*pk）b=arf;arf=（a+arf）/2;x4=x3+arf*pk;f2=subs（f,{x1,x2},{x4（1,1）,x4（2,1）}）;gk2=subs（df,{x1,x2},{x4（1,1）,x4（2,1）}）;
end;
while（1>0）
if（gk2'*pkwhile（f1-f2<=-mu*arf*gk1'*pk）b=arf;arf=（a+arf）/2;x4=x3+arf*pk;f2=subs（f,{x1,x2},{x4（1,1）,x4（2,1）}）;gk2=subs（df,{x1,x2},{x4（1,1）,x4（2,1）}）;
end;
else
break;
endend;result=[arf];运行结果如下：
k=
19x0=
1.0000
1.0000
f0=4.4505e-013
从上述运行结果可以得出：
最优解为x=，最小值约为f*=0。
P229第四章10
用SUM-内点法求解
（1）minf（x）=x1+x2
s.tg1（x）=-x1A2+x2>0
g2（x）=x1>0
用matlab求解如下：
functionsumt=SUMT（x0,e0,max_k0）
%SUM内点法globalxsMx=x0;
e=e0;max_k=max_k0;
trace（1,1）=x
（1）;trace（2,1）=x
（2）;M=100;
c=3;
e_FR=10A-10;
max_FR=200;
fork=0:
max_k
x=FR（x,e_FR,max_FR）;
trace（1,k+2）=x
（1）;
trace（2,k+2）=x
（2）;
iff_pfun（x）break;
end
M=c*M;
end
xf=f_fun（x）ktrace
functionf=FR（xO,e,max_k）
globalxs;
count_k=1;
k=O;
x=xO;
g1=f_dfun（x）;
s=-g1;
whilecount_k<=max_k
ifk==n
g0=f_dfun（x）;s=-g0;
k=0;
else
r_min=fminbnd（@（t）f_fun（x+t*s）,-100,100）;x=x+r_min*s;%
g0=g1;g1=f_dfun（x）ifnorm（g1）end
m=（norm（g1）A2）/（norm（gO）A2）;s=-g1+m*s;
count_k=count_k+1;k=k+1;
end
endcount_k=count_k-1f=x
注意：
以下三个函数体在运行计算习题四第8
（1）题时使用%*****************************************************functionf=f_fun（x）
globalM;
f=x
（1）+x
（2）+M*f_pfun（x）;%题8
（1）
functiong=f_dfun（x）
globalM;
a=-x
（1）A2+x
（2）;
b=x
（1）;
ifa>=O&b>=O
g（1,1）=1;
g（2,1）=1;elseifa<0&b>=0
g（1,1）=1+2*M*（-x
（1）A2+x
（2））*（-2*x
（1））;
g（2,1）=1+2*M*（-x
（1）A2+x
（2））;
elseifa>=0&b<0
g（1,1）=1+2*M*x
（1）;
g（2,1）=1;
else
g（1,1）=1+2*M*（-x
（1）A2+x
（2））*（-2*x
（1））+2*M*x
（1）;g（2,1）=1+2*M*（-x
（1）A2+x
（2））;
end懑8
（1）functionp=f_pfun（x）
%罚函数
a=-x
（1）A2+x
（2）;b=x
（1）;
ifa>=0&b>=0
p=0;
elseifa>=0&b<0
p=b^2;
elseifa<0&b>=0
p=a^2;
else
p=a^2+b^2;
end懑8
（1）运行结果如下：
1.0e-004*
-0.2058
-0.2058
f=-2.0576e-005trace=
100.0000-0.0050-0.0017-0.0006-0.0002-0.0001-0.00003.0000-0.0050-0.0017-0.0006-0.0002-0.0001-0.0000
从上述运行结果可以得出：
最优解为x=，最小值约为f*=0。
（2）minf（x）=x1^2+x2^2
S.t2x1+x2-2w0
-x2+1w0
用matlab求解如下：
主程序与
（1）相同；注意：
以下三个函数体在运行计算习题四第8
（2）题时使用%*****************************************************
functiong=f_dfun（x）globalM;
a=-2*x
（1）-x
（2）+2;b=x
（2）-1;ifa>=0&b>=0g（1,1）=2*x
（1）;g（2,1）=2*x
（2）;
elseifa<0&b>=0
g（1,1）=2*x
（1）+2*M*（-2*x
（1）-x
（2）+2）*（-2）;
g（2,1）=2*x
（2）+2*M*（-2*x
（1）-x
（2）+2）*（-1）;elseifa>=0&b<0
g（1,1）=2*x
（1）;
g（2,1）=2*x
（2）+2*M*（x
（2）-1）;
else
g（1,1）=2*x
（1）+2*M*（-2*x
（1）-x
（2）+2）*（-2）;
g（2,1）=2*x
（2）+2*M*（-2*x
（1）-x
（2）+2）*（-1）+2*M*（x
（2）-1）;
end懑8
（2）%*****************************************************
functionf=f_fun（x）
globalM;
f=x
（1）A2+x
（2）A2+M*1_pfun（x）;%题8
（2）
functionp=f_pfun（x）
a=-2*x
（1）-x
（2）+2;
b=x
（2）-1;
ifa>=0&b>=0
p=0;
elseifa>=0&b<0
p=b^2;
elseifa<0&b>=0
p=a^2;
else
p=a^2+b^2;
end懑8
（2）运行结果如下：
x=
-0.0000
1.0000
1.0000
trace=
100.0000
-0.0000
-0.0000
-0.0000
-0.0000
-0.0000
-0.0000
0.99991.00001.0000
-0.00003.00000.99010.99670.99890.9996
从上述运行结果可以得出：
最优解为x=
，最小值约为f*=1。
P232第四章25
用梯度投影法求解下列线性约束优化问题
minf（x）=x1A2+x1*x2+2*x2A2+2*x3A2+2*x2*x3+4*x1+6*x2+12*x3
s.tx1+x2+x3w6
-x1-x2+2x3>2xi>0,i=1,2,3
取x1=,用matlab求解如下：
f=‘x1^2+x1*x2+2*x2^2+2*x3^2+2*x2*x3+4*x1+6*x2+12*x3a=[111;11-2];
b=[6;-2];
l=zeros（3,1）;
x0=[113];
[x,fval,exitflag,output,lambda,grad,hessian]=fmincon（f,x0,a,b,[],[],l,[]）
运行结果如下：
fval=
14exitflag=
1output=
iterations:
2
funcCount:
14
stepsize:
1
algorithm:
'medium-scale:
SQP,Quasi-Newton,line-search'firstorderopt:
0
cgiterations:
[]
message:
[1x144char]lambda=
lower:
[3x1double]
upper:
[3x1double]eqlin:
[0x1double]eqnonlin:
[0x1double]ineqlin:
[2x1double]ineqnonlin:
[0x1double]grad=
4.0000
8.0000
16.0000
hessian=
1.11460.67710.6042
0.67713.36462.4792
0.60422.47923.4583
从上述运行结果可以得出：
最优解为x=，最小值约为f*=14。
P234第四章34
（1）
用乘子法求解
maxf（x）=10x1+4.4x2^2+2x3
s.tx1+4x2+5x3
x1+3x2+2x3
x3^2/2+x2^2
x1>2,x2
用matlab求解如下：
function[x,minf]=ymh434（l）formatlong;
symsx1x2x3
f=10*（x1）+4.4*（x2）A2+2*（x3）;
h=[-x1-4*x2-5*x3+32,-x1-3*x2-2*x3+29,（（x3）A2）/2+（x2）A2-3,x1-2,x2,x3];
x0=[220];v=[100000];
M=2;alpha=2;gama=0.25;var=[x1x2x3];eps=1.0e-4;ifnargin==8eps=1.0e-4;
endm1=transpose（x0）;m2=inf;
whilel
FE=0;
u=subs（h,{x1,x2,x3},[m1
（1）,m1
（2）,m1（3）]）;
fori=1:
length（h）
if（v（i）+M*u（i））<=0
FE=FE+（（v（i）+M*u（i））A2-（v（i））A2）;
else
FE=FE+（v（i））A2;
end
end
SumF=f+（1/（2*M））*FE;
[m2,minf]=minNT（SumF,transpose（m1）,var,eps）;
Hm2=subs（h,{x1,x2,x3},[m2
（1）,m2
（2）,m2（3）]）;
Hm1=subs（h,{x1,x2,x3},[m1
（1）,m1
（2）,m1（3）]）;
Hx2=Funval（h,var,x2）;
Hx1=Funval（h,var,x1）;
ifnorm（Hx2）x=x2;
break;
elseifHx2/Hx1>=gama
M=alpha*M;x1=x2;
else
v=v-M*transpose（Hx2）;x1=x2;
end
endend
minf=Funval（f,var,x）;formatshort;
运行结果如下：
x=
19.7502.45minf=
-202.42
maxf1=
202.42
从上述运行结果可以得出：
最优解为x=，最大值约为f*=202.42,最小值约为f*=-202.42。
P235第四章35
（2）
用序列二次规划法求解
minf（x）=-5x1-5x2-4x3-x1x3-6x4-5x5/（1+x5）-8x6心+x6）-10（1-2eA-x7+eA-2x7）s.tg1（x）=2x4+x5+0.8x6+x7-5=0
g2（x）=x2A2+x3A2+x5A2+x6A2-5=0
g3（x）=x1+x2+x3+x4+x5+x6+x7
g4（x）=x1+x2+x3+x4w5
g5（x）=x1+x3+x5+x6^2-x7^2-5
xi>0,i=1,2,…,7
用matlab求解如下：
functionf=objfun35（x）
f=-5*x
（1）-5*x
（2）-4*x（3）-x
（1）*x（3）-6*x（4）-5*x（5）/（1+x（5））-8*x（6）/（1+x（6））-10*（1-2*exp（-x（7））+exp（-2*x（7）））;
function[c,ceq]=confun35（x）
ceq=[x
（2）^2+x（3）^2+x（5）^2+x（6）人2-5];
c=[x
（1）+x（3）+x（5）+x（6）A2-x（7）A2-5];
clcclearx0=[1,1,1,1,1,1,1];
aeq=[0,0,0,2,1,0.8,1];
beq=5;
a=[1,1,1,1,1,1,1;1,1,1,1,0,0,0];
b=[10,5]';
lb=[0,0,0,0,0,0,0];
ub=[];
[x,fval,exitflag,output,lambda,grad,hessian]=fmincon（'objfun35',x0,a,b,aeq,beq,lb,ub,'confun35',options）[c,ceq]=confun35（x）
运行结果如下：
derivativeoptimality
IterF-countf（x）constraintStep-size
Procedure
0
8
-31.4958
Infeasiblestartpoint
117-41.7175
0.8642
1
-6.3
1.62
2
26
-43.0309
0.5358
1
-1.09
0.924
3
35
-44.5199
0.8348
1
-1.4
0.438
4
44
-44.4543
0.04621
1
0.126
0.247
5
53
-44.4636
0.004354
1
-0.00456
0.127
6
62
-44.4697
0.005324
1
-0.0007
0.00717
7
71
-44.4687
1.631e-005
1
0.000984
0.000523
989-44.46874.308e-013
1.25e-0094.11e-007
0.88961.24022.8702
x=3.24180.00001.63420.1240fval=-44.4687
hessian=
0.5738
0.3772
-0.2369
0.1036
-0.1250
-0.0757
-0.1739
0.3772
0.6898
0.4769
0.0641
-0.0816
-0.0879
-0.0803
-0.2369
0.4769
1.4590
0.0628
0.0657
0.1229
-0.0055
0.1036
0.0641
0.0628
0.8811
0.0721
0.0120
0.0699
-0.1250
-0.0816
0.0657
0.0721
1.7848
-0.2129
-0.0244
-0.0757
-0.0879
0.1229
0.0120
-0.2129
1.4472
-0.1761
-0.1739
-0.0803
-0.0055
0.0699
-0.0244
-0.1761
1.0268
c=-5.9342ceq=4.3077e-013从上述运行结果可以得出：
升级会员 