梯度投影法-MATLAB程序可执行.docx

《梯度投影法-MATLAB程序可执行.docx》由会员上传分享，免费在线阅读，更多相关内容在工程资料-天天文库。

1、function[x,minf]=minRosen(f,A,b,x0,var,eps)%目标函数：f;%约束矩阵：A;%约束右端力量：b;%初始可行点：x0;%自变量向量：var;%精度：eps;%目标函数取最小值时的自变量值：x;%目标函数的最小值：minf;formatlong;ifnargin==5eps=1.0e-6;endsymsl;x00=transpose(x0);n=length(var);sz=size(A);m=sz(1);gf=jacobian(f,var);bConti=1;

2、whilebContik=0;s=0;A1=A;A2=A;b1=b;b2=b;fori=1:mdfun=A(i,:)*x00-b(i);ifabs(dfun)<0.k=k+1;A1(k,:)=A(i,:);b1(k,1)=b(i);elses=s+1;A2(s,:)=A(i,:);b2(s,1)=b(i);endendifk>0A1=A1(1:k,:);b1=b1(1:k,:);endifs>0A2=A2(1:s,:);b2=b2(1:s,:);endwhile1P=eye(n,n);ifk>0tM

3、=transpose(A1);P=P-tM*inv(A1*tM)*A1;endgv=Funval(gf,var,x0);gv=transpose(gv);d=-P*gv;ifd==0ifk==0x=x0;bConti=0;break;elsew=inv(A1*tM)*A1*gv;ifw>=0x=x0;bConti=0;break;else[u,index]=min(w);sA1=size(A1);ifsA1(1)==1k=0;elsek=sA1(2);A1=[A1(1:(index-1),:);A1

4、((index+1):sA1(2),:)];%去掉A1对应的行endendendelsebreak;endendd1=transpose(d);y1=x0+l*d1;tmpf=Funval(f,var,y1);bb=b2-A2*x00;dd=A2*d;ifdd>=0[tmpI,lm]=minJT(tmpf,0,0.1);elselm=inf;fori=1:length(dd)ifdd(i)<0ifbb(i)/dd(i)

5、inHJ(tmpf,0,lm,1.0e-14);tol=norm(xm*d);iftol

6、~=s2fori=0:((s1-1)/3)k=findstr(varc,var(3*i+1));index=(k-1)/3;varv(i+1)=varval(index+1);%index(i+1);%varv(i+1)=varval(index(i+1));endfv=subs(f,var,varv);elsefv=subs(f,varvec,varval);endfunction[x,minf]=minHJ(f,a,b,eps)formatlong;ifnargin==3eps=1.0e-6;e

7、ndl=a+0.382*(b-a);u=a+0.618*(b-a);k=1;tol=b-a;whiletol>eps&&kfua=l;l=u;u=a+0.618*(b-a);elseb=u;u=l;l=a+0.382*(b-a);endk=k+1;tol=abs(b-a);endifk==disp('找不到最小值！');x=NaN;minf=NaN;return;endx=(a+b)/2;mi

8、nf=subs(f,findsym(f),x);formatshort;function[minx,maxx]=minJT(f,x0,h0,eps)formatlong;ifnargin==3eps=1.0e-6;endx1=x0;k=0;h=h0;while1x4=x1+h;k=k+1;f4=subs(f,findsym(f),x4);f1=subs(f,findsym(f),x1);iff4