matlab怎么实现行随机网络,matlab生成随机行走网格 模拟材料生长 - 计算模拟 - 小木虫 - 学术 科研 互动社区...

CODE:

function growth

%作者：dbb627

%2011.10.27

%======设定初值=============

clc;clear;

rand('state',0);

set(gcf,'DoubleBuffer','on');%消除振动

Cd=0.01;%生长核分布概率数

Vr=0.3;%生长相的体积分数

global kn

kn=500;

p=[1/8 1/8 1/8 1/8 1/4 1/4 1/4 1/4];%方向增长概率数

A=rand(kn,kn);%初始网格kn×kn

A(find(A<=Cd))=0;%黑色，生长核

A(find(A>Cd))=1;%白色

imshow(A,'InitialMagnification','fit')%初始图像

k=0;%迭代步骤

N=kn*kn;%格点数

[x,y]=find(A==0);%生长核脚标

Ixy=sub2ind([kn,kn],x,y); % 转化脚标为索引

gn=length(Ixy);%生长核计数

Vs=gn/N;       %体积分数

ti=title(['time = 0 体积分数=',num2str(Vs)],'Fontsize',14,'Fontname','Times new roman'); % 在图题处实时显示时间

%===========================

while Vs

xL=x-1;xR=x+1;

xL(xL==0)=1;% 对近邻进行边界条件处理

xR(xR==kn+1)=kn;% 对近邻进行边界条件处理

yD=y-1;yU=y+1;

yD(yD==0)=1;% 对近邻进行边界条件处理

yU(yU==kn+1)=kn;% 对近邻进行边界条件处理

InE=sub2ind([kn,kn],xR,y); % 右方转化脚标为索引

InS=sub2ind([kn,kn],x,yD); % 下方转化脚标为索引

InW=sub2ind([kn,kn],xL,y); % 左方转化脚标为索引

InN=sub2ind([kn,kn],x,yU); % 上方转化脚标为索引

InEN=sub2ind([kn,kn],xR,yU); % 右上方转化脚标为索引

InES=sub2ind([kn,kn],xR,yD); % 右下方转化脚标为索引

InWS=sub2ind([kn,kn],xL,yD); % 左下方转化脚标为索引

InWN=sub2ind([kn,kn],xL,yU); % 左上方转化脚标为索引

%===用函数的方法写八个方向的生长迭代=====

[A,Ixy]=growthmod(InE,Ixy,p(5),A);

[A,Ixy]=growthmod(InS,Ixy,p(6),A);

[A,Ixy]=growthmod(InW,Ixy,p(7),A);

[A,Ixy]=growthmod(InN,Ixy,p(8),A);

[A,Ixy]=growthmod(InEN,Ixy,p(1),A);

[A,Ixy]=growthmod(InES,Ixy,p(2),A);

[A,Ixy]=growthmod(InWS,Ixy,p(3),A);

[A,Ixy]=growthmod(InWN,Ixy,p(4),A);

%==========做图========

[x,y]=find(A==0);

Ixy=sub2ind([kn,kn],x,y);

gn=length(Ixy);

Vs=gn/N %体积分数

k=k+1

pause(0.1)

imshow(A,'InitialMagnification','fit')%更新图像

title(['time = ',num2str(k),'体积分数=',num2str(Vs)],'Fontsize',14,'Fontname','Times new roman'); % 更新时间参数

end

function [A,Ixy]=growthmod(Ind1,Indxy,p,A)

%Ind1 当前方向索引

%,Indxy 当前生长核索引的

%p 当前方向的概率值

%A 当前的生长矩阵

global kn

CInE=setdiff(Ind1,Indxy);%提取非生长核边界格点索引

En=rand(1,length(CInE));%各个单元均会被赋值一个随机数

[i,j]=find(En<=p);%大于该方向上的方向增长概率数，则此单元同样成为一个新的生长相

pE=sub2ind([1,length(CInE)],i,j);%新的生长相索引的下标

A(CInE(pE))=0; %更新生长矩阵

[x,y]=find(A==0);%更新生长核坐标

Ixy=sub2ind([kn,kn],x,y);%更新生长核索引