稳定流形和不稳定流形图
大家好,我有一个做稳定流形和不稳定流形图的程序,可是老是运行不出来,请高手指点。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 static/image/common/back.gif
没报错吗!?
若有的话, 建议给齐报错讯息
可是根据报错信息修改出来还是不正确呀,现在都不会修改了 yanxiangya 发表于 2012-12-16 16:55 static/image/common/back.gif
可是根据报错信息修改出来还是不正确呀,现在都不会修改了
请描述具体现象
页:
[1]