优化方法MATLAB编程 大连理工大学.docx
- 文档编号:3045032
- 上传时间:2023-05-05
- 格式:DOCX
- 页数:19
- 大小:278.36KB
优化方法MATLAB编程 大连理工大学.docx
《优化方法MATLAB编程 大连理工大学.docx》由会员分享,可在线阅读,更多相关《优化方法MATLAB编程 大连理工大学.docx(19页珍藏版)》请在冰点文库上搜索。
优化方法MATLAB编程大连理工大学
优化方法上机大作业
学院:
姓名:
学号:
指导老师:
肖现涛
第一题
源程序如下:
functionzy_x=di1ti(x)
%di1ti是用来求解优化作业第一题的函数。
x0=x;yimuxulong=0.000001;
g0=g(x0);s0=-g0;
A=2*ones(100,100);
k=0;
whilek<100
lanmed=-(g0)'*s0/(s0'*A*s0);
x=x0+lanmed*s0;
g=g(x);
k=k+1;
ifnorm(g) zy_x=x; fprintf('After%diterations,obtaintheoptimalsolution.\n\nTheoptimalsolutionis\n%f.\n\nTheoptimal"x"is"ans".',k,f(x)) break; end miu=norm(g)^2/norm(g0)^2; s=-g+miu*s0; g0=g;s0=s;x0=x; end functionf=f(x) f=(x'*ones(100,1))^2-x'*ones(100,1); functiong=g(x) g=(2*x'*ones(100,1))*ones(100,1)-ones(100,1); 代入x0,运行结果如下: >>x=zeros(100,1); >>di1ti(x) After1iterations,obtaintheoptimalsolution. Theoptimalsolutionis -0.250000. Theoptimal"x"is"ans". ans=0.005*ones(100,1). 第二题 1.最速下降法。 源程序如下: functionzy_x=di2titidu(x) %该函数用来解大作业第二题。 x0=x;yimuxulong=1e-5;k=0; g0=g(x0);s0=-g0; whilek>=0 ifnorm(g0) break; else lanmed=10;c=0.1;i=0; whilei>=0&i<100 x=x0+lanmed*s0; iff(x)>(f(x0)+c*lanmed*g0'*s0) lanmed=lanmed/2; i=i+1; else break; end end x=x0+lanmed*s0; x0=x; g0=g(x); s0=-g0; k=k+1; end end zy_x=x; zyj=f(x); fprintf('after%diterations,obtaintheoptimalsolution.\n\nTheoptimalsolutionis%f.\n\nTheoptimal"x"is"ans".\n',k,zyj); functionf=f(x) x1=[1000]*x; x2=[0100]*x; x3=[0010]*x; x4=[0001]*x; f=(x1-1)^2+(x3-1)^2+100*(x2-x1^2)^2+100*(x4-x3^2)^2; functiong=g(x) x1=[1000]*x; x2=[0100]*x; x3=[0010]*x; x4=[0001]*x; g=[2*(x1-1)-400*x1*(x2-x1^2);200*(x2-x1^2);2*(x3-1)-400*x3*(x4-x3^2);200*(x4-x3^2)]; >>x=[-1.21-1.21]'; >>di2titidu(x) after5945iterations,obtaintheoptimalsolution. Theoptimalsolutionis0.000000. Theoptimal"x"is"ans". ans= 1.0000 1.0000 1.0000 1.0000 2.牛顿法 源程序如下: functionzy_x=di2tinewton(x) %该函数用来解大作业第二题。 x0=x;yimuxulong=1e-5;k=0; g0=g(x0);h0=h(x0);s0=-inv(h0)*g0; whilek>=0&k<1000 ifnorm(g0) break; else x=x0+s0; x0=x; g0=g(x); h0=h(x); s0=-inv(h0)*g0; k=k+1; end end zy_x=x; zyj=f(x); fprintf('after%diterations,obtaintheoptimalsolution.\n\nTheoptimalsolutionis%f.\n\nTheoptimal"x"is"ans".\n',k,zyj); functionf=f(x) x1=[1000]*x; x2=[0100]*x; x3=[0010]*x; x4=[0001]*x; f=(x1-1)^2+(x3-1)^2+100*(x2-x1^2)^2+100*(x4-x3^2)^2; functiong=g(x) x1=[1000]*x; x2=[0100]*x; x3=[0010]*x; x4=[0001]*x; g=[2*(x1-1)-400*x1*(x2-x1^2);200*(x2-x1^2);2*(x3-1)-400*x3*(x4-x3^2);200*(x4-x3^2)]; functionh=h(x) x1=[1000]*x; x2=[0100]*x; x3=[0010]*x; x4=[0001]*x; h=[2+1200*x1^2-400*x2-400*x100;-400*x120000;002+1200*x3^2-400*x4-400*x3;00-400*x3200]; 代入初始值,运行结果如下: >>x=[-1.21-1.21]'; >>di2tinewton(x) after6iterations,obtaintheoptimalsolution. Theoptimalsolutionis0.000000. Theoptimal"x"is"ans". ans= 1.0000 1.0000 1.0000 1.0000 可以看出,用Newton法经过6次迭代就能求出最优解。 3.BFGS法 源程序如下: functionzy_x=di2tiBFGS(x) %该函数用来解大作业第二题。 x0=x;yimuxulong=1e-5;k=0; g0=g(x0);H0=eye(4);s0=-H0*g0; whilek>=0&k<100 ifnorm(g0) break; else lanmed=10;c=0.1;i=0; whilei>=0&i<100 x=x0+lanmed*s0; iff(x)>(f(x0)+c*lanmed*g0'*s0) lanmed=lanmed/2; i=i+1; else break; end end x=x0+lanmed*s0; dete_x=x-x0; dete_g=g(x)-g0; miu=1+dete_g'*H0*dete_g/(dete_x'*dete_g); H=H0+(miu*dete_x*dete_x'-H0*dete_g*dete_x'-dete_x*dete_g'*H0)/(dete_x'*dete_g); s=-H*g(x); x0=x; s0=s; H0=H; g0=g(x); k=k+1; end end zy_x=x; zyj=f(x); fprintf('after%diterations,obtaintheoptimalsolution.\n\nTheoptimalsolutionis%f.\n\nTheoptimal"x"is"ans".\n',k,zyj); functionf=f(x) x1=[1000]*x; x2=[0100]*x; x3=[0010]*x; x4=[0001]*x; f=(x1-1)^2+(x3-1)^2+100*(x2-x1^2)^2+100*(x4-x3^2)^2; functiong=g(x) x1=[1000]*x; x2=[0100]*x; x3=[0010]*x; x4=[0001]*x; g=[2*(x1-1)-400*x1*(x2-x1^2);200*(x2-x1^2);2*(x3-1)-400*x3*(x4-x3^2);200*(x4-x3^2)]; 代入初始值,计算结果如下: >>x=[-1.21-1.21]'; >>di2tiBFGS(x) after53iterations,obtaintheoptimalsolution. Theoptimalsolutionis0.000000. Theoptimal"x"is"ans". ans= 1.0000 1.0000 1.0000 1.0000 第三题 1.惩罚函数法 源程序如下: functionzy_x=di3ti(x) %该函数用来解大作业第三题。 x0=x;M=100;c=4;m=1; whilem>0 g0=g(x0,M);yimuxulong=1e-5;k=0;s0=-inv(H(x0,M))*g0; whilek>=0 ifnorm(g0) break; else x=x0+s0;%牛顿法; x0=x; g0=g(x,M); s0=-inv(H(x0,M))*g0; k=k+1; end end ifmax([abs(h(x)),g1(x),g2(x),g3(x)])<0.5 break; else M=M*c; m=m+1; end end zy_x=x; zyj=f(x); fprintf('after%diterations,obtaintheoptimalsolution.\n\nTheoptimalsolutionis%f.\n\nTheoptimal"x"is"ans".\n',m,zyj); functionF=F(x,M) x1=[10]*x; x2=[01]*x; F=4*x1-x2^2-12+M*(h^2+g1^2+g2^2+g3^2); functiong=g(x,M) x1=[10]*x; x2=[01]*x; g=[4+M*(-4*(25-x1^2-x2^2)*x1+2*(10*x1-x1^2+10*x2-x2^2-34)*(10-2*x1)+2*x1);-2*x2+M*(-4*(25-x1^2-x2^2)*x2+2*(10*x1-x1^2+10*x2-x2^2-34)*(10-2*x2)+2*x2)]; functionH=H(x,M) x1=[10]*x; x2=[01]*x; H=[M*(24*x1^2-120*x1+8*x2^2-40*x2+238),M*(16*x1*x2-40*x1-40*x2+200);M*(16*x1*x2-40*x1-40*x2+200),-2+M*(24*x2^2-120*x2+8*x1^2-40*x1+238)]; functionf=f(x) x1=[10]*x; x2=[01]*x; f=4*x1-x2^2-12; functionh=h(x) x1=[10]*x; x2=[01]*x; h=25-x1^2-x2^2; functiong1=g1(x) x1=[10]*x; x2=[01]*x; g=10*x1-x1^2+10*x2-x2^2-34; ifg<0 g1=g; else g1=0; end functiong2=g2(x) x1=[10]*x; x2=[01]*x; ifx1>=0 g2=0; else g2=x1; end functiong3=g3(x) x1=[10]*x; x2=[01]*x; ifx2>=0 g3=0; else g3=x2; end 代入任意初始值,运算结果如下。 >>x=rand(2,1); >>di3ticf(x) after1iterations,obtaintheoptimalsolution. Theoptimalsolutionis-31.490552. Theoptimal"x"is"ans". ans= 1.0024 4.8477 2.乘子法 源程序如下: function[x,mu,lambda,output]=multphr(fun,hf,gf,dfun,dhf,dgf,x0) %功能: 用乘子法解一般约束问题: 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; k=0;ink=0; epsilon=0.00001; x=x0;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; he=feval(hf,x);gi=feval(gf,x); btak=0; fori=1: l btak=btak+he(i)^2; end %更新乘子向量 fori=1: m temp=min(gi(i),lambda(i)/c); btak=btak+temp^2; end btak=sqrt(btak); ifbtak>epsilon ifk>=2&&btak>theta*btaold c=eta*c; end fori=1: l mu(i)=mu(i)-c*he(i); end fori=1: m lambda(i)=max(0,lambda(i)-c*gi(i)); end k=k+1; btaold=btak; x0=x; end end f=feval(fun,x); output.fval=f; output.iter=k; %增广拉格朗日函数 functionpsi=mpsi(x,fun,hf,gf,dfun,dhf,dgf,mu,lambda,c) f=feval(fun,x);he=feval(hf,x);gi=feval(gf,x); l=length(he);m=length(gi); psi=f;s1=0; fori=1: l psi=psi-he(i)*mu(i); s1=s1+he(i)^2; end psi=psi+0.5*c*s1; s2=0; fori=1: m s3=max(0,lambda(i)-c*gi(i)); s2=s2+s3^2-lambda(i)^2; end psi=psi+s2/(2*c); %不等式约束函数文件g1.m functiongi=g1(x) gi=10*x (1)-x (1)^2+10*x (2)-x (2)^2-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;sigma=0.4;epsilon=0.00001; k=0;n=length(x0); Bk=eye(n); while(k gk=feval(gfun,x0,varargin{: }); if(norm(gk) break; end dk=-Bk\gk; m=0;mk=0; while(m<20) newf=feval(fun,x0+rho^m*dk,varargin{: }); oldf=feval(fun,x0,varargin{: }); if(newf mk=m; break; end m=m+1; end x=x0+rho^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{: }); 结果 x=[22]'; [x,mu,lambda,output]=multphr('fun','hf','gf1','df','dh','dg',x0) x= 1.0013 4.8987 mu= 0.7701 lambda= 0 0 0.9434 output= fval: -31.9923 iter: 4
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 优化方法MATLAB编程 大连理工大学 优化 方法 MATLAB 编程