MATLAB常见非线性可视化绘制方法-分岔图与庞加莱截面(混沌可视化、Poincare截面、Logistic、Henon、Lorenz、Rossler、Duffing系统)
本文首发于 matlab爱好者 微信公众号,欢迎关注。
惯例声明:本人没有相关的工程应用经验,只是纯粹对相关算法感兴趣才写此博客。所以如果有错误,欢迎在评论区指正,不胜感激。本文主要关注于算法的实现,对于实际应用等问题本人没有任何经验,所以也不再涉及。
1 引言
本文大概介绍一下常用的非线性可视化的方法,由于篇幅问题,分了两章。第一篇介绍一下相图与相空间的概念。第二篇介绍一下庞加莱截面和分岔图的绘制方法。后面第三篇就不打算写可视化了,看看要不要写一个频谱分析、李雅普诺夫指数分析、维数分析之类的时间序列分析,但这一部分有趣性不高,实在没动力o(╥﹏╥)o。
非线性的数学研究诞生时间还是比较久了,从蝴蝶效应开始一直到现在,常用的分析方法变化不是太大。各领域研究方向有所不同,但在论文中如果涉及到方程的非线性,总会需要定性的分析一下,由于方法比较类似,所以常常也可以参阅别的领域的非线性分析相关的方法或者好看的绘图方式。
这篇博客主要讲怎么绘制分岔图,顺便也提了一下怎么绘制庞加莱截面。分岔一般指的是系统中某一个常数发生变化,导致整个系统都会随之改变的现象。在电路中可能是某个电容改变,在机械系统中可能是某个轴承摩擦阻力的变化,在流场中可能是风速的改变。这些参数微小的变化,就会导致系统突然性情大变,变得震荡发散不收控制。
上一篇介绍了相图的绘制,这个可以看做是非线性可视化的基础,建议先看一眼上一篇文章:
MATLAB常见非线性可视化绘制方法-相图与相空间https://blog.csdn.net/weixin_42943114/article/details/123193855
参考目录
[1]Morris W. Hirsch. 微分方程、动力系统与混沌导论[M]. 人民邮电出版社, 2008.
[2]刘秉正. 非线性动力学与混沌基础[M]. 东北师范大学出版社, 1994.
[3] Khalil H K . 非线性系统(第3版)[M].电子工业出版社,2005.
[4]Morris W. Hirsch. Differential equations, dynamical systems, and an introduction to chaos [M] Academic Press
[5] 微分方程、动力系统与混沌导论[M]
[6] 计算物理基础-第10章第77讲(北京师范大学)(中国大学MOOC)计算物理基础_北京师范大学_中国大学MOOC(慕课) (icourse163.org)
[7]Computing accurate Poincaré maps[J]. PHYSICA D, 2002, 171(3):127-137.
[8]Santo Banerjee. Applications of Chaos and Nonlinear Dynamics in Engineering[M]. Springer, 2011
[9]Ali.H.Nayfeh. Applied Nonlinear Dynamics Analytical, Computational, and Experimental Methods[M]. 1995.(用到了P284里面的Henon系统)
[10]詹姆斯·格雷克. 混沌:开创新科学[M].
2 离散非线性系统的分岔图绘制
2.1 一维Logistic系统分岔图
Logistic系统是一种非常简单的二次多项式形式的映射。该映射是1976 年生物学家罗伯特·梅 (Robert May) 发扬光大的,用来说明即使很简单的系统也可以产生混沌现象。其迭代方程如下:
x n + 1 = a ( 1 x n ) x n x_{n+1}=a(1-x_n)x_n xn+1=a(1xn)xn
这个迭代方程用来表述某个物种的数量增长变化。其中x的含义为现有物种数量除以理论最大数量,a为增长率。方程认为,物种数量的变化,与物种现存人口x有关,也与人口过多带来的环境压力(1-x)有关。
当增长率a小于1时,物种逐渐灭亡;当增长率在1-3之间,则物种能够稳定在一个范围内;当增长率更大,则会逐渐出现非线性现象。
绘制其分岔图有两种方法:
第一种,均匀分布初始点,然后进行逐次迭代,观察其运动规律。
比如下图选取a=2.1和a=3.3的两个参数,初始分布各个初值,进行迭代,观察这些初值的变化:
可以看到当a=2.1时,最终所有种群收敛到了1个点。当a=3.3时,最终所有种群点收敛到了2个点。
继续细化a的值,列出一个图,可以看到不同参数a变化下,系统的收敛规律:
最终迭代200次的结果如下图:
这种体现随着某个参数(比如a)变化,导致系统行为发生改变的图,叫做分岔图。
当然为了好看,也可以根据点的密度来加一点颜色,比如下面这个样子。
上面几个图的代码如下:
clear
clc
close all
%分岔图
%离散系统
%% 1Logistic系统
%% 1.1Logistic系统 两个典型的a,进行迭代的效果
d=0.004;
x=d:d:1-d;
Nx=length(x);
BF=zeros(1,Nx);
a1=2.1;
a_k1=a1;%第一个a=2.1
a2=3.3;
a_k2=a2;%第二个a=3.3
x1=x;
x2=x;
%在系统中迭代50次
for m=1:50
%把结果绘图
f1=figure(1);
set(f1,'Color',[1,1,1])
clf
subplot(1,2,1)
plot(a1*ones(size(x1)),x1,'.')
ylim([0,1])
text(a1-0.8,0.5,['a1=',num2str(a1)])
text(a1-0.8,0.4,['迭代次数',num2str(m)])
subplot(1,2,2)
plot(a2*ones(size(x1)),x2,'.')
ylim([0,1])
text(a2-0.8,0.5,['a2=',num2str(a2)])
text(a2-0.8,0.4,['迭代次数',num2str(m)])
pause(0.1)
%更新下一次迭代
x1=Logistic(x1,a_k1);
x2=Logistic(x2,a_k2);
end
%% 1.2 不同a,绘制分岔图
%初始种群采用均匀分布
d=0.005;
x=d:d:1-d;
a=0:0.004:4;%把a采集的足够密,就可以绘制随参数a变化的分岔图
Nx=length(x);
Na=length(a);
BF=zeros(Na,Nx);%初始化最终储存的矩阵
for k=1:Na
a_k=a(k);
x1=x;
%在系统中迭代200次
for m=1:200
x1=Logistic(x1,a_k);
end
%把结果保存
BF(k,:)=x1;
end
%画图
figure()
hold on
for k=1:Na
a_k=a(k);
plot(a_k*ones(1,Nx),BF(k,:),...
'LineStyle','none','Marker','.','MarkerFaceColor','k','MarkerEdgeColor','k',...
'MarkerSize',1)
end
hold off
xlabel('a')
%% 1.3 不同a 动图,把上面那个分岔图的过程展示出来
d=0.005;
x=d:d:1-d;
a=0:0.004:4;
Nx=length(x);
Na=length(a);
BF=ones(Na,1)*x;
BFx=a'*ones(1,Nx);
%在系统中迭代100次
for m=1:90
%绘制图片
f3=figure(3);
clf
scatter(BFx(:),BF(:),0.5,'k','MarkerEdgeAlpha',0.5)
text(0.7,0.6,['迭代次数',num2str(m)])
xlabel('a')
ylabel('x')
set(f3,'Color',[1,1,1])
pause(0.1)
%迭代更新
for k=1:Na
a_k=a(k);
x1=BF(k,:);
x1=Logistic(x1,a_k);
%把结果保存
BF(k,:)=x1;
end
end
%% 1.4 不同a上颜色
d=0.002;
x=d:d:1-d;
a=0:0.002:4;
Nx=length(x);
Na=length(a);
BF=zeros(Na,Nx);
for k=1:Na
a_k=a(k);
x1=x;
%在系统中迭代250次
for m=1:250
x1=Logistic(x1,a_k);
end
%把结果保存
BF(