分形几何中一些经典图形的Matlab画法说课材料Word下载.docx
《分形几何中一些经典图形的Matlab画法说课材料Word下载.docx》由会员分享，可在线阅读，更多相关《分形几何中一些经典图形的Matlab画法说课材料Word下载.docx（17页珍藏版）》请在冰豆网上搜索。
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）<
alpha=alpha+pi;
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条线段两端点的坐标
functionw=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];
图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）;
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））;
X=XX;
Y=YY;
plot（X,Y）
%由以（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）;
a21,b21,a22,b22];
图2Levy曲线
（3）分形树程序tree.h
functiontree（n,a,b）
%tree（8,pi/8,pi/8）,n为分形树迭代次数
%a,b为分枝与竖直方向夹角
n=8;
a=pi/8;
b=pi/8;
x2=0;
y2=1;
plot（[x1,x2],[y1,y2]）
[X,Y]=tree1（x1,y1,x2,y2,a,b）;
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:
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）;
w1=W（:
w2=W（:
w=[w1;
%由每个分枝两端点坐标（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）<
a=pi+atan（（y2-y1）/（x2-x1））;
else
a=atan（（y2-y1）/（x2-x1））;
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）
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<
r（i）<
1/3时,进行M1对应的压缩映射;
%当1/3=<
2/3时,进行M2对应的压缩映射;
%当2/3=<
1时,进行M3对应的压缩映射;
n
ifr（i）<
w1
a=M1
（1）;
b=M1
（2）;
e=M1（3）;
c=M1（4）;
d=M1（5）;
f=M1（6）;
elseifr（i）<
w1+w2
a=M2
（1）;
b=M2
（2）;
e=M2（3）;
c=M2（4）;
d=M2（5）;
f=M2（6）;
w1+w2+w3
a=M3
（1）;
b=M3
（2）;
e=M3（3）;
c=M3（4）;
d=M3（5）;
f=M3（6）;
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）;
cx=-0.77;
cy=0.08;
%z^2+c=z0,x=real（z0）;
y=image（z0）;
x=1;
y=1;
%A为产生的服从标准正态分布的n维随机数组
A=randn（1,n）;
wx=x-cx;
wy=y-cy;
ifwx>
alpha=atan（wy/wx）;
ifwx<
alpha=pi+atan（wy/wx）;
ifwx==0
alpha=pi/2;
alpha=alpha/2;
r=sqrt（wx^2+wy^2）;
ifA（i）<
r=-sqrt（r）;
r=sqrt（r）;
x=r*cos（alpha）;
y=r*sin（alpha）;
B（1,k）=x;
B（2,k）=y;
k=k+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;
m
x0=a+（c-a）*（i-1）/m;
y0=b+（d-b）*（j-1）/m;
x=x0;
y=y0;
ify>
0.5
x=2*x;
y=2*y-1;
elseifx>
=0.5
x=2*x-1;
y=2*y;
ifx^2+y^2>
r
break;
ifk==n
B（1,w）=i;
B（2,w）=j;
w=w+1;
图6Sierpinski三角形垫片
（7）元胞自动机算法画Sierpinski三角形程序
✧一维元胞自动机sierpinski_ca1.h
functionsierpinski_ca1（m,n）
%sierpinski_ca1（1000,3000）
m=1000;
n=3000;
t=1;
w=zeros（2,m*n）;
s=zeros（m,n）;
s（1,fix（n/3））=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+1）==1）
s（i+1,j）=1;
w（1,t）=x+3+3*j;
w（2,t）=y+5*i;
t=t+1;
plot（w（1,:
）,w（2,:
1）
图7.1一维元胞自动机画Sierpinski三角形
✧二维元胞自动机sierpinski_ca2.h
functionsierpinski_ca2（m,n）
%sierpinski_ca2（400,400）
m=400;
n=400;
s（m/2,n/2）=1;
fori=[m/2:
-1:
2,m/2:
m-1]
forj=[n/2:
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;
图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];
图8Helix曲线
升级会员 