DLA模型的Matlab程序

大学里写了一些关于DLA模型的Matlab程序，由于不懂编程，而且也不懂matlab，在程序算法的编制上走了很多弯路，这里将陆续贴出关于DLA、DLCA、先DLA后DLCA模型的Matlab下实现和分维分析，请大家指导。程序按gpl发布（:(不明白gpl什么意思，只知道共享），如果你能有所改进也希望你能与我分享。二维dla模型程序：文件名：dla2d.m功能：生成所需粒子数的二维dla絮体，绘制回转半径与粒子数的对数关系图，空隙率与絮体半径的关系图。代码：%(1)设置初始凝结半径,总凝集粒子数nsum,绘图区域;%(2)根据绘图区域确定凝集中心在屏幕上的坐标，并在中心位置画点，置n=0;%(3)根据凝结半径确定自由离子释放半径r，，最大可凝集半径;%(4)在半径r0的圆周的任意位置释放一粒子，让其作布朗运动;%(5)粒子每移动一次，都要检测其上下左右的紧邻位置是否为凝集体，若是则凝集下来，n=n+1.执行第(6)步,否则，执行第(7)步;%(6)每凝结一个粒子，都要计算此粒子与凝集中心的距离是否大于凝结半径r0,%若是则将凝结半径置为当前距离，然后执行第4步;%(7)若粒子与凝集中心的距离小于rmax，继续作布朗运动。否则抛弃此粒子，执行第4步;%循环执行第3、4、5、6、7步，直到n=nsum为止.functiondla2d(nsum)r0=2;r=8;Rr=0.5;L=2000;Xhalf=floor(L/2);Yhalf=Xhalf;szPoints=sparse(L,L);szPoints(Xhalf,Xhalf)=1;nx=Xhalf;ny=Yhalf;n=0;theta=2*pi*rand;x=nx+r*cos(theta);y=ny+r*sin(theta);nx=round(x);ny=round(y);i=0;j=1;p=-1:1;N2(Xhalf/2)=0;N2(1)=1;R=2;%%用于计算生长维的一个变量whilennsumi=i+1;theta=2*pi*rand;x=nx+r0*cos(theta);y=ny+r0*sin(theta);nx=round(x);ny=round(y);ifsqrt((nx-Xhalf)^2+(ny-Yhalf)^2)r+20%theta=2*pi*rand;nx=Xhalf+r*cos(theta);ny=Yhalf+r*sin(theta);elseifszPoints(nx-1,ny)+szPoints(nx+1,ny)+szPoints(nx,ny+1)+szPoints(nx,ny-1)0&szPoints(nx,ny)==0szPoints(nx,ny)=1;n=n+1;s=sqrt((nx-Xhalf)^2+(ny-Yhalf)^2);k=round(s)+1;N2(k)=N2(k)+1;ifsRR=s;R1(j)=log(R/Rr);N1(j)=n;I(j)=i;j=j+1;endifsrr=s+5;endelseifszPoints(nx,ny)==1%theta=2*pi*rand;nx=Xhalf+r*cos(theta);ny=Yhalf+r*sin(theta);endend[u,v]=find(szPoints);figure(1);plot(u,v,'*');%%%%%%%%%%%%%%%%%%%%%%%%%%%分维采样数据分析整理%%%%%%%%%%%%%%%%%%%%%%%%%num=find(N2);b=size(num,2);fori=2:bN2(i)=N2(i)+N2(i-1);num(i)=N2(i);R2(i)=log(2*i);endR2(1)=[];num(1)=[];savenum2num;saveR2R2;figure(2);plot(R2,log(num),'*');%%%%%%%%%%%%%%%%%%%%%%%%%生长维曲线绘制%%%%%%%%%%%%%%%%%%%%%%%figure(3);%画回旋半径与其内粒子数的对数关系图plot(R1,log(N1),'*');%%%%%%%%%%%%%%%%%%%%%%%%%粒子数量与运行步数的关系图绘制%%%%%%%%%%%%%%%%%%%%%%%%figure(4);%画凝聚粒子数量与运行步数的关系图%%plot(N1,I,'o');