大连理工优化方法大作业MATLAB编程Word文档下载推荐.docx
- 文档编号:20999968
- 上传时间:2023-01-26
- 格式:DOCX
- 页数:19
- 大小:36.33KB
大连理工优化方法大作业MATLAB编程Word文档下载推荐.docx
《大连理工优化方法大作业MATLAB编程Word文档下载推荐.docx》由会员分享,可在线阅读,更多相关《大连理工优化方法大作业MATLAB编程Word文档下载推荐.docx(19页珍藏版)》请在冰豆网上搜索。
注:
通过比习题验证共範梯度法求辉门无二次西数极小点至多需要“次迭代.
functionf=fun(X)
%所求问题目标函数
f=X
(1)A2-2*X
(1)*X
(2)+2*X
(2)A2+X(3)A2+X(4)A2-
X
(2)*X(3)+2*X
(1)+3*X
(2)-X(3);
endfunctiong=gfun(X)
%所求问题目标函数梯度
g=[2*X
(1)-2*X
(2)+2,-2*X
(1)+4*X
(2)-X(3)+3,2*X(3)-X
(2)-1,2*X(4)];
endfunction[x,val,k]=frcg(fun,gfun,xO)
%功能:
用FR共轭梯度法求无约束问题最小值
%输入:
x0是初始点,fun和gfun分别是目标函数和梯度
%输出:
x、val分别是最优点和最优值,k是迭代次数
maxk=5000;
%最大迭代次数
rho=0.5;
sigma=0.4;
eps=10e-6;
n=length(x0);
while(k<
maxk)
g=feval(gfun,x0);
%计算梯度itern=k-(n+1)*floor(k/(n+1));
itern=itern+1;
%计算搜索方向
if(itern==1)
d=-g;
else
beta=(g*g'
)/(g0*g0'
);
d=-g+beta*d0;
gd=g'
*d;
if(gd>
=0.0)
if(norm(g)<
eps)
break;
m=0;
mk=0;
while(m<
20)
if(feval(fun,xO+rhoAm*d)<
feval(fun,xO)+sigma*rhoAm*g'
*d)
mk=m;
break;
m=m+1;
x0=x0+rhoAmk*d;
val=feval(fun,x0);
g0=g;
d0=d;
x=x0;
>
x0=[0,0,0,0];
[x,val,k]=frcg('
fun'
'
gfun'
x0)
x=
-4.0000-3.0000-1.00000
val=
-8.0000
k=
或者
function[x,f,k]=second(x)
dk=dfun(x);
g0=gfun(x);
s=-g0;
x=x+dk*s;
g1=gfun(x);
while(norm(g1)>
=0.02)
if(k==3)
elseif(k<
3)
u=((norm(g1))A2)/(norm(gO)A2);
s=-g1+u*s;
g0=g1;
f=fun(x);
functionf=fun(x)
f=x(1F2-2*x
(1)*x
(2)+2*x
(2)A2+x(3)A2+x(4)A2-x
(2)*x(3)+2*x
(1)+3*x
(2)-x(3);
gf=[2*x
(1)-2*x
(2)+2,-2*x
(1)+4*x
(2)-x(3)+3,2*x(3)-x
(2)-1,2*x(4)];
function[p,q]=con(x,d)
ss=-gfun(x);
p=fun(x)-fun(x+d*ss)+0.2*d*gfun(x)*(ss)'
q=gfun(x+d*ss)*(ss)'
-0.6*gfun(x)*(ss)'
functiondk=dfun(x)
[p,q]=con(x,d);
if(p<
(q>
End
x=[0,0,0,0];
[x,f,k]=second(x)
-1.00900
x=-4.0147-3.0132
f=-7.9999
k=1
取初始点3=(0」)二考虑下面无约東优化问题:
minf(x)二冷+2x2+exp(xf+天孑),
其中歩长Qk的选取可別用习题1或精确一维搜索•搜索方向为一HNW
♦取垃=b
•取皿=R2f防)]"
9耳丈啟为BFG5公式亠
通过此习题体会上述三种算法的收敛速度.
function[f,x,k]=third_1(x)k=0;
g=gfun(x);
while(norm(g)>
=0.001)s=-g;
dk=dfun(x,s);
f=fun(x);
functionf=fun(x)
f=x
(1)+2*x
(2)A2+exp(x
(1)A2+x
(2)A2);
functiongf=gfun(x)
gf=[1+2*x
(1)*exp(x
(1)A2+x
(2)A2),4*x
(2)+2*x
(2)*(x
(1)A2+x
(2)A2)];
function[j_1,j_2]=con(x,d,s)
j_1=fun(x)-fun(x+d*s)+0.1*d*gfun(x)*(s)'
j_2=gfun(x+d*s)*(s)'
-0.5*gfun(x)*(s)'
functiondk=dfun(x,s)%获取步长flag=0;
[p,q]=con(x,d,s);
[f,x,k]=third_1(x)
f=0.7729
x=-0.41960.0001
k=8
(1)程序:
function[f,x,k]=third_2(x)
H=inv(ggfun(x));
g=gfun(x);
=0.001)
s=(-H*g'
)'
dk=dfun(x,s);
f=x
(1)+2*x
(2)A2+exp(x(1F2+x
(2)A2);
functiongf=gfun(x)gf=[1+2*x
(1)*exp(x(1F2+x
(2)A2),4*x
(2)+2*x
(2)*(x(1F2+x
(2)A2)];
functionggf=ggfun(x)
ggf=[(4*x
(1)A2+2)*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+(4*x
(2)A2+2)*exp(x
(1)A2+x
(2)A2)];
function[j_1,j_2]=con(x,d,s)
j_2=gfun(x+d*s)*(s)'
functiondk=dfun(x,s)%步长获取
b=10000;
if2*d>
=(d+b)/2
d=(d+b)/2;
elsed=2*d;
[f,x,k]=third_2(x)
x=-0.41930.0001
(2)程序:
function[f,x,k]=third_3(x)
X=cell
(2);
g=cell
(2);
X{1}=x;
H=eye
(2);
g{1}=gfun(X{1});
s=(-H*g{1}'
dk=dfun(X{1},s);
X{2}=X{1}+dk*s;
g{2}=gfun(X{2});
while(norm(g{2})>
dx=X{2}-X{1};
dg=g{2}-g{1};
v=dx/(dx*dg'
)-(H*dg'
/(dg*H*dg'
h1=H*dg'
*dg*H/(dg*H*dg'
h2=dx'
*dx/(dx*dx'
h3=dg*H*dg'
*v'
*v;
H=H-h1+h2+h3;
X{1}=X{2};
norm(g{2});
x=X{2};
ggf=[(4*x
(1)A2+2)*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+(4*x
(2)A2+2)*exp(x
(1)A2+x
(2)A2);
function[p,q]=con(x,d,s)
p=fun(x)-fun(x+d*s)+0.1*d*gfun(x)*(s)'
q=gfun(x+d*s)*(s)'
functiondk=dfun(x,s)
[f,x,k]=third_3(x)
0.0000
x=-0.4195
k=6
习题4
*U用有效集法求解下面勺勺二次规划问题:
mm
Si.
(XI一I)2+(x2一2.5)2
X1-2X2+2>
0
-Xi—2>
(2+6>
-Xi+2X2+2>
xi,x2>
functioncallqpact
H=[20;
02];
c=[-2-5]'
Ae=[];
be=[];
Ai=[1-2;
-1-2;
-12;
10;
01];
bi=[-2-6-200]'
x0=[00]'
[x,lambda,exitflag,output]=qpact(H,c,Ae,be,Ai,bi,xO)
function[x,lamk,exitflag,output]=qpact(H,c,Ae,be,Ai,bi,x0)
epsilon=1.0e-9;
err=1.0e-6;
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;
endwhile(k<
=kmax)
Aee=[];
if(ne>
0),Aee=Ae;
end
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>
exitflag=0;
exitflag=1;
if(index(i)&
(ne+sum(index(1:
i)))==jk)
index(i)=0;
alpha=1.0;
tm=1.0;
if((index(i)==0)&
(Ai(i,:
)*dk<
0))
tm1=(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;
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>
rb=Ae*ginvH*c+be;
lambda=pinv(Ae*ginvH*Ae'
)*rb;
x=ginvH*(Ae'
*lambda-c);
x=-ginvH*c;
lambda=0;
结果
callqpact
1.4000
1.7000
lambda=
0.8000
exitflag=
output=
fval:
-6.4500
iter:
7
function[x,mu,lambda,output]=multphr(fun,hf,gf,dfun,dhf,dgf,xO)
用乘子法解一般约束问题:
minf(x),s.t.h(x)=0,g(x).=0
x0是初始点,fun,dfun分别是目标函数及其梯度;
%hf,dhf分别是等式约束(向量)函数及其Jacobi矩阵的转置;
%gf,dgf分别是不等式约束(向量)函数及其Jacobi矩阵的转置;
x是近似最优点,mu,lambda分别是相应于等式约束和不等式约束的乘子向量
%output是结构变量,输出近似极小值f,迭代次数,内迭代次数等
maxk=500;
c=2.0;
eta=2.0;
theta=0.8;
ink=0;
epsilon=0.00001;
x=xO;
he=feval(hf,x);
gi=feval(gf,x);
n=length(x);
l=length(he);
m=length(gi);
mu=zeros(l,1);
lambda=zeros(m,1);
btak=10;
btaold=10;
while(btak>
epsilon&
k<
%调用BFGS算法程序求解无约束子问题
[x,ival,ik]=bfgs('
mpsi'
'
dmpsi'
x0,fun,hf,gf,dfun,dhf,dgf,mu,lambda,c);
ink=ink+ik;
btak=0;
fori=1:
l
btak=btak+he(y2;
%更新乘子向量
m
temp=min(gi(i),lambda(i)/c);
btak=btak+tempA2;
btak=sqrt(btak);
ifbtak>
epsilon
ifk>
=2&
btak>
theta*btaold
c=eta*c;
mu(i)=mu(i)-c*he(i);
lambda(i)=max(0,lambda(i)-c*gi(i));
btaold=btak;
x0=x;
f=feval(fun,x);
output.fval=f;
%增广拉格朗日函数
functionpsi=mpsi(x,fun,hf,gf,dfun,dhf,dgf,mu,lambda,c)
l=length(he);
m=length(gi);
psi=f;
s1=0;
psi=psi-he(i)*mu(i);
s仁s1+he(y2;
psi=psi+0.5*c*s1;
s2=0;
s3=max(0,lambda(i)-c*gi(i));
s2=s2+s3A2-lambda(i)A2;
psi=psi+s2/(2*c);
%不等式约束函数文件g1.m
functiongi=g1(x)
gi=10*x
(1)-x
(1)A2+10*x
(2)-x
(2)A2-34;
%目标函数的梯度文件df1.m
functiong=df1(x)
g=[4,-2*x
(2)]'
%等式约束(向量)函数的Jacobi矩阵(转置)文件dh1.m
functiondhe=dh1(x)
dhe=[-2*x
(1),-2*x
(2)]'
%不等式约束(向量)函数的Jacobi矩阵(转置)文件dg1.m
functiondgi=dg1(x)
dgi=[10-2*x
(1),10-2*x
(2)]'
function[x,val,k]=bfgs(fun,gfun,x0,varargin)maxk=500;
rho=0.55;
epsilon=0.00001;
Bk=eye(n);
gk=feval(gfun,x0,varargin{:
});
if(norm(gk)<
epsilon)
dk=-Bk\gk;
newf=feval(fun,x0+rhoAm*dk,varargin{:
oldf=feval(fun,x0,varargin{:
if(newf<
oldf+sigma*rhoAm*gk'
*dk)mk=m;
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);
val=feval(fun,x0,varargin{:
x=[22]'
[x,mu,lambda,output]=multphr('
hf'
gf1'
df'
dh'
dg'
x0)
1.0013
4.8987
mu=
0.7701
0.9434
-31.9923
4
习题6
利用序列二次规划方法求解习题5中的约束优化问题:
min4xi一好一12
s.t.25-x?
—x孑=Q
10xt一召+10旳-xj-34>
0X1,X2>
f=[3,1,1];
A=[2,1,1;
1,-1,-1];
b=[2;
-1];
lb=[0,0,0];
x=linprog(f,A,b,zeros(3),[0,0,0]'
lb)
Optimizationterminated.
0.5000
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 大连理工 优化 方法 作业 MATLAB 编程
![提示](https://static.bdocx.com/images/bang_tan.gif)