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);%更新生长核索引

4522277楼: Originally posted by 张曙亮 at 2012-04-12 11:32:18:
??? Error using ==> imshow>ParseInputs
Unknown option string "fit"

Error in ==> imshow at 85
[cdata, cdatamapping, clim, map, xdata, ydata, filename, ...

Error in ==>  ...

1楼: Originally posted by dbb627 at 2011-10-27 09:20:06:
这个板块帖子太少，发一个增加点人气，原始问题见
http://muchong.com/bbs/viewthread.php?tid=3723739
算法如下：
随机生长的步骤如下：
1) 建立材料生长的初始网格，并设定所有网格材料皆为主相； ...

4楼: Originally posted by cxguo at 2011-10-28 12:41:57:
楼主高手啊，请教以下可以用Matalab画原子生长模拟图不？三位岛状生长的。