优化方法上机大作业Word文档下载推荐.docx
《优化方法上机大作业Word文档下载推荐.docx》由会员分享，可在线阅读，更多相关《优化方法上机大作业Word文档下载推荐.docx（24页珍藏版）》请在冰点文库上搜索。
x1=
-0.0000
1.0000d=
1.1102e-016fl=
2
习题2
取初始点X°
=（0,0,0,0）7,利）+]共純梯度法求解下而无约束优化问题：
minf（x）=好一2xiX2—2%孑+x孑+爲一x2x3一2xi+3x2-x3?
其中步长廊的逸取可利用习題1或精确一维捜索+
注:
通过此习题验证共梔梯度法求解门元二次禹数极小点至多需要门次迭代，
functionf=fun（x）
%UNTITLED3Summaryofthisfunctiongoeshere
%Detailedexplanationgoeshere
f=x（1F2-2*x
（1）*x
（2）+2*x
（2）A2+x（3F2+X⑷A2-x
（2）*x（3）+2*x
（1）+3*x
（2）-x（3）;
End
functiong=fun（x）
%UNTITLED4Summaryofthisfunctiongoeshere
g=[2-200;
-24-10;
0-120;
0002]*x+[2;
3;
-1;
0];
O;
O];
%初始值
eps=1.0e-4;
%精度
gO=gfun（xO）;
sO=-gO;
n=O;
symsd1;
whilenorm（gO）>
eps
if*3
g=gfun（x0+d1*s0）;
d=double（solve（s0'
*g））;
x1=x0+d*s0;
g1=gfun（x1）;
ifnorm（g1）<
n=n+1;
x0=x1;
break
sO=-g1+（norm（g1F2/norm（g0F2）*s0;
x0=x1;
gO=g1;
elseifn=3
n=0;
xO
n
fun（xO）
xO=
-4
-3
-1
ans=
-8
习题3
取初始点0=（0,1）1考虑下的无约束优化问题：
minf（x）=Xi+2迓+expfx^+於），
-Jt中步长®
的选取可利用习题1鼓精确一维捜索.搜索方向为-佔心），
°
取Hk=I2
。
取皿=[V2"
xk）]T
»
取总为BFGS公扎
通过此习题休会上述三种算法的收效速度.
functionf=fun3_1（x）
%FUN3Summaryofthisfunctiongoeshere
f=x
（1）+2*x
（2）A2+exp（x
（1）A2+x
（2）A2）;
functiong=gfun3_1（x）
%GFUN3_1Summaryofthisfunctiongoeshere
g=[1+2*x
（1）*exp（x
（1）A2+x（2F2）;
4*x
（2）+2*x
（2）*exp（x（1F2+x
⑵A2）];
（1）最速下降法
x0=[0;
1];
eps=1.0e-5;
%精度
g0=gfun3_1（x0）;
whilenorm（g0）>
=eps
s0=-g0;
g=gfun3_1（x0+d1*s0）;
d=double（solve（s0'
g1=gfun3_1（x1）;
if（norm（gl）veps）n=n+1;
break;
elsex0=x1;
f0=fun3_1（x0）x0nf0=
0.7729x0=
-0.4194
0.0000
（2）牛顿法
functiong2=hesse（x）
%HESSESummaryofthisfunctiongoeshere
g2=[2*exp（x
（1）A2+x
（2）A2）+4*x
（1）A2*exp（x
（1）A2+x
（2）A2）,4*x
（1）*x
（2）*exp（
x
（1）A2+x
（2）A2）
4*x
（1）*x
（2）*exp（x
（1）A2+x
（2）A2）,4+2*exp（x
（1）A2+x
（2）A2）+4*x（2f2*exp（x（
1）A2+x
（2）A2）];
gO=gfun3_1（x0）;
g20=hesse（x0）;
whilenorm（g0）>
d=-g20\g0;
x仁x0+d;
g1=gfun3_1（x1）;
if（norm（g1）<
eps）n=n+1;
g0=gfun3_1（x0）;
end
endf0=fun3_1（x0）x0
f0=
0.0000n=
（3）BFGS
gO=gfun3_1（x0）;
hO=eye
（2）;
s0=-h0*g0;
if（norm（g1）<
eps）
break;
h0=h0-（h0*（g1-g0）*（g1-g0）'
*h0）/（（g1-g0）'
*h0*（g1-g0））+...（（x1-x0）*（x1-x0）'
）/（（x1-x0）'
*（g1-g0））+...
（（g1-g0）'
*h0*（g1-g0））*（（x1-x0）*（x1-x0）'
）;
f0=fun3_1（x0）
x0
nf0=
function[x,lamk,exitflag,output]=qpact（H,c,Ae,be,Ai,bi,x0
）
epsilon=1.0e-9;
err=1.0e-6;
k=0;
x=x0;
n=length（x）;
kmax=1.0e3;
ne=length（be）;
ni=length（bi）;
lamk=zeros（ne+ni,1）;
index=ones（ni,1）;
for（i=1:
ni）
if（Ai（i,:
）*x>
bi（i）+epsilon）,index（i）=0;
while（k<
=kmax）
Aee=[];
if（ne>
0）,Aee=Ae;
for（j=1:
if（index（j）>
0）,Aee=[Aee;
Ai（j,:
）];
gk=H*x+c;
[m1,n1]=size（Aee）;
[dk,lamk]=qsubp（H,gk,Aee,zeros（m1,1））;
if（norm（dk）<
=err）
y=0.0;
if（length（lamk）>
ne）
[y,jk]=min（lamk（ne+1:
length（lamk）））;
if（y>
=0）
exitflag=0;
exitflag=1;
for（i=1:
if（index（i）&
（ne+sum（index（1:
i）））==jk）index（i）=0;
k=k+1;
alpha=1.0;
tm=1.0;
if（（index（i）==0）&
（Ai（i,:
）*dk<
0））tm仁（bi（i）-Ai（i,:
）*x）/（Ai（i,:
）*dk）;
if（tm1<
tm）
tm=tm1;
ti=i;
alpha=min（alpha,tm）;
x=x+alpha*dk;
if（tm<
1）,index（ti）=1;
endend
if（exitflag==0）,break;
output.fval=0.5*x'
*H*x+c'
*x;
output.iter=k;
function[x,lambda]=qsubp（H,c,Ae,be）
ginvH=pinv（H）;
[m,n]=size（Ae）;
if（m>
0）
rb=Ae*ginvH*c+be;
lambda=pinv（Ae*ginvH*Ae'
）*rb;
x=ginvH*（Ae'
*lambda-c）;
x=-ginvH*c;
lambda=zeros（m,1）;
callqpact.m文件
functioncallqpactH=[20;
02];
c=[-2-5]'
;
Ae=[];
be=[];
Ai=[1-2;
-1-2;
-12;
10;
01];
bi=[-2-6-200]'
xO=[O0]'
[x,lambda,exitflag,output]=qpact（H,c,Ae,be,Ai,bi,xO）
运行结果：
callqpact
x=
1.4000
1.7000
lambda=
0.8000
exitflag=
output=
fval:
-6.4500
iter:
7
习题5
利用乘予法求解下面的约束优化问题：
min4xi—x孑一12
s,t.25——X2=0
10x1—xf+10x2—X2—34>
0
Xr,X2>
其中每次迭代需求解一个无约束极小化问题」可利用习题2或习
题3中的程序.
Function
[x,mu,lambda,output]=multphr（fun,hf,gf,dfun,dhf,dgf,x0）
function
psi=mpsi（x,fun,hf,gf,dfun,dhf,dgf,mu,lambda,sigma）
functionhe=h1（x）
he=-x
（1）A2-x
（2）A2+25.0;
functionf=f1（x）
f=4*x
（1）-x
（2）A2-12;
functiongi=g1（x）
gi=10*x
（1）-x
（1）A2+10*x
（2）-x
（2）八2-34;
functiondhe=dh1（x）
dhe=[-1*x
（1）,-1*x
（2）]'
functiondgi=dg1（x）
dgi=[10-2*x
（1）,10-2*x
（2）]'
functiong=df1（x）
g=[4,-2.0*x
（2）]'
function[x,val,k]=bfgs（fun,gfun,x0,varargin）maxk=500;
sigma=2.0;
eta=2.0;
theta=0.8;
ink=0;
epsilon=1e-5;
x=x0;
he=feval（hf,x）;
gi=feval（gf,x）;
n=length（x）;
l=length（he）;
m=length（gi）;
mu=0.1*ones（l,1）;
lambda=0.1*ones（m,1）;
btak=10;
btaold=10;
while（btak>
epsilon&
k<
maxk）
[x,ival,ik]=bfgs（'
mpsi'
'
dmpsi'
x0,fun,hf,gf,dfun,dhf,dgf,mu,lambda,sigma）;
ink=ink+ik;
he=feval（hf,x）;
btak=0.0;
l）,btak=btak+he（i）A2;
fori=1:
m
temp=min（gi（i）,lambda（i）/sigma）;
btak=btak+tempA2;
btak=sqrt（btak）;
ifbtak>
epsilon
if（k>
=2&
btak>
theta*btaold）sigma=eta*sigma;
l）,mu（i）=mu（i）-sigma*he（i）;
m）lambda（i）=max（0.0,lambda（i）-sigma*gi（i））;
btaold=btak;
x0=x;
f=feval（fun,x）;
output.fval=f;
output.iter=k;
output.inner_iter=ink;
output.bta=btak;
f=feval（fun,x）;
l=length（he）;
psi=f;
s1=0.0;
l）
psi=psi-he（i）*mu（i）;
s仁s1+he（i）A2;
psi=psi+0.5*sigma*s1;
s2=0.0;
m）
s3=max（0.0,lambda（i）-sigma*gi（i））;
s2=s2+s3A2-lambda（i）A2;
psi=psi+s2/（2.0*sigma）;
maxk=500;
rho=0.55;
sigma1=0.4;
epsilon1=1e-5;
n=length（x0）;
Bk=eye（n）;
while（kvmaxk）
gk=feval（gfun,x0,varargin{:
}）;
if（norm（gk）<
epsilon1）,break;
dk=-Bk\gk;
m=0;
mk=0;
while（m<
20）
newf=feval（fun,x0+rh°
Am*dk,varargin{:
oldf=feval（fun,x0,varargin{:
if（newf<
oldf+sigma1*rhoAm*gk'
*dk）
mk=m;
m=m+1;
x=x0+rhoAmk*dk;
sk=x-x0;
yk=feval（gfun,x,varargin{:
}）-gk;
if（yk'
*sk>
Bk=Bk-（Bk*sk*sk'
*Bk）/（sk'
*Bk*sk）+（yk*yk'
）/（yk'
*sk）;
endk=k+1;
x0=x;
val=feval（fun,x0,varargin{:
dpsi=dmpsi（x,fun,hf,gf,dfun,dhf,dgf,mu,lambda,sigma）dpsi=feval（dfun,x）;
dhe=feval（dhf,x）;
dgi=feval（dgf,x）;
dpsi二dpsi+（sigma*he（i）-mu（i））*dhe（:
i）;
m）dpsi=dpsi+（sigma*gi（i）-lambda（i））*dgi（:
>
xO=[1,1]'
[x,mu,lambda,output]=multphr（'
f1'
h1'
g1'
df1'
dh1'
dg1'
x0）x=
1.0013
4.8987mu=
2.0312lambda=0.7545output=
-31.9923iter:
5
inner_iter:
58
bta:
4.3187e-07
习题6
利fl序列二次规划方法求解习题5中的约束优化问题:
min
s.t.
4xi—xf—12
25-好一烤=0
10xi—xf+10x2—x孑一34二0
Xi,X2>
function[x,mu,lam,val,k]=sqpm（x0,mu0,lam0）maxk=100;
n=length（x0）;
l=length（mu0）;
m=length（lam0）;
rho=0.5;
eta=0.1;
B0=eye（n）;
mu=mu0;
lam=lam0;
Bk=BO;
sigma=0.8;
epsilon1=1e-6;
epsilon2=1e-5;
[hk,gk]=cons（x）;
dfk=df1（x）;
[Ae,Ai]=dcons（x）;
Ak=[Ae;
Ai];
k=0;
[dk,mu,lam]=qpsubp（dfk,Bk,Ae,hk,Ai,gk）;
mp1=norm（hk,1）+norm（max（-gk,0）,1）;
if（norm（dk,1）<
epsilon1）&
（mp1<
epsilon2）break;
deta=0.05;
tau=max（norm（mu,inf）,norm（lam,inf））;
if（sigma*（tau+deta）<
1）sigma=sigma;
sigma=1.0/（tau+2*deta）;
im=0;
while（im<
=20）
if（phi1（x+rhoAim*dk,sigma）-phi1（x,sigma）<
eta*rhoAim*dphi1（x,sigma,dk））mk=im;
im=im+1;
if（im==20）,mk=10;
alpha=rhoAmk;
x1=x+alpha*dk;
[hk,gk]=cons（x1）;
dfk=df1（x1）;
[Ae,Ai]=dcons（x1）;
lamu=pinv（Ak）'
*dfk;
if（l>
0&
m>
mu=lamu（1:
l）;
lam=lamu（l+1:
l+m）;
if（l==0）,mu=[];
lam=lamu;
if（m==0）,mu=lamu;
lam=[];
sk=alpha*dk;
yk=dlax（x1,mu,lam）-dlax（x,mu,lam）;
if（sk'
*yk>
0.2*sk'
*Bk*sk）
theta=1;
theta=0.8*sk'
*Bk*sk/（sk'
*Bk*sk-sk'
*yk）;
zk=theta*yk+（1-theta）*Bk*sk;
Bk=Bk+zk*zk'
/（sk'
*zk）-（Bk*sk）*（Bk*sk）'
*Bk*sk）;
x=x1;
val=f1（x）;
functionp=phi1（x,sigma）
f=f1（x）;
[h,g]=cons（x）;
gn=max（-g,0）;
IO=length（h）;
mO=length（g）;
if（IO==O）,p=f+1.0/sigma*norm（gn,1）;
endif（m0==0）,p=f+1.0/sigma*norm（h,1）;
endif（l0>
m0>
p=f+1.0/sigma*（norm（h,1）+norm（gn,1））;
functiondp=dphi1（x,sigma,d）df=df1（x）;
[h,g]=cons（x）;
l0=length（h）;
m0=length（g）;
if（l0==0）,dp=df*
升级会员 