资源描述:
《运筹学实验报告共轭梯度内外点法》由会员上传分享,免费在线阅读,更多相关内容在应用文档-天天文库。
1、运筹与优化实验报告姓名:学号:学院:专业:指导老师:实验一实验名称:MATLAB实现黄金分割法实验要求:用MATLAB实现黄金分割法()试验程序:function[x,minf]=minHJ(f,a,b,eps)formatlong;eps=0.16;symsx;f=2*x^2-x-1;a=-1;b=1;l=a+0.382*(b-a);u=a+0.618*(b-a);k=1;tol=b-a;whiletol>epsf1=2*l^2-l-1;fu=2*u^2-u-1;iff1>fua=1;l=u;u=a+0
2、.618*(b-a);elseb=u;u=l;l=a+0.382*(b-a);endk=k+1;tol=abs(b-a);endkifk==10000disp('找不到最小值!');x=NuN;minf=NuN;returnendx=(a+b)/2;min=subs(f,findsym(f),x);formatshort;运行结果:黄金分割法的结果:K=7,ans=-1.0341实验二实验名称:MATLAB实现共轭梯度法实验要求:用MATLAB实现共轭梯度法求解实例()试验程序:functionf=con
3、jugate(x0,t)x=x0;symsxi,yi,a;f=xi^2-xi*yi+yi^2+2*xi-4*yi;fx=diff(f,xi);fy=diff(f,yi);fx=subs(fx,{xi,yi},x0);fy=subs(fy,{xi,yi},x0);fi=[fx,fy];count=0;whiledouble(sqrt(fx^2+fy^2))>ts=-fi;ifcount<=0s=-fielsex=x+a*s;f=subs(f,{xi,yi},x);f1=diff(f);f1=solve(f)
4、;iff1~=0ai=double(f1);elsebreakx.f=subs(f,{xi,yi},x),countendx=subs(x,a,ai);f=xi^2-xi*yi+yi^2+2*xi-4*yi;fxi=diff(f,xi);fyi=diff(f,yi);fxi=subs(fxi,{xi,yi},x);fyi=subs(fyi,{xi,yi},x);fii=[fxi,fyi];运行结果:共轭梯度法:conjugate([0,0],0.000001),f=-4,count=2,ans=-4实验三
5、实验名称:MATLAB实现外点法实验要求:用MATLAB实现外点法和内点法求解实例:mins.t.,实验程序:m=zeros(1,50);a=zeros(1,50);b=zeros(1,50);f0=zeros(1,50);symsdx1x2e;m(1)=1;c=10;a(1)=0;b(1)=0;f=x1^2+x2^2+e*(1-x1)^2;f0(1)=1;fx1=diff(f,'x1');fx2=diff(f,'x2');fork=1:100x1=a(k);x2=b(k);e=m(k);forn=1:1
6、00f1=subs(fx1);f2=subs(fx2);if(double(sqrt(f1^2+f2^2))<=0.002)a(k+1)=double(x1);b(k+1)=double(x2);f0(k+1)=double(subs(f));break;elseD=(x1-d*f1)^2+(x2-d*f2)^2+e*(1-(x1-d*f1))^2;Dd=diff(D,'d');dd=solve(Dd);x1=x1-dd*f1;x2=x2-dd*f2;endendif(double(sqrt((a(k+1
7、)-a(k))^2+(b(k+1)-b(k))^2))<=0.001)&&(double(abs((f0(k+1)-f0(k))/f0(k)))<=0.001)a(k+1)b(k+1)kf0(k+1)break;elsem(k+1)=c*m(k);endend外点法运行结果:实验四实验名称:MATLAB实现内点法实验要求:用MATLAB实现内点法求解实例:mins.t.,实验程序:m=zeros(1,50);a=zeros(1,50);b=zeros(1,50);f0=zeros(1,50);symsx1
8、x2e;m(1)=1;c=0.2;a(1)=2;b(1)=-3;f=x1^2+x2^2-e*(1/(2*x1+x2-2)+1/(1-x1));f0(1)=15;fx1=diff(f,'x1');fx2=diff(f,'x2');fx1x1=diff(fx1,'x1');fx1x2=diff(fx1,'x2');fx2x1=diff(fx2,'x1');fx2x2=diff(fx2,'x2');fork=1:100x1=a(k)