非线性动力学学习报告
- 格式：doc
- 大小：377.50 KB
- 文档页数：15
非线性动力学学习报告
在课堂上老师以生动活泼的方式介绍了分形的相关知识，特别是展现了一些美丽的分形图案，我对此十分感兴趣，所以课后找了一些相关资料，学会了用仿射变换的循环迭代方法，在MATLAB平台下，实现了一些简单的飞行图案的绘制。具体内容见项目一。其中的数学原理由于我还不是特别清楚，所以在此进仅做一简要汇报，下面会具体叙述用MATLAB绘制分形图案的过程。
在项目二中，探讨了对于一根细长压杆，端部的压力大小与杆件变形之间的关系。这里的端部压力是较大的载荷（即大于临界力），那么经典的材料力学理论便束手无策，这里构建了一个压杆变形的微段迭代模型，把一个大变形非线性问题转化为有限个小变形的迭加，用MATLAB编程迭代计算的结果较好的吻合了铁木辛哥弹性稳定理论中有关压杆弹性屈曲中的一些成果。
项目一：用MATLAB绘制美丽的分形图案
上个世纪60年代，B.Mandelbrot对一个具有复杂几何性质但局部看起来仍然一样的几何对象提出了分形概念。在很多非线性动力学系统等血多领域都会看到分形的例子，随着电子计算机的发展，我们绘制出了很多分形图案。
在这个项目中，实现了用MATLAB来绘制蕨类植物枝叶和著名的Sierpinski三角形；另外还给出了一个通过编程绘制树枝的例子没有用到仿射变换，只是复杂的循环。
经过翻阅相关资料（考文献[1]），我了解到数学中的仿射变换的定义如下：设x是一个n维向量，A是n*n的矩阵，b是与x同维的向量，那么变换bAxx称作仿射变换，去不同的A，b就会得到不同的变换结果。如果打印前k次（k应该取较大的值）迭代过程中向量x在坐标系中所表示的所有点，那么就可以得到一幅漂亮的分形图案。其中矩阵A和向量b的取法涉及到很复杂的数学理论，在这里不做详细介绍。
基于前面的理论分析很容易得到MATLAB绘图程序代码及其运行结果。
1.、使用数学中的仿射变换理论，绘制蕨类植物枝叶
程序：%fenxing_juelei.m
%蕨类植物模拟
x = [.5; .5]; %初值
h = plot(x(1),x(2),'.'); %绘制初值点
%设置用于后面随机数的判别向量
p = [ .85 .92 .99 1.00];
b1 = [0; 1.6];
b2 = [0; 1.6];
b3 = [0; .44];
%------仿射变换矩阵
A1 = [ .85 .04; -.04 .85];
A2 = [ .20 -.26; .23 .22];
A3 = [-.15 .28; .26 .24];
A4 = [ 0 0 ; 0 .16];
for i=1:20000
r = rand; %产生随机数
if r < p(1)
x = A1*x + b1;
elseif r < p(2)
x = A2*x + b2;
elseif r < p(3)
x = A3*x + b3;
else
x = A4*x;
end
plot(x(1),x(2),'g'),hold on %采用绿色绘制
end
axis off %取消坐标轴
把该m文件放置到Matlab的当前工作目录下在命令行中输入fenxing_juelei，变得到了下面的运行结果。
在具体操作中，迭代了20000次，电脑计算了十多秒钟才计算绘制完毕，着看起来好像迭代次数，可是得到的图形靠下的部分还是不够细密。而实际上再多尺度上分型图形保持了其固有的形态，就是说当把分形图形的一部分无穷放大后仍然保留了其原来的形状特性。如果增加迭代次数显然可以得到更为细密的蕨类植物枝叶图形。
2.、绘制Sierpinski三角形分形图案
程序：fenxing_sie.m
%Sierpinski三角形分形图案
%初值
x = [0; 0];
h = plot(x(1),x(2),'.'); %绘制初值
A=[1/2 0;0 1/2]; %仿射变换矩阵
b1=[0 0]';
b2=[1/2 0]';
b3=[1/4 sqrt(3)/4]';
for i=1:10000
r = rand; %产生随机数
if r >0.66666666667
x = A*x + b1;
elseif r >0.333333333333
x=A*x+b2;
else
x = A*x+b3;
end
plot(x(1),x(2),'-'),hold on %保持图形
end
axis off %取消坐标轴
程序编制完成后将其放在MATLAB当前工作目录下，在命令行中输入fenxing_sie。得到下面的运行结果。
3.、编程绘制树枝图案的形成过程。
○1绘图原理：
树形分形图，将一条线段三等分，在等分点上个画一条长度为原线段长度的1/3长的线段，该线段与原线段成一定夹角。依次做下去就会得到树形。
○2程序代码%tree.m
clf;
clear;
theta=pi/6;
u(:,:,1)=[0,0;0,1];
n=1;
for k=1:6;
subplot(2,3,k)
hold on;
o=u;
for i=0:n-1;
diff=(o(:,2,i+1)-o(:,1,i+1))/3;
u(:,1,5*i+1)=o(:,1,i+1);
u(:,2,5*i+1)=u(:,1,5*i+1)+diff;
u(:,1,5*i+2)=u(:,2,5*i+1);
lp=[cos(theta),-sin(theta);sin(theta),cos(theta)]*diff;
u(:,2,5*i+2)=u(:,1,5*i+2)+lp;
u(:,1,5*i+3)=u(:,2,5*i+1);
u(:,2,5*i+3)=u(:,2,5*i+1)+diff;
u(:,1,5*i+4)=u(:,2,5*i+3);
rp=[cos(theta),sin(theta);-sin(theta),cos(theta)]*diff;
u(:,2,5*i+4)=u(:,1,5*i+4)+rp;
u(:,1,5*i+5)=u(:,1,5*i+4);
u(:,2,5*i+5)=u(:,1,5*i+5)+diff;
end
n=length(u);
for j=1:n;
plot(u(1,:,j),u(2,:,j));
end
axis off
end
○3程序运行结果
4.、编程绘制雪花图案的形成过程。
○1绘图原理：
雪花的分形图，将等边三角形的一条边三等分，将中间的一段去掉，代之以一个更小的等边三角形的两条边，依次做下去，将得到雪花图。
○2程序简单分析
在程序中用复数的实部和虚部分别表示点的坐标x,y。plot是作图指令，它对复数作图相当于用复数的实部和虚部分别表示曲线的x,y坐标。因此只要计算出各点的坐标，依次连接起来，就是所求的结果。
程序中的第一句是列出第一个分区图中a,b,c,d共4点的坐标，其中a,d两点是重合的，然后画出一个等边三角形，接下来对k的循环画出了后面的5各分区图，而对jj的循环，是对每条边计算新图形中4个点的坐标。
○3程序代码%xuehua.m
clf;
clear;
new=[.5+(sqrt(3)/2*i),-0.5+(sqrt(3)/2*i),0+0*i,0.5+(sqrt(3)/2*i)];%%给分区图1中a,b,c,d点的坐标
subplot(2,3,1);plot(new),axis equal%%画出分区图2
for k=1:5%%画出分区图2~6
old=new;%%保存原有个点的坐标
n=length(old)-1;%%计算点的数目
subplot(2,3,k+1);%
axis([-1 1 -0.5 1.5])%
grid%
hold on%
for j=1:n+1%
plot(old(j),'r.','marker','.','markersize',5)%
end%
for jj=0:n-1; %%对每条边计算4个新的点的坐标
diff=(old(jj+2)-old(jj+1))/3;
new(4*jj+1)=old(jj+1);
plot(new(4*jj+1),'marker','*','markersize',5)%
new(4*jj+2)=old(jj+1)+diff;
plot(new(4*jj+2),'marker','*','markersize',5)%
new(4*jj+3)=new(4*jj+2)+diff*((1-sqrt(3)*i)/2);
plot(new(4*jj+3),'marker','*','markersize',5)%
new(4*jj+4)=old(jj+1)+2*diff;subplot(2,3,k+1);
plot(new(4*jj+4),'marker','*','markersize',5)%
end
new(4*n+1)=old(n+1);
subplot(2,3,k+1);
plot(new),axis equal;
end
○4程序运行结果：
