论坛 › Matlab讨论区 › 稳定流形和不稳定流形图

yanxiangya 发表于 2012-12-16 14:30

稳定流形和不稳定流形图

大家好，我有一个做稳定流形和不稳定流形图的程序，可是老是运行不出来，请高手指点。
            
function M = Umanifold(p0,Arc,alphamax, alphamin,deltamax,deltamin,deltaalphamax,deltaalphamin,deltak,pp,epsb,lamb,del,v,n)
%Compute 1-D unstable manifold         
global M            
%pp is vector of parameters
p1=p0+del*v;
M=;p = n2c(pp);
preimage=p0-del/lamb;
pleft=p0;lindex=1;
pright=p1;rindex=2;
arcl=norm(p1-p0);
while( arcl < Arc)
    =Search_Line(M, pleft,pright,preimage,lindex,rindex, deltak,deltamin,alphamax,pp,epsb,deltamax,n);
   
    pbar=M(:,end)+((norm(M(:,end)-M(:,end-1))/norm(M(:,end)-pcan)))*(M(:,end)-pcan);
    alphak=norm(pbar-M(:,end-1))/norm(M(:,end)-M(end-1));
   
    deltak=norm(pcan-M(:,end));
   

if ((alphak<alphamax)&& (deltak*alphak<deltaalphamax))   
      
      M(:,end+1)=pcan;
      arcl=arcl+deltak;
      pleft=preimage;
      pright=pj;      
   
   elseif ((alphak<alphamax)&& (deltak*alphak<deltamax))
   
      deltak=2*deltak;
elseif deltak > deltamin
      deltak=deltak/2;
         
end %if   

end %while
%-----------------------------------------------------------------------
function = Search_Line(M, pleft,pright,preimage,lindex,rindex, deltak,deltamin,alphamax,pp,epsb,deltamax,n);
=Find_intersection_On_circle(M, pleft,pright,preimage,lindex,rindex,deltak,deltamin,alphamax,pp,epsb,deltamax,n);
=Bisection(M,pleft,pright,preimage,deltak,epsb,pp,alphamax,n) ;
pi=pleft;pj=pright;
%--------------------------------------------------------------------------
function =Find_intersection_On_circle(M, pleft,pright,preimage,lindex,rindex,deltak,deltamin,alphamax,pp,epsb,deltamax,n)
                  
=Side_of_line(M, pleft,pright,preimage,lindex,rindex,deltak,deltamin,alphamax,pp,epsb,deltamax,n);
j=0;
while (S==0)&&(j<100)   
      j=j+1;
      if M(:,rindex)~=M(:,end)      
         lindex=rindex;rindex=rindex+1;
         pleft=M(:,lindex); pright=M(:,rindex);
      elseif deltak>deltamin
          deltak=deltak/2;
      end
      =Side_of_line(M, pleft,pright,preimage,lindex,rindex,deltak,deltamin,alphamax,pp,epsb,deltamax,n);
end   
% -------------------------------------------------------------------------
function =Side_of_line(M, pleft,pright,preimage,lindex,rindex,deltak,deltamin,alphamax,pp,epsb,deltamax,n)
p = n2c(pp);
fp=norm(fun_eval(0,preimage,p{:},n)-M(:,end));
fr=norm(fun_eval(0,pright,p{:},n)-M(:,end));
if ((fp<deltak)&&(fr>deltak))|((fp>deltak)&&(fr<deltak))
%The line segment L that is mapped by f,intersects the circle with center pk and radius deltak
    S=1;
else
    S=0;
end
    %--------------------------------------------------------------------
function =Bisection(M,pleft,pright,preimage,deltak,epsb,pp,alphamax,n)
pleft=preimage;
ptry=(pleft+pright)/2;
p=n2c(pp);
ft=norm(fun_eval(0,ptry,p{:},n)-M(:,end));
fl=norm(fun_eval(0,preimage,p{:},n)-M(:,end));
fr=norm(fun_eval(0,pright,p{:},n)-M(:,end));
if (abs(fl)> abs(fr))
      temp=fl;
      fl=fr;
      fr=temp;
end   
i=0;
while (ft>(1+epsb)*deltak) && (ft<(1-epsb)*deltak)&&(i<10)
   i=i+1;
if (ft>= deltak)
   % f(ptry) is on the same side as f(pend)   
      pright=ptry;
      ptry=(pleft+pright)/2;
   
else   
   pleft=ptry;
   ptry=(pleft+pright)/2;
end

end %while
p=n2c(pp);ftry=fun_eval(0,ptry,p{:},n);
pcan=ftry;preimage=ptry;
%------------------------------------------------------------------------
function map = fun_eval(t,kmrgd,alpha,beta,R,S,n)

   for i=1:n      
    map=;
    kmrgd(1)=map(1);
    kmrgd(2)=map(2);
   end
%--------------------------------------------------------------------------

ChaChing 发表于 2012-12-16 15:35

没报错吗!?
若有的话, 建议给齐报错讯息

yanxiangya 发表于 2012-12-16 16:55

ChaChing 发表于 2012-12-16 15:35 static/image/common/back.gif
没报错吗!?
若有的话, 建议给齐报错讯息

可是根据报错信息修改出来还是不正确呀，现在都不会修改了

happy 发表于 2012-12-17 10:33

yanxiangya 发表于 2012-12-16 16:55 static/image/common/back.gif
可是根据报错信息修改出来还是不正确呀，现在都不会修改了

请描述具体现象

查看完整版本: 稳定流形和不稳定流形图