分形几何

一、 问题描述：

从一条直线段开始，将线段中间的三分之一部分用一个等边三角形的两边代替，形成山丘形图形如下

图1

在新的图形中，又将图中每一直线段中间的三分之一部分都用一个等边三角形的两条边代替，再次形成新的图形如此迭代，形成Koch分形曲线。 二、 算法分析：

考虑由直线段（2个点）产生第一个图形（5个点）的过程。图1中，设P和P分别为原始直线段的两个端点，现需要在直线

15段的中间依次插入三个点P，P，P。显然P位于线段三分之一

2342处，P位于线段三分之二处，P点的位置可看成是由P点以P点

4342为轴心，逆时针旋转60而得。旋转由正交矩阵

cos()3Asin()3sin()3cos()30



实现。

算法根据初始数据（P和P点的坐标），产生图1中5个结点的

15坐标。结点的坐标数组形成一个52矩阵，矩阵的第一行为P的

1坐标，第二行为P的坐标……，第五行为P的坐标。矩阵的第一

25列元素分别为5个结点的x坐标，第二列元素分别为5个结点的y坐标。

进一步考虑Koch曲线形成过程中结点数目的变化规律。设第

k次迭代产生的结点数为n，第k1次迭代产生的结点数为n，

kk1kk1k1则n和n中间的递推关系为n三、实验程序及注释：

4nk3。

p=[0 0;10 0]; %P为初始两个点的坐标,第一列为x坐标,第二列为y坐标

n=2; %n为结点数

A=[cos(pi/3) -sin(pi/3);sin(pi/3) cos(pi/3)]; %旋转矩阵 for k=1:4

d=diff(p)/3; %diff计算相邻两个点的坐标之差,得到相邻两点确定的向量

%则d就计算出每个向量长度的三分之一,与题中将线段三等分对应 m=4*n-3; %迭代公式

q=p(1:n-1,:); %以原点为起点,前n-1个点的坐标为终点形成向量

p(5:4:m,:)=p(2:n,:); %迭代后处于4k+1位置上的点的坐

标为迭代前的相应坐标

p(2:4:m,:)=q+d; %用向量方法计算迭代后处于4k+2位置上的点的坐标

p(3:4:m,:)=q+d+d*A'; %用向量方法计算迭代后处于4k+3位置上的点的坐标

p(4:4:m,:)=q+2*d; %用向量方法计算迭代后处于4k位置上的点的坐标

n=m; %迭代后新的结点数目 end

plot(p(:,1),p(:,2)) %绘出每相邻两个点的连线 axis([0 10 0 10]) 四、实验数据记录：

由第三部分的程序，可得到如下的Koch分形曲线：

图2

五、注记：

１．参照实验方法，可绘制如下生成元的Koch 分形曲线：

图3

此时，旋转矩阵为：

cos()2Asin()2sin()021cos()210

程序和曲线如下：

p=[0 0;10 0]; %P为初始两个点的坐标,第一列为x坐标,第二列为y坐标

n=2; %n为结点数 A=[0 -1;1 0]; %旋转矩阵

for k=1:4

d=diff(p)/3; %diff计算相邻两个点的坐标之差,得到相邻两点确定的向量

%则d就计算出每个向量长度的三分之一,与题中将线段三等分对应 m=5*n-4; %迭代公式

q=p(1:n-1,:); %以原点为起点,前n-1个点的坐标为终点形成向量

p(6:5:m,:)=p(2:n,:); %迭代后处于5k+1位置上的点的坐标为迭代前的相应坐标

p(2:5:m,:)=q+d; %用向量方法计算迭代后处于5k+2位置上的点的坐标

p(3:5:m,:)=q+d+d*A'; %用向量方法计算迭代后处于5k+3位置上的点的坐标

p(4:5:m,:)=q+2*d+d*A'; %用向量方法计算迭代后处于5k+4位置上的点的坐标

p(5:5:m,:)=q+2*d; %用向量方法计算迭代后处于5k位置上的点的坐标

n=m; %迭代后新的结点数目 end

plot(p(:,1),p(:,2)) %绘出每相邻两个点的连线

axis([0 10 0 10])

图4

由于中间三分之一部分是一个正方形时，有很多连接的部分。所以我们将高度压缩到原来的0.7倍，即中间部分为一个长与宽之比为1:0.7的矩形时，得到程序和曲线如下：

p=[0 0;10 0]; %P为初始两个点的坐标,第一列为x坐标,第二列为y坐标

n=2; %n为结点数 A=[0 -1;1 0]; %旋转矩阵 for k=1:4

d=diff(p)/3; %diff计算相邻两个点的坐标之差,得到相邻两点确定的向量

%则d就计算出每个向量长度的三分之一,与题中将线段三等分对应

m=5*n-4; %迭代公式

q=p(1:n-1,:); %以原点为起点,前n-1个点的坐标为终点形成向量

p(6:5:m,:)=p(2:n,:); %迭代后处于5k+1位置上的点的坐标为迭代前的相应坐标

p(2:5:m,:)=q+d; %用向量方法计算迭代后处于5k+2位置上的点的坐标

p(3:5:m,:)=q+d+0.7*d*A'; %用向量方法计算迭代后处于5k+3位置上的点的坐标

p(4:5:m,:)=q+2*d+0.7*d*A'; %用向量方法计算迭代后处于5k+4位置上的点的坐标

p(5:5:m,:)=q+2*d; %用向量方法计算迭代后处于5k位置上的点的坐标

n=m; %迭代后新的结点数目 end

plot(p(:,1),p(:,2)) %绘出每相邻两个点的连线 axis([0 10 0 10])

图５

２．参照实验方法，我们由四边形的四个初始点出发，对于四边形的每条边，生成元如下：

图6

可得到火焰般的图形。 程序和曲线如下：

p=[0 10;10 0;0 -10;-10 0;0 10];

%P为四边形四个顶点的坐标,其中第五个点与第一个点重合,以

便于绘图

%第一列为x坐标,第二列为y坐标 n=5; %n为结点数

A=[cos(-pi/3) -sin(-pi/3);sin(-pi/3) cos(-pi/3)]; %旋转矩阵,顺时针旋转60度 for k=1:5

d=diff(p)/3;m=4*n-3; %迭代公式 q=p(1:n-1,:); p(5:4:m,:)=p(2:n,:); p(2:4:m,:)=q+d; p(3:4:m,:)=q+2*d+d*A'; p(4:4:m,:)=q+2*d; n=m; end

plot(p(:,1),p(:,2)) axis([-10 10 -10 10])

图7

３．参照实验方法，由下列的生成元，绘制Koch分形曲线：

图8

分析：为了绘图方便，我们将结点数处理一下，把第一次迭代产生的六个点看成十个点，即图中有五条线段（１－２，３－４，５－６，７－８，９－１０），我们将每条线段的每个端点看成新的两个结点，这样我们就可以很方便地用plot绘图了。 程序和曲线如下：

p=[0 0;10 10]; %P为初始两个点的坐标,第一列为x坐标,第二列

为y坐标

n=2; %n为结点数

A=[cos(pi/3) -sin(pi/3);sin(pi/3) cos(pi/3)]; B=[cos(-pi/3) -sin(-pi/3);sin(-pi/3) cos(-pi/3)]; %旋转矩阵A对应于第一次逆时针旋转60度,旋转矩阵B对应于第二次顺时针旋转60度 for k=1:4 d=diff(p)/3;

d1=d(1:2:n,:);%取每条线段对应的向量 m=5*n; %迭代公式 q1=p(1:2:n-1,:); p(10:10:m,:)=p(2:2:n,:);

p(1:10:m,:)=p(1:2:n,:); %迭代后处于10k与10k+1位置上的点的坐标为迭代前的相应坐标 p(2:10:m,:)=q1+d1;

%用向量方法计算迭代后处于10k+2,10k+3,10k+5位置上的点的坐标,都相同

p(3:10:m,:)=p(2:10:m,:);

p(4:10:m,:)=q1+d1+d1*A'; %用向量方法计算迭代后处于10k+4位置上的点的坐标 p(5:10:m,:)=p(2:10:m,:); p(6:10:m,:)=q1+2*d1;

%用向量方法计算迭代后处于10k+6,10k+7,10k+9位置上的点的坐标,都相同

p(7:10:m,:)=p(6:10:m,:); p(8:10:m,:)=q1+2*d1+d1*B'; p(9:10:m,:)=p(6:10:m,:);

n=m; %迭代后新的结点数目 end

plot(p(:,1),p(:,2)) %绘出每相邻两个点的连线 axis([0 10 0 10])

另一种方法、function koch(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); for i=1:n

for j=1:length(A)/5;

w=sub_koch2(A(1+5*(j-1):5*j),B(1+5*(j-1):5*j)); for k=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));

end end A=AA; B=BB; end plot(A,B) hold on axis equal

%由以(ax,ay),(bx,by)为端点的线段生成新的中间三点坐标并把这五点横、纵坐标依次分别存%储在数组A,B中

function [A,B]=sub_koch1(ax,ay,bx,by) cx=ax+(bx-ax)/3; cy=ay+(by-ay)/3; ex=bx-(bx-ax)/3; ey=by-(by-ay)/3;

L=sqrt((ex-cx).^2+(ey-cy).^2); alpha=atan((ey-cy)./(ex-cx)); if (ex-cx)<0 alpha=alpha+pi; end

dx=cx+cos(alpha+pi/3)*L; dy=cy+sin(alpha+pi/3)*L; A=[ax,cx,dx,ex,bx]; B=[ay,cy,dy,ey,by];

%把由函数sub_koch1生成的五点横、纵坐标A,B顺次划分为四组,分别对应四条折线段中

%每条线段两端点的坐标,并依次分别存储在4*4阶矩阵k中,k中第i(i=1,2,3,4)行数字代表第%i条线段两端点的坐标 function w=sub_koch2(A,B) a11=A(1);b11=B(1); a12=A(2);b12=B(2);

a21=A(2);b21=B(2); a22=A(3);b22=B(3); a31=A(3);b31=B(3); a32=A(4);b32=B(4); a41=A(4);b41=B(4); a42=A(5);b42=B(5);

w=[a11,b11,a12,b12;a21,b21,a22,b22;a31,b31,a32,b32;a41,b41,a42,b42];

图1 Von Koch曲线

（2）、Levy 曲线程序levy.m function levy(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); for i=1:n

for j=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)

hold on axis equal

%由以(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条线段%两端点的坐标

function w=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];

图2 Levy 曲线

（3） 分形树程序tree.h function tree(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]) hold on

[X,Y]=tree1(x1,y1,x2,y2,a,b); hold on 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]; for k=1:n

for i=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; else if (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) axis equal hold on X=[x2,x3,x4]; Y=[y2,y3,y4];

%把由函数tree1生成的X,Y顺次划分为两组,分别对应两分枝两

个端点的坐标,并存储在一维%数组w中 function w=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 function sierpinski_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.5 0 0 0 0.5 0]; M2=[0.5 0 0.5 0 0.5 0]; M3=[0.5 0 0.25 0 0.5 0.5]; x=0;y=0;

% r为[0,1]区间内产生的n维随机数组 r=rand(1,n); B=zeros(2,n); k=1;

% 当0a=M1(1);b=M1(2);e=M1(3);c=M1(4);d=M1(5);f=M1(6); else if r(i)a=M2(1);b=M2(2);e=M2(3);c=M2(4);d=M2(5);f=M2(6); else if r(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)

图4 Sierpinski三角形

（5）IFS算法画Julia集程序julia_ifs.h function julia_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); for i=1:n wx=x-cx; wy=y-cy; if wx>0

alpha=atan(wy/wx); end if wx<0

alpha=pi+atan(wy/wx);

end if wx==0 alpha=pi/2; end

alpha=alpha/2; r=sqrt(wx^2+wy^2); if A(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)

图5 Julia 集

(6)

逃逸时间算法画Sierpinski垫片程序sierpinski.h

function sierpinski(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; for i=1:m

x0=a+(c-a)*(i-1)/m;

for j=1:m

y0=b+(d-b)*(j-1)/m; x=x0; y=y0; for k=1:n if y>0.5 x=2*x; y=2*y-1; else if x>=0.5 x=2*x-1; y=2*y; else x=2*x; y=2*y; end

if x^2+y^2>r break; end end if k==n B(1,w)=i; B(2,w)=j;

w=w+1; end end end

plot(B(1,:),B(2,:),'.','markersize',0.1)

图6 Sierpinski三角形垫片

(7)

元胞自动机算法画Sierpinski三角形程序

 一维元胞自动机sierpinski_ca1.h function sierpinski_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; for i=1:m-1 for j=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 function sierpinski_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;

for i=[m/2:-1:2,m/2:m-1] for j=[n/2:-1:2,n/2:n-1] if

mod(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

function helix_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.424242 1.758647 0.242424 0.859848 1.408065];

M2=[-0.121212 0.257576 -6.721654 0.05303 0.05303 1.377236];

M3=[0.181818 -0.136364 6.086107 0.090909 0.181818 1.568035]; x=0;y=0;

% r为[0,1]区间内产生的n维随机数组 r=rand(1,n); B=zeros(2,n); k=1;

% 当0a=M1(1);b=M1(2);e=M1(3);c=M1(4);d=M1(5);f=M1(6);

else if r(i)a=M2(1);b=M2(2);e=M2(3);c=M2(4);d=M2(5);f=M2(6); else if r(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)

图8 Helix曲线

三、 问题描述：

从一条直线段开始，将线段中间的三分之一部分用一个等边三角形的两边代替，形成山丘形图形如下

图1

在新的图形中，又将图中每一直线段中间的三分之一部分都用一个等边三角形的两条边代替，再次形成新的图形如此迭代，形成Koch分形曲线。 四、 算法分析：

考虑由直线段（2个点）产生第一个图形（5个点）的过程。图1中，设P和P分别为原始直线段的两个端点，现需要在直线

15段的中间依次插入三个点P，P，P。显然P位于线段三分之一

2342处，P位于线段三分之二处，P点的位置可看成是由P点以P点

4342为轴心，逆时针旋转60而得。旋转由正交矩阵

cos()3Asin()3sin()3cos()30



实现。

算法根据初始数据（P和P点的坐标），产生图1中5个结点的

15坐标。结点的坐标数组形成一个52矩阵，矩阵的第一行为P的

1坐标，第二行为P的坐标……，第五行为P的坐标。矩阵的第一

25列元素分别为5个结点的x坐标，第二列元素分别为5个结点的y坐标。

进一步考虑Koch曲线形成过程中结点数目的变化规律。设第

k次迭代产生的结点数为n，第k1次迭代产生的结点数为n，

kk1kk1k1则n和n中间的递推关系为n三、实验程序及注释：

4nk3。

p=[0 0;10 0]; %P为初始两个点的坐标,第一列为x坐标,第二列

为y坐标

n=2; %n为结点数

A=[cos(pi/3) -sin(pi/3);sin(pi/3) cos(pi/3)]; %旋转矩阵 for k=1:4

d=diff(p)/3; %diff计算相邻两个点的坐标之差,得到相邻两点确定的向量

%则d就计算出每个向量长度的三分之一,与题中将线段三等分对应 m=4*n-3; %迭代公式

q=p(1:n-1,:); %以原点为起点,前n-1个点的坐标为终点形成向量

p(5:4:m,:)=p(2:n,:); %迭代后处于4k+1位置上的点的坐标为迭代前的相应坐标

p(2:4:m,:)=q+d; %用向量方法计算迭代后处于4k+2位置上的点的坐标

p(3:4:m,:)=q+d+d*A'; %用向量方法计算迭代后处于4k+3位置上的点的坐标

p(4:4:m,:)=q+2*d; %用向量方法计算迭代后处于4k位置上的点的坐标

n=m; %迭代后新的结点数目 end

plot(p(:,1),p(:,2)) %绘出每相邻两个点的连线 axis([0 10 0 10]) 四、实验数据记录：

由第三部分的程序，可得到如下的Koch分形曲线：

图2

五、注记：

１．参照实验方法，可绘制如下生成元的Koch 分形曲线：

图3

此时，旋转矩阵为：

cos()2Asin()2sin()021cos()210

程序和曲线如下：

p=[0 0;10 0]; %P为初始两个点的坐标,第一列为x坐标,第二列为y坐标

n=2; %n为结点数 A=[0 -1;1 0]; %旋转矩阵 for k=1:4

d=diff(p)/3; %diff计算相邻两个点的坐标之差,得到相邻两点确定的向量

%则d就计算出每个向量长度的三分之一,与题中将线段三等分对应 m=5*n-4; %迭代公式

q=p(1:n-1,:); %以原点为起点,前n-1个点的坐标为终点形成向量

p(6:5:m,:)=p(2:n,:); %迭代后处于5k+1位置上的点的坐标为迭代前的相应坐标

p(2:5:m,:)=q+d; %用向量方法计算迭代后处于5k+2位置上的点的坐标

p(3:5:m,:)=q+d+d*A'; %用向量方法计算迭代后处于5k+3位置上的点的坐标

p(4:5:m,:)=q+2*d+d*A'; %用向量方法计算迭代后处于5k+4位置上的点的坐标

p(5:5:m,:)=q+2*d; %用向量方法计算迭代后处于5k

位置上的点的坐标

n=m; %迭代后新的结点数目 end

plot(p(:,1),p(:,2)) %绘出每相邻两个点的连线 axis([0 10 0 10])

图4

由于中间三分之一部分是一个正方形时，有很多连接的部分。所以我们将高度压缩到原来的0.7倍，即中间部分为一个长与宽之比为1:0.7的矩形时，得到程序和曲线如下：

p=[0 0;10 0]; %P为初始两个点的坐标,第一列为x坐标,第二列为y坐标

n=2; %n为结点数 A=[0 -1;1 0]; %旋转矩阵

for k=1:4

d=diff(p)/3; %diff计算相邻两个点的坐标之差,得到相邻两点确定的向量

%则d就计算出每个向量长度的三分之一,与题中将线段三等分对应 m=5*n-4; %迭代公式

q=p(1:n-1,:); %以原点为起点,前n-1个点的坐标为终点形成向量

p(6:5:m,:)=p(2:n,:); %迭代后处于5k+1位置上的点的坐标为迭代前的相应坐标

p(2:5:m,:)=q+d; %用向量方法计算迭代后处于5k+2位置上的点的坐标

p(3:5:m,:)=q+d+0.7*d*A'; %用向量方法计算迭代后处于5k+3位置上的点的坐标

p(4:5:m,:)=q+2*d+0.7*d*A'; %用向量方法计算迭代后处于5k+4位置上的点的坐标

p(5:5:m,:)=q+2*d; %用向量方法计算迭代后处于5k位置上的点的坐标

n=m; %迭代后新的结点数目 end

plot(p(:,1),p(:,2)) %绘出每相邻两个点的连线

axis([0 10 0 10])

图５

２．参照实验方法，我们由四边形的四个初始点出发，对于四边形的每条边，生成元如下：

图6

可得到火焰般的图形。 程序和曲线如下：

p=[0 10;10 0;0 -10;-10 0;0 10];

%P为四边形四个顶点的坐标,其中第五个点与第一个点重合,以便于绘图

%第一列为x坐标,第二列为y坐标 n=5; %n为结点数

A=[cos(-pi/3) -sin(-pi/3);sin(-pi/3) cos(-pi/3)]; %旋转矩阵,顺时针旋转60度 for k=1:5

d=diff(p)/3;m=4*n-3; %迭代公式 q=p(1:n-1,:); p(5:4:m,:)=p(2:n,:); p(2:4:m,:)=q+d; p(3:4:m,:)=q+2*d+d*A'; p(4:4:m,:)=q+2*d; n=m; end

plot(p(:,1),p(:,2)) axis([-10 10 -10 10])

图7

３．参照实验方法，由下列的生成元，绘制Koch分形曲线：

图8

分析：为了绘图方便，我们将结点数处理一下，把第一次迭代产生的六个点看成十个点，即图中有五条线段（１－２，３－４，５－６，７－８，９－１０），我们将每条线段的每个端点看成新的两个结点，这样我们就可以很方便地用plot绘图了。 程序和曲线如下：

p=[0 0;10 10]; %P为初始两个点的坐标,第一列为x坐标,第二列

为y坐标

n=2; %n为结点数

A=[cos(pi/3) -sin(pi/3);sin(pi/3) cos(pi/3)]; B=[cos(-pi/3) -sin(-pi/3);sin(-pi/3) cos(-pi/3)]; %旋转矩阵A对应于第一次逆时针旋转60度,旋转矩阵B对应于第二次顺时针旋转60度 for k=1:4 d=diff(p)/3;

d1=d(1:2:n,:);%取每条线段对应的向量 m=5*n; %迭代公式 q1=p(1:2:n-1,:); p(10:10:m,:)=p(2:2:n,:);

p(1:10:m,:)=p(1:2:n,:); %迭代后处于10k与10k+1位置上的点的坐标为迭代前的相应坐标 p(2:10:m,:)=q1+d1;

%用向量方法计算迭代后处于10k+2,10k+3,10k+5位置上的点的坐标,都相同

p(3:10:m,:)=p(2:10:m,:);

p(4:10:m,:)=q1+d1+d1*A'; %用向量方法计算迭代后处于10k+4位置上的点的坐标 p(5:10:m,:)=p(2:10:m,:); p(6:10:m,:)=q1+2*d1;

%用向量方法计算迭代后处于10k+6,10k+7,10k+9位置上的点的坐标,都相同

p(7:10:m,:)=p(6:10:m,:); p(8:10:m,:)=q1+2*d1+d1*B'; p(9:10:m,:)=p(6:10:m,:);

n=m; %迭代后新的结点数目 end

plot(p(:,1),p(:,2)) %绘出每相邻两个点的连线 axis([0 10 0 10])

六，结束语

通过图形显示，更好的理解的分形同时也也加深对分形概

念的进一步掌握

参考文献：<>