2022年完整word版,Matlab下三维DLA模型模拟 .pdf

《2022年完整word版,Matlab下三维DLA模型模拟 .pdf》由会员分享，可在线阅读，更多相关《2022年完整word版,Matlab下三维DLA模型模拟 .pdf（6页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、Matlab下三维 DLA 模型模拟2007-01-11 19:18分类： science字号：大中小function dla3dv5(Nsum,Wstep) %定义 dla 函数， Nsum 为所生成絮体包含的颗粒数， Wstep 为计算过程中所采取的步长% % 本程序内变量的定义% %radius为颗粒半径， release 为起始释放半径% % radius=0.5; % 颗粒半径release=2; % 起始释放半径L=200; Xhalf=floor(L/2); Yhalf=Xhalf; Zhalf=Xhalf; n=1; %粒子计数N(1)=1; N(Xhalf)=0; p(1,:

2、)=Xhalf Yhalf Zhalf; szpoints=zeros(L,L,L); % 网格点阵szpoints(Xhalf,Yhalf,Zhalf)=1; % 种子位置标志% 释放初始粒子theta=2*pi*rand; gama=pi*rand; 名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 1 页，共 6 页 - - - - - - - - - M=p(1,:)+Wstep*cos(theta) sin(theta) cos(gama); while n(release+

3、15)2 % 判断是否逃逸theta=2*pi*rand; gama=pi*rand; M=p(1,:)+release*cos(theta) sin(theta) cos(gama); elseif szpoints(T(1)-1),T(2),T(3)+szpoints(T(1)+1),T(2),T(3)+szpoints(T(1),(T(2)-1),T(3)+szpoints(T(1),(T(2)+1),T(3)+szpoints(T(1),T(2),(T(3)-1)+szpoints(T(1),T(2),(T(3)+1)0&szpoints(T(1),T(2),T(3)=1 % 判断是否

4、凝结n=n+1; szpoints(T(1),T(2),T(3)=1; p(n,:)=T; %存储凝聚颗粒的球心坐标。s=sqrt(M(1)-Xhalf)2+(M(2)-Yhalf)2+(M(3)-Zhalf)2); k=round(s)+1; N(k)=N(k)+1; if srelease % 调整释放半径release=s+3; end elseif szpoints(T(1),T(2),T(3)=1 % 检查是否出现漏检,即运动一步后进入粒子内部的情况. 名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 -

5、- - - - - - 第 2 页，共 6 页 - - - - - - - - - theta=2*pi*rand(1); gama=pi*rand(1); M=p(1,:)+Wstep*cos(theta) sin(theta) cos(gama); end end nmax=size(p,1); j1=0; j2=0; j3=0; for i=1:nmax if p(i,1)=Xhalf j1=j1+1; X(j1,:)=p(i,:); end if p(i,2)=Yhalf j2=j2+1; Y(j2,:)=p(i,:); end if p(i,3)=Zhalf j3=j3+1; Z(j

6、3,:)=p(i,:); end end %绘制立体图figure(1); for i=1:nmax ssphere(p(i,:),radius); hold on 名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 3 页，共 6 页 - - - - - - - - - end shading interp; colormap(gray); title( 絮凝分形仿真模拟结果); %沿轴线切割图形绘制，分别为垂直于x,y,z 轴的切割面图figure(2); for i=1:j1 ss

7、phere(X(i,:),radius); hold on end shading interp; colormap(gray); title( 过中心垂直于X 轴的切割面 ); %y=X(:,2); %z=X(:,3); %plot(y,z,+); figure(3); for i=1:j2 ssphere(Y(i,:),radius); hold on end shading interp; colormap(gray); title( 过中心垂直于Y 轴的切割面 ); %x=Y(:,1); %z=Y(:,3); %plot(x,z,+); figure(4); 名师资料总结 - - -精

8、品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 4 页，共 6 页 - - - - - - - - - for i=1:j3 ssphere(Z(i,:),radius); hold on end title( 过中心垂直于Z 轴的切割面 ); %x=Z(:,1); %y=Z(:,2); %plot(x,y,+); shading interp; colormap(gray); % 为了与 dlacon 连用而进行存储数据，不用时可以注销%save datap p; %save datarelease rel

9、ease; %save dataszp szpoints; %save dataM M; %save dataNsum Nsum; %save datan n; %save datahalf Xhalf; % % % 分维分析计算% num=find(N); R(1)=0.5; for i=2:size(num,2) 名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 5 页，共 6 页 - - - - - - - - - N(i)=N(i)+N(i-1); num(i)=N(i); R

10、(i)=(i-0.5)/0.5; end num(1)=; R(1)=; %save dataN num; %save dataR R; figure(5); plot(log(R),log(num),*); % %figure(3); %kxl=1-num.*(radius./R).3; %plot(R,kxl,*); % % 图形绘制% %figure; %绘制回转半径与其内粒子数的对数关系图% %plot(R, N,*); % %figure(3); %绘制凝聚粒子数目与步长的关系图% %plot(Nn,I,o); % % 名师资料总结 - - -精品资料欢迎下载 - - - - - - - - - - - - - - - - - - 名师精心整理 - - - - - - - 第 6 页，共 6 页 - - - - - - - - -