急求:用matlab来画出Lorenz的关于参数r的分叉图程序
急求:用matlab来画出Lorenz的关于参数r的分叉图程序老板要求做出来谢谢各位大虾 clear all
global r
t0=;%积分时间
y0=;
%bifurcation
for r=20:0.05:30 %r的变化精度
=ode45('Lorenz',t0,y0);
=getmax(y);
plot(r,Xmax,'b','markersize',1)
hold on
clear Xmax
end
xlabel('r')
ylabel('Xmax')
function = getmax(y)
a=length(y);
j=1;
for i=(a-1)/2:a
b=(y(i,1)-y(i-2,1))/2;
c=(y(i,1)+y(i-2,1))/2-y(i-1,1);
if y(i-2,1)<=y(i-1,1)&y(i-1,1)>=y(i,1)&c==0
Xmax(j)=y(i-1,1);
j=j+1;
elseif y(i-2,1)<=y(i-1,1)&y(i-1,1)>=y(i,1)
Xmax(j)=y(i-1,1)-b^2/(4*c);
j=j+1;
end
end
function dy = Lorenz(t,y)
global r
dy=zeros(3,1);
dy(1)=-10*(y(1)-y(2));
dy(2)=-y(1)*y(3)+r*y(1)-y(2);
dy(3)=y(1)*y(2)-8*y(3)/3;
仅供参考!要想图漂亮一点,可以把积分时间和r的变化精度增大,但是计算时间会变长。 jumpwolf谢谢你
小问题
两点建议1 该程序没有考虑过渡过程的剔除,建议每次积分都忽略前面一段足够长时间的解。因为如果是周期解或不动点,这样会得到错误的分岔图。
2 该程序是用某个分量的局部极大值描绘分岔图,其实没有必要使用两两比较的办法,这样如果数据量大会运行很长的时间。以x分量为例,x的局部极大值就是x的导数为零的点,只需考虑导数值(就是第一个方程的右边)的符号即可,另外应该考虑插值以得到足够的精度。 楼上的说的很好
页:
[1]