大连理工大学优化方法上机作业.docx
- 文档编号:2224810
- 上传时间:2022-10-28
- 格式:DOCX
- 页数:12
- 大小:360.06KB
大连理工大学优化方法上机作业.docx
《大连理工大学优化方法上机作业.docx》由会员分享,可在线阅读,更多相关《大连理工大学优化方法上机作业.docx(12页珍藏版)》请在冰豆网上搜索。
大连理工大学优化方法上机作业
大连理工大学优化方法上机作业
优化方法上机大作业
学院:
电子信息与电气工程学部
姓名:
学号:
指导老师:
上机大作业
(一)
%目标函数
functionf=fun(x)
f=100*(x
(2)-x
(1)^2)^2+(1-x
(1))^2;
end
%目标函数梯度
functiongf=gfun(x)
gf=[-400*x
(1)*(x
(2)-x
(1)^2)-2*(1-x
(1));200*(x
(2)-x
(1)^2)];
End
%目标函数Hess矩阵
functionHe=Hess(x)
He=[1200*x
(1)^2-400*x
(2)+2,-400*x
(1);
-400*x
(1),200;];
end
%线搜索步长
functionmk=armijo(xk,dk)
beta=0.5;sigma=0.2;
m=0;maxm=20;
while(m<=maxm)
if(fun(xk+beta^m*dk)<=fun(xk)+sigma*beta^m*gfun(xk)'*dk)
mk=m;break;
end
m=m+1;
end
alpha=beta^mk
newxk=xk+alpha*dk
fk=fun(xk)
newfk=fun(newxk)
%最速下降法
function[k,x,val]=grad(fun,gfun,x0,epsilon)
%功能:
梯度法求解无约束优化问题:
minf(x)
%输入:
fun,gfun分别是目标函数及其梯度,x0是初始点,
%epsilon为容许误差
%输出:
k是迭代次数,x,val分别是近似最优点和最优值
maxk=5000;%最大迭代次数
beta=0.5;sigma=0.4;
k=0;
while(k gk=feval(gfun,x0);%计算梯度 dk=-gk;%计算搜索方向 if(norm(gk) m=0;mk=0; while(m<20)%用Armijo搜索步长 if(feval(fun,x0+beta^m*dk)<=feval(fun,x0)+sigma*beta^m*gk'*dk) mk=m;break; end m=m+1; end x0=x0+beta^mk*dk; k=k+1; end x=x0; val=feval(fun,x0); >>x0=[0;0]; >>[k,x,val]=grad('fun','gfun',x0,1e-4) 迭代次数: k= 1033 x= 0.9999 0.9998 val= 1.2390e-008 %牛顿法 x0=[0;0];ep=1e-4;maxk=10;k=0; while(k gk=gfun(x0); if(norm(gk) x=x0 miny=fun(x) k0=k break; else H=inv(Hess(x0)); x0=x0-H*gk; k=k+1; end end x= 1.0000 1.0000 miny= 4.9304e-030 迭代次数 k0= 2 %BFGS方法 function[k,x,val]=bfgs(fun,gfun,x0,varargin) %功能: 梯度法求解无约束优化问题: minf(x) %输入: fun,gfun分别是目标函数及其梯度,x0是初始点, %epsilon为容许误差 %输出: k是迭代次数,x,val分别是近似最优点和最优值 N=1000; epsilon=1e-4; beta=0.55;sigma=0.4; n=length(x0);Bk=eye(n); k=0; while(k gk=feval(gfun,x0,varargin{: }); if(norm(gk) dk=-Bk\gk; m=0;mk=0; while(m<20) newf=feval(fun,x0+beta^m*dk,varargin{: }); oldf=feval(fun,x0,varargin{: }); if(newf<=oldf+sigma*beta^m*gk'*dk) mk=m;break; end m=m+1; end x=x0+beta^mk*dk; sk=x-x0; yk=feval(gfun,x,varargin{: })-gk; if(yk'*sk>0) Bk=Bk-(Bk*sk*sk'*Bk)/(sk'*Bk*sk)+(yk*yk')/(yk'*sk); end k=k+1; x0=x; end val=feval(fun,x0,varargin{: }); >>x0=[0;0]; >>[k,x,val]=bfgs('fun','gfun',x0) k= 20 x= 1.0000 1.0000 val= 2.2005e-011 %共轭梯度法 function[k,x,val]=frcg(fun,gfun,x0,epsilon,N) ifnargin<5,N=1000;end ifnargin<4,epsilon=1e-4;end beta=0.6;sigma=0.4; n=length(x0);k=0; while(k gk=feval(gfun,x0); itern=k-(n+1)*floor(k/(n+1)); itern=itern+1; if(itern==1) dk=-gk; else betak=(gk'*gk)/(g0'*g0); dk=-gk+betak*d0;gd=gk'*dk; if(gd>=0),dk=-gk;end end if(norm(gk) m=0;mk=0; while(m<20) if(feval(fun,x0+beta^m*dk)<=feval(fun,x0)+sigma*beta^m*gk'*dk) mk=m;break; end m=m+1; end x=x0+beta^m*dk; g0=gk;d0=dk; x0=x;k=k+1; end val=feval(fun,x); >>x0=[0;0]; [k,x,val]=frcg('fun','gfun',x0,1e-4,1000) k= 122 x= 1.0001 1.0002 val= 7.2372e-009 上机大作业 (二) %目标函数 functionf_x=fun(x) f_x=4*x (1)-x (2)^2-12; %等式约束条件 functionhe=hf(x) he=25-x (1)^2-x (2)^2; end %不等式约束条件 functiongi_x=gi(x,i) switchi case1 gi_x=10*x (1)-x (1)^2+10*x (2)-x (2)^2-34; case2 gi_x=x (1); case3 gi_x=x (2); otherwise end %求目标函数的梯度 functionL_grad=grad(x,lambda,cigma) d_f=[4;2*x (2)]; d_g(: 1)=[-2*x (1);-2*x (2)]; d_g(: 2)=[10-2*x (1);10-2*x (2)]; d_g(: 3)=[1;0]; d_g(: 4)=[0;1]; L_grad=d_f+(lambda (1)+cigma*hf(x))*d_g(: 1); fori=1: 3 iflambda(i+1)+cigma*gi(x,i)<0 L_grad=L_grad+(lambda(i+1)+cigma*gi(x,i))*d_g(: i+1); continue end end %增广拉格朗日函数 functionLA=lag(x,lambda,cee) LA=fun(x)+lambda (1)*hf(x)+0.5*cee*hf(x)^2; fori=1: 3 LA=LA+1/(2*cee)*(min(0,lambda(i+1)+cee*gi(x,i))^2-lambda(i+1)^2); end functionxk=BFGS(x0,eps,lambda,cigma) gk=grad(x0,lambda,cigma);res_B=norm(gk);k_B=0; a_=1e-4;rho=0.5;c=1e-4; length_x=length(x0);I=eye(length_x);Hk=I; whileres_B>eps&&k_B<=10000 dk=-Hk*gk; m=0; whilem<=5000 iflag(x0+a_*rho^m*dk,lambda,cigma)-lag(x0,lambda,cigma)<=c*a_*rho^m*gk'*dk mk=m; break; end m=m+1; end ak=a_*rho^mk; xk=x0+ak*dk; delta=xk-x0; y=grad(xk,lambda,cigma)-gk; Hk=(I-(delta*y')/(delta'*y))*Hk*(I-(y*delta')/(delta'*y))+(delta*delta')/(delta'*y); k_B=k_B+1;x0=xk;gk=y+gk;res_B=norm(gk); end %增广拉格朗日法 functionval_min=ALM(x0,eps) lambda=zeros(4,1);cigma=5;alpha=10;k=1; res=[abs(hf(x0)),0,0,0]; fori=1: 3 res(1,i+1)=norm(min(gi(x0,i),-lambda(i+1)/cigma)); end res=max(res); whileres>eps&&k<1000 xk=BFGS(x0,eps,lambda,cigma); lambda (1)=lambda (1)+cigma*hf(xk); fori=1: 3 lambda(i+1)=lambda(i+1)+min(0,lambda(i+1)+gi(x0,1)); end k=k+1;cigma=alpha*cigma;x0=xk; res=[norm(hf(x0)),0,0,0]; fori=1: 3 res(1,i+1)=norm(min(gi(x0,i),-lambda(i+1)/cigma)); end res=max(res); end val_min=fun(xk); fprintf('k=%d\n',k); fprintf('fmin=%.4f\n',val_min); fprintf('x=[%.4f;%.4f]\n',xk (1),xk (2)); >>x0=[0;0]; >>val_min=ALM(x0,1e-4) k=10 fmin=-31.4003
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 大连理工大学 优化 方法 上机 作业