《梯度投影法 MATLAB程式可執行》由會員分享，可線上閱讀，更多相關《梯度投影法 MATLAB程式可執行(10頁珍藏版)》請在人人文庫網上搜尋。
1、function x,minf=minRosen(f,A,b,x0,var,eps)%目标函數：f;%限制矩陣：A;%限制右端力量：b;%初始可行點：x0;%自變量向量：var;%精度：eps;%目标函數取最小值時的自變量值：x;%目标函數的最小值：minf;format long;if nargin = 5eps=1.0e-6;endsyms l;x00=transpose(x0);n=length(var);sz=size(A);m=sz(1);gf=jacobian(f,var);bConti=1;while bContik=0;s=0;A1=A;A2=A;b1=b;b2=b;for i。
2、=1:mdfun=A(i,:)*x00-b(i);if abs(dfun)0A1=A1(1:k,:);b1=b1(1:k,:);endif s0A2=A2(1:s,:);b2=b2(1:s,:);endwhile 1P=eye(n,n);if k0tM=transpose(A1);P=P-tM*inv(A1*tM)*A1;endgv=Funval(gf,var,x0);gv=transpose(gv); d=-P*gv ;if d=0if k=0x=x0;bConti=0;break;elsew=inv(A1*tM)*A1*gv;if w=0 x=x0;bConti=0;break;elseu。
3、,index=min(w);sA1=size(A1);if sA1(1)=1k=0;elsek=sA1(2); A1=A1(1:(index-1),:);A1(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;if dd=0tmpI,lm=minJT(tmpf,0,0.1); elselm=inf;for i=1:length(dd)if dd(i)eps & k fua = l;l = u;u =。
4、 a + 0.618*(b - a);elseb = u;u = l;l = a + 0.382*(b-a);endk = k+1;tol = abs(b - a);endif k = 100000disp(找不到最小值！);x = NaN;minf = NaN;return;endx = (a+b)/2;minf = subs(f, findsym(f),x);format short;function minx,maxx = minJT(f,x0,h0,eps)format long;if nargin = 3eps = 1.0e-6;endx1 = x0;k = 0;h = h0;whi。
5、le 1x4 = x1 + h;k = k+1;f4 = subs(f, findsym(f),x4);f1 = subs(f, findsym(f),x1);if f4 0 A1=A1(1:k,:); b1=b1(1:k,:); end if s0 A2=A2(1:s,:); b2=b2(1:s,:); end while 1 P=eye(n,n); if k0 tM=transpose(A1); P=P-tM*inv(A1*tM)*A1; %calculate P; end gv=Funval(gf,var,x0); gv=transpose(gv); d=-P*gv; %calculat。
6、e the searching direction 計算搜尋方向% flg=1; % if(P=zeros(n) % flg =0; % end % if flg=1 % d=d/norm(d); %normorlize the searching direction 搜尋方向% end % 加入這部分會無止境的循環if d=0 if k=0 x=x0; bConti=0; break; else w=-inv(A1*tM)*A1*gv; if w=0 x=x0; bConti=0; break; else u,index=min(w);%find the negative component。
7、 in w sA1=size(A1); if sA1(1)=1 k=0; else k=sA1(2); A1=A1(1:(index-1),:);A1(index+1):sA1(1),:); %delete corresponding row in A1 删除對應的行A1end end end else break; end end d1=transpose(d);y1=x0+l*d1; %new iteration variable 新的疊代變量tmpf=Funval(f,var,y1); bb=b2-A2*x00; dd=A2*d; if dd=0 tmpI,lm= ForwardBack。
8、Method(tmpf,0,0.001); %find the searching interval 找到搜尋區間 else lm=inf; %find lambda_max for i=1:length(dd) % if(dd(i)0) % % if dd(i)0 保證 0%find the minimizer by one dimension searching method 找到由一維搜尋方法得到目标tol=norm(xm*d); if toleps&kfu a=l; l=u; u=a+0.618*(b-a); else b=u; u=l; l=a+0.382*(b-a); end k=。
9、k+1; tol=abs(b-a); end if k=100000 disp(找不到最小值！); x=NaN; minf=NaN; return; end x=(a+b)/2; %return the minimizer 傳回最小值minf=subs(f, findsym(f),x); format short; ForwardBackMethod.m 進退法确定搜尋區間function left,right=ForwardBackMethod(f,x0,step) if nargin=2 step = 0.01 end if nargin=1 x0=0; step = 0.01 end f。
10、0 =subs(f,findsym(f),x0); x1=x0+step; f1=subs(f,findsym(f),x1); if(f1f2) x0=x1; x1=x2; f0=f1; f1=f2; step=2*step; x2=x1+step; f2=subs(f,findsym(f),x2); end left=min(x0, x2); right=max(x0, x2); else step=2*step; x2=x1-step; f2=subs(f,findsym(f),x2); while(f0f2) x1=x0; x0=x2; f1=f0; f0=f2; step=2*step。
11、; x2=x1-step; f2=subs(f,findsym(f),x2); end left=min(x1,x2);%left end point right=max(x1,x2);%right end point end Funval.mfunction fv=Funval(f,varvec,varval) var=findsym(f); %找出表達式包含的變量t,s f=t2+s+1varc=findsym(varvec); %找出傳遞參數的變量t s中的 t ss1=length(var); %函數的個數s2=length(varc); %變量的個數m=floor(s1-1)/3+1。
12、); %靠近左邊的整數varv=zeros(1,m); if s1 = s2%if the number of variable is different, deal with it specially 如果變量的數量是不同的,專門處理它for i=0: (s1-1)/3) k=findstr(varc,var(3*i+1); index=(k-1)/3; varv(i+1) = varval(index+1); end fv=subs(f,var,varv); else fv=subs(f,varvec,varval); %return the value of the function 如果原來函數變量個數與傳遞參數中變量個數一緻，調用subs函數，計算給定點處函數值end % disp(here) % f % varvec % varval %fv = subs(f,varvec,varval);。