Josephon系统(Poincare图象和Lyapunov指数)
Josephon.m
function dx=Josephon(t,x);
% Josephon方程
% dx=Josephon(t,[x;y;b])
% t-时间，x,y-为自变量，b-为如下方程所示的参数
% eg: dx=Josephon(10,[0;0;.33803])
%
% 方程如下:
% θ''+G*θ'+sinθ=I+A*sin(ωt)+αsin(βωt)
% 变化：
% dx=y
% dy=-G*y-sin(x)+I+A*sin(w*t)+a*sin(b*w*t)
%
% Example(Poincare图象):
% [T,Y]=ode45('Josephon',[0,1000],[0;0;.33803]);
% plot(mod(Y(:,1),6.283),Y(:,2),'.','markersize',2);
%
% Author's email: ustb03-07@yahoo.com.cn
%
G=.7;A=.4;w=.25;
a=0.0125;I=.905;
% b=.33803;
b=x(3);
dx(1,1)=x(2);
dx(2,1)=-G*x(2)-sin(x(1))+I+A*sin(w*t)+a*sin(b*w*t);
dx(3,1)=0;
jose_ly.m
function ly=jose_ly(b,k)
% the largest lyapunov exponent of josephson
% k 迭代步数，b 参数
% 方程如下:
% θ''+G*θ'+sinθ=I+A*sin(ωt)+αsin(βωt)
% 变化：
% dx=y
% dy=-G*y-sin(x)+I+A*sin(w*t)+a*sin(b*w*t)
%
% Example:
% ly=jose_ly(0,800)
%
% Author:LDYU
% Author's email: ustb03-07@yahoo.com.cn
%
d0=1e-8;
ly=0;
lsum=0;
x=[0;2;b];
x1=[d0;2;b];
for t=1:k
[T1,Y1]=ode45('Josephon',[t-1,t],x);
[T2,Y2]=ode45('Josephon',[t-1,t],x1);
x=Y1(end,:);
x1=Y2(end,:);
d1=norm(x-x1);
x1=x+(d0/d1)*(x1-x);
lsum=lsum+log(d1/d0);
end
ly=lsum/k;
>> help jose_ly
the largest lyapunov exponent of josephson
k 迭代步数，b 参数
方程如下:
θ''+G*θ'+sinθ=I+A*sin(ωt)+αsin(βωt)
变化：
dx=y
dy=-G*y-sin(x)+I+A*sin(w*t)+a*sin(b*w*t)
Example:
ly=jose_ly(0,800)
Author:LDYU
Author's email: ustb03-07@yahoo.com.cn
>> ly=jose_ly(0,800)
ly =
0.0522
josel.m
% the largest lyapunov exponents of josephson
%
% Author's email: ustb03-07@yahoo.com.cn
%
clear
Z=[];
for b=linspace(0,1.5,50)
Z=[Z b+i*jose_ly(b,300)];
end
plot(Z,'-')
title('Josephon最大Lyapunov指数图'),xlabel('b'),ylabel('lyapunov')
grid on
>> josel
