分形几何中一些经典图形的Matlab画法.docx
- 文档编号:27686119
- 上传时间:2023-07-04
- 格式:DOCX
- 页数:19
- 大小:180.21KB
分形几何中一些经典图形的Matlab画法.docx
《分形几何中一些经典图形的Matlab画法.docx》由会员分享,可在线阅读,更多相关《分形几何中一些经典图形的Matlab画法.docx(19页珍藏版)》请在冰豆网上搜索。
分形几何中一些经典图形的Matlab画法
分形几何中一些经典图形的Matlab画法
分形几何中一些经典图形的Matlab画法
(1)Koch曲线程序koch.m
functionkoch(a1,b1,a2,b2,n)
%koch(0,0,9,0,3)
%a1,b1,a2,b2为初始线段两端点坐标,n为迭代次数
a1=0;b1=0;a2=9;b2=0;n=3;
%第i-1次迭代时由各条线段产生的新四条线段的五点横、纵坐标存储在数组A、B中
[A,B]=sub_koch1(a1,b1,a2,b2);
fori=1:
n
forj=1:
length(A)/5;
w=sub_koch2(A(1+5*(j-1):
5*j),B(1+5*(j-1):
5*j));
fork=1:
4
[AA(5*4*(j-1)+5*(k-1)+1:
5*4*(j-1)+5*(k-1)+5),BB(5*4*(j-1)+5*(k-1)+1:
5*4*(j-1)+5*(k-1)+5)]=sub_koch1(w(k,1),w(k,2),w(k,3),w(k,4));
图1VonKoch曲线
(2)Levy曲线程序levy.m
functionlevy(n)
%levy(16),n为levy曲线迭代次数
%x1,y1,x2,y2为初始线段两端点坐标,nn为迭代次数
n=16;
x1=0;y1=0;
x2=1;y2=0;
%第i-1次迭代时由各条线段产生的新两条线段的三端点横、纵坐标存储在数组X、Y中
[X,Y]=levy1(x1,y1,x2,y2);
fori=1:
n
forj=1:
length(X)/3
w=levy2(X(1+3*(j-1):
3*j),Y(1+3*(j-1):
3*j));
[XX(3*2*(j-1)+1:
3*2*(j-1)+3),YY(3*2*(j-1)+1:
3*2*(j-1)+3)]=levy1(w(1,1),w(1,2),w(1,3),w(1,4));
[XX(3*2*(j-1)+3+1:
3*2*(j-1)+3+3),YY(3*2*(j-1)+3+1:
3*2*(j-1)+3+3)]=levy1(w(2,1),w(2,2),w(2,3),w(2,4));
end
X=XX;
Y=YY;
end
plot(X,Y)
holdon
axisequal
%由以(x1,y1),(x2,y2)为端点的线段生成新的中间点坐标并把(x1,y1),(x2,y2)连同新点横、纵坐%标依次分别存储在数组X,Y中
function[X,Y]=levy1(x1,y1,x2,y2)
x3=1/2*(x1+x2+y1-y2);
y3=1/2*(-x1+x2+y1+y2);
X=[x1,x3,x2];
Y=[y1,y3,y2];
%把由函数levy1生成的三点横、纵坐标X,Y顺次划分为两组,分别对应两条折线段中每条线%段两端点的坐标,并依次分别存储在2*4阶矩阵w中,w中第i(i=1,2)行数字代表第i条线段%两端点的坐标
functionw=levy2(X,Y)
a11=X
(1);b11=Y
(1);
a12=X
(2);b12=Y
(2);
a21=X
(2);b21=Y
(2);
a22=X(3);b22=Y(3);
w=[a11,b11,a12,b12;a21,b21,a22,b22];
图2Levy曲线
(3)分形树程序tree.h
functiontree(n,a,b)
%tree(8,pi/8,pi/8),n为分形树迭代次数
%a,b为分枝与竖直方向夹角
%x1,y1,x2,y2为初始线段两端点坐标,nn为迭代次数
n=8;a=pi/8;b=pi/8;
x1=0;y1=0;
x2=0;y2=1;
plot([x1,x2],[y1,y2])
holdon
[X,Y]=tree1(x1,y1,x2,y2,a,b);
holdon
W=tree2(X,Y);
w1=W(:
1:
4);
w2=W(:
5:
8);
%w为2^k*4维矩阵,存储第k次迭代产生的分枝两端点的坐标,
%w的第i(i=1,2,…,2^k)行数字对应第i个分枝两端点的坐标
w=[w1;w2];
fork=1:
n
fori=1:
2^k
[X,Y]=tree1(w(i,1),w(i,2),w(i,3),w(i,4),a,b);
W(i,:
)=tree2(X,Y);
end
w1=W(:
1:
4);
w2=W(:
5:
8);
w=[w1;w2];
end
%由每个分枝两端点坐标(x1,y1),(x2,y2)产生两新点的坐标(x3,y3),(x4,y4),画两分枝图形,并把%(x2,y2)连同新点横、纵坐标分别存储在数组X,Y中
function[X,Y]=tree1(x1,y1,x2,y2,a,b)
L=sqrt((x2-x1)^2+(y2-y1)^2);
if(x2-x1)==0
a=pi/2;
elseif(x2-x1)<0
a=pi+atan((y2-y1)/(x2-x1));
else
a=atan((y2-y1)/(x2-x1));
end
end
x3=x2+L*2/3*cos(a+b);
y3=y2+L*2/3*sin(a+b);
x4=x2+L*2/3*cos(a-b);
y4=y2+L*2/3*sin(a-b);
a=[x3,x2,x4];
b=[y3,y2,y4];
plot(a,b)
axisequal
holdon
X=[x2,x3,x4];
Y=[y2,y3,y4];
%把由函数tree1生成的X,Y顺次划分为两组,分别对应两分枝两个端点的坐标,并存储在一维%数组w中
functionw=tree2(X,Y)
a1=X
(1);b1=Y
(1);
a2=X
(2);b2=Y
(2);
a3=X
(1);b3=Y
(1);
a4=X(3);b4=Y(3);
w=[a1,b1,a2,b2,a3,b3,a4,b4];
图3分形树
(4)IFS算法画Sierpinski三角形程序sierpinski_ifs.h
functionsierpinski_ifs(n,w1,w2,w3)
%sierpinski_ifs(10000,1/3,1/3,1/3)
%w1,w2,w3出现频率
n=10000;
w1=1/3;
w2=1/3;
w3=1/3;
M1=[0.50000.50];
M2=[0.500.500.50];
M3=[0.500.2500.50.5];
x=0;y=0;
%r为[0,1]区间内产生的n维随机数组
r=rand(1,n);
B=zeros(2,n);
k=1;
%当0 %当1/3= %当2/3= fori=1: n ifr(i) a=M1 (1);b=M1 (2);e=M1(3);c=M1(4);d=M1(5);f=M1(6); elseifr(i) a=M2 (1);b=M2 (2);e=M2(3);c=M2(4);d=M2(5);f=M2(6); elseifr(i) a=M3 (1);b=M3 (2);e=M3(3);c=M3(4);d=M3(5);f=M3(6); end end end x=a*x+b*y+e; y=c*x+d*y+f; B(1,k)=x; B(2,k)=y; k=k+1; end plot(B(1,: ),B(2,: ),'.','markersize',0.1) 图4Sierpinski三角形 (5)IFS算法画Julia集程序julia_ifs.h functionjulia_ifs(n,cx,cy) %julia_ifs(100000,-0.77,0.08) %f(z)=z^2+c,cx=real(c);cy=image(c); n=10000; cx=-0.77; cy=0.08; %z^2+c=z0,x=real(z0);y=image(z0); x=1;y=1; B=zeros(2,n); k=1; %A为产生的服从标准正态分布的n维随机数组 A=randn(1,n); fori=1: n wx=x-cx; wy=y-cy; ifwx>0 alpha=atan(wy/wx); end ifwx<0 alpha=pi+atan(wy/wx); end ifwx==0 alpha=pi/2; end alpha=alpha/2; r=sqrt(wx^2+wy^2); ifA(i)<0 r=-sqrt(r); else r=sqrt(r); end x=r*cos(alpha); y=r*sin(alpha); B(1,k)=x; B(2,k)=y; k=k+1; end plot(B(1,: ),B(2,: ),'.','markersize',0.1) 图5Julia集 (6)逃逸时间算法画Sierpinski垫片程序sierpinski.h functionsierpinski(a,b,c,d,n,m,r) %sierpinski(0,0,1,1,12,200,200) %(a,b),(c,d)收敛区域左上角和右下角坐标,m为分辨率 %n为逃逸时间,需要反复试探,r逃逸半径 a=0;b=0;c=1;d=1;n=12;m=200;r=200; B=zeros(2,m*m); w=1; fori=1: m x0=a+(c-a)*(i-1)/m; forj=1: m y0=b+(d-b)*(j-1)/m; x=x0; y=y0; fork=1: n ify>0.5 x=2*x; y=2*y-1; elseifx>=0.5 x=2*x-1; y=2*y; else x=2*x; y=2*y; end ifx^2+y^2>r break; end end ifk==n B(1,w)=i; B(2,w)=j; w=w+1; end end end plot(B(1,: ),B(2,: ),'.','markersize',0.1) 图6Sierpinski三角形垫片 (7)元胞自动机算法画Sierpinski三角形程序 ✧一维元胞自动机sierpinski_ca1.h functionsierpinski_ca1(m,n) %sierpinski_ca1(1000,3000) m=1000;n=3000; x=1;y=1; t=1;w=zeros(2,m*n); s=zeros(m,n); s(1,fix(n/3))=1; fori=1: m-1 forj=2: n-1 if(s(i,j-1)==1&s(i,j)==0&s(i,j+1)==0)|(s(i,j-1)==0&s(i,j)==0&s(i,j+1)==1) s(i+1,j)=1; w(1,t)=x+3+3*j; w(2,t)=y+5*i; t=t+1; end end end plot(w(1,: ),w(2,: ),'.','markersize',1) 图7.1一维元胞自动机画Sierpinski三角形 ✧二维元胞自动机sierpinski_ca2.h functionsierpinski_ca2(m,n) %sierpinski_ca2(400,400) m=400;n=400; t=1;w=zeros(2,m*n); s=zeros(m,n); s(m/2,n/2)=1; fori=[m/2: -1: 2,m/2: m-1] forj=[n/2: -1: 2,n/2: n-1] ifmod(s(i-1,j-1)+s(i,j-1)+s(i+1,j-1)+s(i-1,j)+s(i+1,j)+s(i-1,j+1)+s(i,j+1)+s(i+1,j+1),2)==1 s(i,j)=1; w(1,t)=i; w(2,t)=j; t=t+1; end end end plot(w(1,: ),w(2,: ),'.','markersize',0.1) 图7.2二维元胞自动机画Sierpinski三角形 (8)IFS算法画Helix曲线程序helix_ifs.h functionhelix_ifs(n,w1,w2,w3) %helix_ifs(20000,0.9,0.05,0.05) %w1,w2,w3为出现频率 n=20000;w1=0.9;w2=0.05;w3=0.05; M1=[0.787879-0.4242421.7586470.2424240.8598481.408065]; M2=[-0.1212120.257576-6.7216540.053030.053031.377236]; M3=[0.181818-0.1363646.0861070.0909090.1818181.568035]; x=0;y=0; %r为[0,1]区间内产生的n维随机数组 r=rand(1,n); B=zeros(2,n); k=1; %当0 %当1/3= %当2/3= fori=1: n ifr(i) a=M1 (1);b=M1 (2);e=M1(3);c=M1(4);d=M1(5);f=M1(6); elseifr(i) a=M2 (1);b=M2 (2);e=M2(3);c=M2(4);d=M2(5);f=M2(6); elseifr(i) a=M3 (1);b=M3 (2);e=M3(3);c=M3(4);d=M3(5);f=M3(6); end end end x=a*x+b*y+e; y=c*x+d*y+f; B(1,k)=x; B(2,k)=y; k=k+1; end plot(B(1,: ),B(2,: ),'.','markersize',0.1) 图8Helix曲线
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 几何 一些 经典 图形 Matlab 画法