机械优化设计MATLAB程序.docx

1、机械优化设计MATLAB程序机械优化设计MATLAB程序机械优化设计作业1.用二次插值法求函数tt1t22极小值，精度e=。 在MATLAB的M文件xx器中编写的M文件，如下： f=inline(t+1)*(t-2)2,t) a=0;b=3;epsilon=; t1=a;f1=f(t1); t3=b;f3=f(t3);t2=*(t1+t3);f2=f(t2); c1=(f3-f1)/(t3-t1);c2=(f2-f1)/(t2-t1)-c1)/(t2-t3); t4=*(t1+t3-c1/c2);f4=f(t4); k=0;while(abs(t4-t2)=epsilon)if t2f4f1=

2、f2;t1=t2;t2=t4;f2=f4;elsef3=f4;t3=t4;endelseif f2f4f3=f2;t3=t2;t2=t4;f2=f4;elsef1=f4;t2=t4;endendc1=(f3-f1)/(t3-t1);c2=(f2-f1)/(t2-t1)-c1)/(t2-t3);t4=*(t1+t3-c1/c2);f4=f(t4);k=k+1; end%输出最优解 if f2f4t=t4;f=f(t4); elset=t2;f=f(t2); endfprintf(1,迭代计算k=%n,k) fprintf(1,极小点坐标t=%n,t) fprintf(1,函数值f=%n,f)运行

3、结果如下： 迭代计算k= 7 极小点坐标t= 2 函数值f=2.用黄金分割法求函数t3t23t21的极小值，精度e=。在MATLAB的M文件xx器中编写的M文件，如下： f=inline(t(2/3)-(t2+1)(1/3),t); a=0;b=3;epsilon=; t1=*(b-a);f1=f(t1); t2=a+*(b-a);f2=f(t2); k=1;while abs(b-a)=epsilonif f1epson)p=-G1g1;x0=x0+p;g1=subs(df,x1,x2,x0(1,1),x0(2,1);G1=subs(G,x1,x2,x0(1,1),x0(2,1);k=k+1

4、; end x=x0;fx=subs(f,x1,x2,x(1,1),x(2,1); 运行结果如下： x,fx,k=niudunfa(1;1)x =fx = k =23(2)用阻尼牛顿法在MATLAB的M文件xx器中编写的M文件，如下：function x,fx,k=zuniniudunfa(x0)%阻尼牛顿法 syms x1 x2f=(x1-2)4+(x1-2*x2)2; fx=0;v=x1,x2;df=jacobian(f,v); df=df.;G=jacobian(df,v);epson=1e-12;%停机原则g1=subs(df,x1,x2,x0(1,1),x0(2,1);G1=subs

5、(G,x1,x2,x0(1,1),x0(2,1); k=0;%迭代次数 p=-G1g1;a0=-p*g1/(p*G1*p); x0=x0+a0*p;while(norm(a0*p)epson)p=-G1g1;a0=-p*g1/(p*G1*p);x0=x0+a0*p;g1=subs(df,x1,x2,x0(1,1),x0(2,1);G1=subs(G,x1,x2,x0(1,1),x0(2,1);k=k+1; end x=x0;fx=subs(f,x1,x2,x0(1,1),x0(2,1);运行结果如下：x,fx,k=zuniniudunfa(1;1)x=fx= k=23(3)用变尺度法在MATL

6、AB的M文件xx器中编写的M文件，如下：32124.用共轭梯度法求函数fx1,x2x1x2x1x22x1的极小点22用共轭梯度法在MATLAB的M文件xx器中编写的M文件，如下：functiony,x,k=CG(A,b,c,x0)%共轭梯度法解minf=*X*A*X+bx+c eps=1e-6;%迭代停机原则%fx=*x0.*A.*x0+b.*x0+c; r0=A*x0+b; if norm(r0)epsbeta=(r1*r1)/(r0*r0);p1=-r1+beta*p0;alpha=-(r1*p1)/(p1*A*p1);x1=x1+alpha*p1;r2=A*x1+b;p0=p1;r0=r

7、1;r1=r2;k=k+1; end x=x1;y=*x*A*x+b*x+c;运行结果如下：y,x,k=CG(3 -1;-1 1,-2;0,0,2;1) y = -1x =k = 1(2)用变尺度法在MATLAB的M文件xx器中编写的M文件，如下：function x,fx,k=bianchidufa(A,b,c,x0) %用变尺度法求fx=*x*A*x+b*x+c; epson=1e-12; g0=A*x0+b; G0=A; H0=eye(2); k=0; d0=-H0*g0;a0=-d0*g0/(d0*G0*d0); s0=a0*d0; %x(k+1)-x(k); y0=A*a0*d0;

8、%g(k+1)-g(k); x1=x0+a0*d0;while (norm(s0)=epson) switch kcase10 x0=x1; g0=A*x0+b; H0=eye(2); k=0; d0=-H0*g0;a0=-d0*g0/(d0*G0*d0); s0=a0*d0; x1=x0+a0*d0;breakotherwiseg1=A*x1+b;y0=A*a0*d0;s0=a0*d0;% H1=H0+s0*s0/(s0*y0)-H0*y0*y0*H0/(y0*H0*y0);H1=H0+(1+y0*H0*y0/(s0*y0)*s0*s0-H0*y0*s0-s0*y0*H0)/(s0*y0);

9、k=k+1;d1=-H1*g1;a1=-d1*g1/(d1*G0*d1);a0=a1;d0=d1;H0=H1;s0=a0*d0;x1=x1+a0*d0;break end end x=x1;fx=*x1*A*x1+b*x1+c;运行结果如下： x,fx,k=bianchidufa(3 -1;-1 1,-2;0,0,2;1) H1 =fx =-1 x =fx =-1 k =1故函数极小点是点25.用鲍威尔法求函数fx1,x2x122x24x12x1x2的极小点。用鲍威尔法在MATLAB的M文件xx器中编写的M文件，如下：function x,fx,k=bowell(A,b,c,x0)%鲍威尔法

10、d01=1;0; d02=0;1; x02=0;0; esp=1e-12;%停机原则 k=0;%迭代次数while norm(x0-x02)=espk=k+1;g01=A*x0+b;a01=-d01*g01/(d01*A*d01);x01=x0+a01*d01;g02=A*x01+b;a02=-d02*g02/(d02*A*d02);x02=x01+a02*d02;d10=x02-x0;g10=A*x02+b;a10=-d10*g10/(d10*A*d10);x10=x0+a01*d01;d01=d02;d02=d10;x0=x10; end x=x0;fx=*x*A*x+b*x+c;运行结果

11、如下：x,fx,k=bowell(2 -2;-2 4,-4;0,0,2;1) fx =-8 x =42 fx =-8 k =36.用单纯形法求线性规划问题minf(x)x1x2x34x12x23x45x0(j1,2,3,4)j%单纯形法matlab程序-danchunxingfa% 求解标准型线性规划:max c*x; A*x=b; x=0% 本函数中的A是单纯初始表，包括:最后一行是初始的检验数，最后一列是资源向量b % N是初始的基变量的下标% 输出变量sol是最优解, 其中松弛变量可能不为0 % 输出变量val是最优目标值，kk是迭代次数function sol,val,kk=danch

12、unxingfa(A,N) mA,nA=size(A); kk=0; % 迭代次数 flag=1; while flag用单纯形法在MATLAB的M文件xx器中编写的M文件，如下：kk=kk+1;if A(mA,:)0&A(1:mA-1,i)temptemp=A(mA,i);inb=i; % 进基变量的下标 endendsita=zeros(1,mA-1);for i=1:mA-1if A(i,inb)0 sita(i)=A(i,nA)/A(i,inb);endendtemp=inf; for i=1:mA-1if sita(i)0&sita(i) A=1 1 1 0 4;1 23 5;-0;

13、N=3;4;sol,val,kk=danchunxingfa(A,N) sol =00val =kk =2 所以。经两次转轴运算，得到的最优解为x1x20,x3,x4,minf(x)7. 求解线性规划问题minz7x112x29x14x2x33604x5xx2003x110x2x5300xj0(j1,2,3,4,5)用单纯形法在MATLAB的M文件xx器中编写的M文件，如下：%单纯形法matlab程序-danchunxingfa% 求解标准型线性规划:max c*x; A*x=b; x=0% 本函数中的A是单纯初始表，包括:最后一行是初始的检验数，最后一列是资源向量b % N是初始的基变量的下

14、标% 输出变量sol是最优解, 其中松弛变量可能不为0 % 输出变量val是最优目标值，kk是迭代次数function sol,val,kk=danchunxingfa(A,N) mA,nA=size(A);kk=0; % 迭代次数 flag=1; while flagkk=kk+1;if A(mA,:)0&A(1:mA-1,i)temptemp=A(mA,i);inb=i; % 进基变量的下标 endendsita=zeros(1,mA-1);for i=1:mA-1if A(i,inb)0 sita(i)=A(i,nA)/A(i,inb);endendtemp=inf; for i=1:m

15、A-1if sita(i)0&sita(i) A=9 4 1 0 0 360;4 5 0 1 0 200;3 10 0 0 1 300; 7 12 0 0 0 0;N=3;4;5;sol,val,kk=danchunxingfa(A,N) sol =20248400 val =420 kk =3所以，经3次转轴运算，得到的最优解为x120,x224,x384,x4x50,minz420.endend% 以下进行转轴运算A(outb,:)=A(outb,:)/A(outb,inb);for i=1:mAif i=outb A(i,:)=A(i,:)-A(outb,:)*A(i,inb); end

16、endendend end ;运行结果如下：令Y7x112x2,则求minz7x112x2可转变为求maxY,即minz-maxY. A=9 4 1 0 0 360;4 5 0 1 0 200;3 10 0 0 1 300; 7 12 0 0 0 0;N=3;4;5;sol,val,kk=danchunxingfa(A,N) sol =20248400 val =420 kk =3所以，经3次转轴运算，得到的最优解为x120,x224,x384,x4x50,minz420.机械优化设计作业1.用二次插值法求函数tt1t22极小值，精度e=。 在MATLAB的M文件xx器中编写的M文件，如下：

17、f=inline(t+1)*(t-2)2,t) a=0;b=3;epsilon=; t1=a;f1=f(t1); t3=b;f3=f(t3);t2=*(t1+t3);f2=f(t2); c1=(f3-f1)/(t3-t1);c2=(f2-f1)/(t2-t1)-c1)/(t2-t3); t4=*(t1+t3-c1/c2);f4=f(t4); k=0;while(abs(t4-t2)=epsilon)if t2f4f1=f2;t1=t2;t2=t4;f2=f4;elsef3=f4;t3=t4;endelseif f2f4f3=f2;t3=t2;t2=t4;f2=f4;elsef1=f4;t2=t

18、4;endendc1=(f3-f1)/(t3-t1);c2=(f2-f1)/(t2-t1)-c1)/(t2-t3);t4=*(t1+t3-c1/c2);f4=f(t4);k=k+1; end%输出最优解 if f2f4t=t4;f=f(t4); elset=t2;f=f(t2); endfprintf(1,迭代计算k=%n,k) fprintf(1,极小点坐标t=%n,t) fprintf(1,函数值f=%n,f)运行结果如下： 迭代计算k= 7 极小点坐标t= 2 函数值f=2.用黄金分割法求函数t3t23t21的极小值，精度e=。在MATLAB的M文件xx器中编写的M文件，如下： f=in

19、line(t(2/3)-(t2+1)(1/3),t); a=0;b=3;epsilon=; t1=*(b-a);f1=f(t1); t2=a+*(b-a);f2=f(t2); k=1;while abs(b-a)=epsilonif f1epson)p=-G1g1;x0=x0+p;g1=subs(df,x1,x2,x0(1,1),x0(2,1);G1=subs(G,x1,x2,x0(1,1),x0(2,1);k=k+1; end x=x0;fx=subs(f,x1,x2,x(1,1),x(2,1); 运行结果如下： x,fx,k=niudunfa(1;1)x =fx = k =23(2)用阻尼

20、牛顿法在MATLAB的M文件xx器中编写的M文件，如下：function x,fx,k=zuniniudunfa(x0)%阻尼牛顿法 syms x1 x2f=(x1-2)4+(x1-2*x2)2; fx=0;v=x1,x2;df=jacobian(f,v); df=df.;G=jacobian(df,v);epson=1e-12;%停机原则g1=subs(df,x1,x2,x0(1,1),x0(2,1);G1=subs(G,x1,x2,x0(1,1),x0(2,1); k=0;%迭代次数 p=-G1g1;a0=-p*g1/(p*G1*p); x0=x0+a0*p;while(norm(a0*p

21、)epson)p=-G1g1;a0=-p*g1/(p*G1*p);x0=x0+a0*p;g1=subs(df,x1,x2,x0(1,1),x0(2,1);G1=subs(G,x1,x2,x0(1,1),x0(2,1);k=k+1; end x=x0;fx=subs(f,x1,x2,x0(1,1),x0(2,1);运行结果如下：x,fx,k=zuniniudunfa(1;1)x=fx= k=23(3)用变尺度法在MATLAB的M文件xx器中编写的M文件，如下：32124.用共轭梯度法求函数fx1,x2x1x2x1x22x1的极小点22用共轭梯度法在MATLAB的M文件xx器中编写的M文件，如下：

22、functiony,x,k=CG(A,b,c,x0)%共轭梯度法解minf=*X*A*X+bx+c eps=1e-6;%迭代停机原则%fx=*x0.*A.*x0+b.*x0+c; r0=A*x0+b; if norm(r0)epsbeta=(r1*r1)/(r0*r0);p1=-r1+beta*p0;alpha=-(r1*p1)/(p1*A*p1);x1=x1+alpha*p1;r2=A*x1+b;p0=p1;r0=r1;r1=r2;k=k+1; end x=x1;y=*x*A*x+b*x+c;运行结果如下：y,x,k=CG(3 -1;-1 1,-2;0,0,2;1) y = -1x =k =

23、 1(2)用变尺度法在MATLAB的M文件xx器中编写的M文件，如下：function x,fx,k=bianchidufa(A,b,c,x0) %用变尺度法求fx=*x*A*x+b*x+c; epson=1e-12; g0=A*x0+b; G0=A; H0=eye(2); k=0; d0=-H0*g0;a0=-d0*g0/(d0*G0*d0); s0=a0*d0; %x(k+1)-x(k); y0=A*a0*d0; %g(k+1)-g(k); x1=x0+a0*d0;while (norm(s0)=epson) switch kcase10 x0=x1; g0=A*x0+b; H0=eye(2); k=0; d0=-H0*g0;a0=-d0*g0/(d0*G0*d0); s0=a0*d0; x1=x0+a0*d0;breakotherwise