機(jī)械優(yōu)化設(shè)計(jì)MATLAB程序

1、機(jī)械優(yōu)化設(shè)計(jì)作業(yè)e=0.01 。M 文件，如下：1. 用二次插值法求函數(shù)極小值，精度 在 MATLAB 的 M 文件編輯器中編寫的f=in li ne('(t+1)*(t-2)A2','t') a=0;b=3;epsilon=0.01;t1=a;f1=f(t1);t3=b;f3=f(t3); t2=0.5*(t1+t3);f2=f(t2); c1=(f3-f1)/(t3-t1); c2=(f2-f1)/(t2-t1)-c1)/(t2-t3); t4=0.5*(t1+t3-c1/c2);f4=f(t4);k=0;while(abs(t4-t2)>=epsil

2、on)if t2<t4if f2>f4f1=f2;t1=t2;t2=t4;f2=f4;elsef3=f4;t3=t4;end elseif f2>f4f3=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=0.5*(t1+t3-c1/c2);f4=f(t4); k=k+1;end%輸出最優(yōu)解if f2>f4t=t4;f=f(t4);elset=t2;f=f(t2);endfprintf(1,' 迭代計(jì)算 k=%3.0f

3、n',k) fprintf(1,' 極小點(diǎn)坐標(biāo) t=%3.0fn',t) fprintf(1,' 函數(shù)值 f=%3.4fn',f) 運(yùn)行結(jié)果如下：迭代計(jì)算 k= 7極小點(diǎn)坐標(biāo) t= 2函數(shù)值 f=0.00012. 用黃金分割法求函數(shù)的極小值，精度e=0.01 。在 MATLAB 的 M 文件編輯器中編寫的 M 文件，如下：f=i nlin e('tA(2/3)-(tA2+1)A(1/3)','t');a=0;b=3;epsilon=0.01;t1=b-0.618*(b-a);f1=f(t1);t2=a+0.618*(b-a

4、);f2=f(t2);k=1;while abs(b-a)>=epsilonif f1<f2b=t2;t2=t1;f2=f1;t1=b-0.618*(b-a);f1=f(t1);elsea=t1;t1=t2;f1=f2;t2=a+0.618*(b-a);f2=f(t2);endt=0.5*(b+a);k=k+1;f0=f(t);endfprintf(1,' 迭代次數(shù) k=% 3.0fn',k)fprintf(1,' 迭代區(qū)間 -左端 a=%3.4fn',a)fprintf(1,' 試點(diǎn) 1 坐標(biāo)值 t1=%3.4fn',t1)fpri

5、ntf(1,' 函數(shù)值 f1=%3.4fn',f(t1)fprintf(1,' 迭代區(qū)間 -右端 b=%3.4fn',b)fprintf(1,' 試點(diǎn) 2坐標(biāo)值 t2=%3.4fn',t2)fprintf(1,' 函數(shù)值 f2=%3.4fn',f(t2)fprintf(1,' 區(qū)間中點(diǎn) t=%3.4fn',t)fprintf(1,' 函數(shù)值 f0=%3.4fn',f(t)運(yùn)行結(jié)果如下：迭代次數(shù) k= 13迭代區(qū)間 -左端 a=0.0000試點(diǎn) 1 坐標(biāo)值 t1=0.0036函數(shù)值 f1=-0.976

6、7迭代區(qū)間 -右端 b=0.0093試點(diǎn) 2 坐標(biāo)值 t2=0.0058函數(shù)值 f2=-0.9679區(qū)間中點(diǎn) t=0.0047函數(shù)值 f0=-0.9721由黃金分割法在初始區(qū)間 0，3求得的極小值點(diǎn)為 t=0.0047 ，極小值為 -0.9721 。3. 用牛頓法、阻尼牛頓法及變尺度法求函數(shù)的極小點(diǎn)。(1) 在用牛頓法在 MATLAB的M文件編輯器中編寫的 M文件，如下： function x,fx,k=niudunfa(x0)syms x1 x2f=(x1-2)M+(x1-2*x2)A2;fx=0;v=x1,x2;df=jacobian(f,v);df=df.'G=jacobian(

7、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;x0=x0+p;while(norm(g1)>epson)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;endx=x0;fx=subs(f,x1,x2,x(1,1),x(2,1);運(yùn)行結(jié)果如下：>> x,fx,k=niudunfa(1;1)x =1.137789

8、70.fx =0.000000000000000007M 文件，如下：k =23(2) 用阻尼牛頓法在 MATLAB 的 M 文件編輯器中編寫的 function x,fx,k=zuniniudunfa(x0)% 阻尼牛頓法 syms x1 x2f=(x1-2)A4+(x1-2*x2)A2;fx=0;v=x1,x2;df=jacobian(f,v);df=df.'G=jacobian(df,v);epso n=1e-12;% 停機(jī)原則g1=subs(df,x1,x2,x0(1,1),x0(2,1);G1=subs(G,x1,x2,x0(1,1),x0(2,1);k=0;%迭代次數(shù)p=-

9、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;endx=x0;fx=subs(f,x1,x2,x0(1,1),x0(2,1);運(yùn)行結(jié)果如下：>>x,fx,k=zuniniudunfa(1;1)x=1.13778970.fx=0.00000000

10、0000000007k=23(3) 用變尺度法在 MATLAB的M文件編輯器中編寫的M文件，如下：4. 用共軛梯度法求函數(shù)的極小點(diǎn)(1 )用共軛梯度法在 MATLAB的M文件編輯器中編寫的M文件，如下:functiony,x,k=CG(A,b,c,x0)%共軛梯度法解 minf ( x) =0.5*X'*A*X+b'x+ceps=1e-6;%迭代停機(jī)原則%fx=0.5*x0'.*A.*x0+b'.*x0+c;r0=A*x0+b;if norm(r0)<=epsx=x0;y=0.5*x'*A*x+b'*x+c;k=0;endp0=-r0;a=

11、-r0'*p0/(p0'*A*p0);x1=x0+a*p0;r1=A*x1+b;k=0;while norm(r1)>eps beta=(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;endx=x1;y=0.5*x'*A*x+b'*x+c; 運(yùn)行結(jié)果如下： y,x,k=CG(3 -1;-1 1,-2;0,0,2;1)y = -1x = 1.00001.00

12、00k = 1(2) 用變尺度法在 MATLAB 的 M 文件編輯器中編寫的 M 文件，如下： function x,fx,k=bianchidufa(A,b,c,x0)%用變尺度法求 fx=0.5*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 k

13、case10x0=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 g1=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);

14、 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;breakend end x=x1;fx=0.5*x1'*A*x1+b'*x1+c; 運(yùn)行結(jié)果如下： x,fx,k=bianchidufa(3 -1;-1 1,-2;0,0,2;1) H1 =0.4031 0.25780.2578 0.8945 fx = -1 x =1.00001.0000 fx = -1 k = 1故函數(shù)極小點(diǎn)是點(diǎn)( 1,1)5. 用鮑威爾法求函數(shù)的極小點(diǎn)。M 文件，如下：用鮑威爾

15、法在MATLAB的M文件編輯器中編寫的 function x,fx,k=bowell(A,b,c,x0)% 鮑威爾法 d01=1;0;d02=0;1; x02=0;0;esp=1e-12;%停機(jī)原則 k=0;%迭代次數(shù) 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=-d

16、10'*g10/(d10'*A*d10); x10=x0+a01*d01;d01=d02; d02=d10;x0=x10;endx=x0;fx=0.5*x'*A*x+b'*x+c; 運(yùn)行結(jié)果如下：x,fx,k=bowell(2 -2;-2 4,-4;0,0,2;1)fx =-8x =42fx =-8k =36. 用單純形法求線性規(guī)劃問題用單純形法在 MATLAB的M文件編輯器中編寫的M文件，如下：%單純形法 matlab 程序 -danchunxingfa% 求解標(biāo)準(zhǔn)型線性規(guī)劃 :max c*x; s.t. A*x=b; x>=0% 本函數(shù)中的 A 是單純

17、初始表，包括 :最后一行是初始的檢驗(yàn)數(shù)，最后一列是資源向量b% N 是初始的基變量的下標(biāo)% 輸出變量 sol 是最優(yōu)解 , 其中松弛變量(或剩余變量)可能不為 0% 輸出變量 val 是最優(yōu)目標(biāo)值， kk 是迭代次數(shù) function sol,val,kk=danchunxingfa(A,N) mA,nA=size(A);kk=0; % 迭代次數(shù)flag=1;while flagkk=kk+1;if A(mA,:)<=0 % 已找到最優(yōu)解flag=0; sol=zeros(1,nA-1); for i=1:mA-1sol(N(i)=A(i,nA);end val=-A(mA,nA);el

18、sefor i=1:nA-1if A(mA,i)>0&A(1:mA-1,i)<=0 % 問題有無(wú)界解 disp('have infinite solution!');flag=0; break;endendif flag % 還不是最優(yōu)表，進(jìn)行轉(zhuǎn)軸運(yùn)算temp=0;for i=1:nA-1if A(mA,i)>temp temp=A(mA,i); inb=i; % 進(jìn)基變量的下標(biāo) end end sita=zeros(1,mA-1); for i=1:mA-1if A(i,inb)>0 sita(i)=A(i,nA)/A(i,inb);enden

19、dtemp=inf;for i=1:mA-1if sita(i)>0&sita(i)<temp temp=sita(i); outb=i; % 出基變量下標(biāo) end end % 以下更新 N for i=1:mA-1if i=outbN(i)=inb;endend% 以下進(jìn)行轉(zhuǎn)軸運(yùn)算A(outb,:)=A(outb,:)/A(outb,inb);for i=1:mAif i=outbA(i,:)=A(i,:)-A(outb,:)*A(i,inb);endendendendend ;運(yùn)行結(jié)果如下：>> A=1 11 0 4;1 21.1 2.22.5 3-3.3

20、4.45;0;N=3;4;sol,val,kk=danchunxingfa(A,N)sol=004.0000 1.6667val=7.3333kk =2所以，7. 求解線性規(guī)劃問題 用單純形法在 MATLAB 的 M 文件編輯器中編寫的 M 文件，如下：%單純形法 matlab 程序 -danchunxingfa% 求解標(biāo)準(zhǔn)型線性規(guī)劃 :max c*x; s.t. A*x=b; x>=0% 本函數(shù)中的 A 是單純初始表，包括 : 最后一行是初始的檢驗(yàn)數(shù)，最后一列是資源向量b% N 是初始的基變量的下標(biāo)% 輸出變量 sol 是最優(yōu)解 , 其中松弛變量(或剩余變量)可能不為 0% 輸出變量 val 是最優(yōu)目標(biāo)值， kk 是迭代次數(shù) function sol,val,kk=danchunxingfa(A,N)mA,nA=size(A);kk=0; % 迭代次數(shù)flag=1;while flagkk=kk+1;if A(mA,:)<=0 % 已找到最優(yōu)解flag=0; sol=zeros(1,nA-1);for i=1:mA-1sol(N(i)=A(i,nA);end val=-A(mA,nA);elsefor i=1:nA-1i