元胞自动机CA代码及应用.doc

元胞自动机 （CA） 代码及应用引言元胞自动机（ CA ）是一种用来仿真局部规则和局部联系的方法典型的元 胞自动机是定义在网格上的，每一个点上的网格代表一个元胞与一种有限的状 态变化规则适用于每一个元胞并且同时进行 典型的变化规则， 决定于元胞的 状态，以及其（ 4 或 8 ）邻居的状态元胞自动机已被应用于物理模拟，生物 模拟等领域本文就一些有趣的规则，考虑如何编写有效的 MATLAB 的程序来 实现这些元胞自动机MATLAB 的编程考虑 元胞自动机需要考虑到下列因素， 下面分别说明如何用 MATLAB 实现这些部分 并以 Conway 的生命游戏机的程序为例，说明怎样实现一个元胞自动机矩阵和图像可以相互转化，所以矩阵的显示是可以真接实现的如果矩阵cells的所有元素只包含两种状态且矩阵 Z含有零，那么用image函数来显示 cat命令建的RGB图像，并且能够返回句柄imh = image(cat(3,cells,z,z)); set(imh, 'erasemode', 'none') axis equal axis tight矩阵和图像可以相互转化，所以初始条件可以是矩阵，也可以是图形。

以下 代码生成一个零矩阵，初始化元胞状态为零，然后使得中心十字形的元胞状 态= 1z = zeros(n,n); cells = z; cells(n/2,.25*n:.75*n) = 1; cells(.25*n:.75*n,n/2) = 1;Matlab 的代码应尽量简洁以减小运算量以下程序计算了最近邻居总和，并 按照 CA 规则进行了计算本段 Matlab 代码非常灵活的表示了相邻邻居x = 2:n-1;y = 2:n-1;sum(x,y) = cells(x,y-1) + cells(x,y+1) + ... cells(x-1, y) + cells(x+1,y) + ... cells(x-1,y-1) + cells(x-1,y+1) + ... cells(x+1,y-1) + cells(x+1,y+1);cells = (sum==3) | (sum==2 & cells);加入一个简单的图形用户界面是很容易的在下面这个例子中，应用了三个 按钮和一个文本框三个按钮，作用分别是运行，停止，程序退出按钮文 框是用来显示的仿真运算的次数build the GUI %define the plot button plotbutton=uicontrol('style','pushbutton',...'string','Run', ...'fontsize',12, ... 'position',[100,400,50,20], ...'callback', 'run=1;');%define the stop button erasebutton=uicontrol('style','pushbutton',...'string','Stop', ...'fontsize',12, ... 'position',[200,400,50,20], ...'callback','freeze=1;');%define the Quit button quitbutton=uicontrol('style','pushbutton',...'string','Quit', ...'fontsize',12, ... 'position',[300,400,50,20], ...'callback','stop=1;close;');number = uicontrol('style','text', ...'string','1', ... 'fontsize',12, ...'position',[20,400,50,20]);经过对控件(和CA)初始化，程序进入一个循环，该循环测试由回调函数的每 个按钮控制的变量。

刚开始运行时，只在嵌套的 while 循环和 if 语句中运行 直到退出按钮按下时，循环停止另外两个按钮按下时执行相应的 if 语句stop= 0; %wait for a quit button pushrun = 0; %wait for a drawfreeze = 0; %wait for a freezewhile (stop==0)if (run==1)%nearest neighbor sum sum(x,y) = cells(x,y-1) + cells(x,y+1) + ...cells(x-1, y) + cells(x+1,y) + ... cells(x-1,y-1) + cells(x-1,y+1) + ... cells(3:n,y-1) + cells(x+1,y+1);% The CA rulecells = (sum==3) | (sum==2 & cells);%draw the new image set(imh, 'cdata', cat(3,cells,z,z) ) %update the step number diaplay stepnumber = 1 + str2num(get(number,'string')); set(number,'string',num2str(stepnumber)) end if (freeze==1)run = 0;freeze = 0;enddrawnow %need this in the loop for controls to work end例子1 .Conway 的生命游戏机。

规则是：对周围的 8 个近邻的元胞状态求和如果总和为 2 的话，则下一时刻的状态不改变如果总和为 3 ，则下一时刻的状态为 1 否则状态 = 0 核心代码：x = 2:n-1;y = 2:n-1;%nearest neighbor sumsum(x,y) = cells(x,y-1) + cells(x,y+1) + ...cells(x-1, y) + cells(x+1,y) + ... cells(x-1,y-1) + cells(x-1,y+1) + ... cells(3:n,y-1) + cells(x+1,y+1);% The CA rulecells = (sum==3) | (sum==2 & cells);2 . 表面张力规则是：对周围的 8 近邻的元胞以及自身的状态求和 如果总和 < 4 或= 5 ，下一时刻的状态 = 0 否则状态= 1核心代码：x = 2:n-1;y = 2:n-1; %nearest neighbor sumsum(x,y) = cells(x,y-1) + cells(x,y+1) + ...cells(x-1, y) + cells(x+1,y) + ... cells(x-1,y-1) + cells(x-1,y+1) + ... cells(3:n,y-1) + cells(x+1,y+1)+...cells(x,y);% The CA rulecells = ~((sum< 4) | (sum==5));3 . 渗流集群规则： 对周围相邻的 8 邻居求和(元胞只有两种状态， 0 或 1 )。

元胞也有一个单 独的状态参量(所谓 '记录 ' )记录它们之前是否有非零状态的邻居 在 0 与 1 之间产生一个随机数 r 如果总和 > 0 (至少一个邻居)并且 r >阈值，或者元胞从未有过一个邻居， 则元胞 = 1 如果总和 > 0则设置"记录"的标志，记录这些元胞有一个非零的邻居 核心代码：sum(2:a-1,2:b-1) = cells(2:a-1,1:b-2) + cells(2:a-1,3:b) + ...cells(1:a-2, 2:b-1) + cells(3:a,2:b-1) + ... cells(1:a-2,1:b-2) + cells(1:a-2,3:b) + ... cells(3:a,1:b-2) + cells(3:a,3:b);pick = rand(a,b); cells = cells | ((sum>=1) & (pick>=threshold) & (visit==O)) visit = (sum>=1);变量a和b是图像的尺寸最初的图形是由图形操作决定的 以下程序设定坐标 系为一个固定的尺寸，在坐标系里写入文本，然后获得并返回坐标内的内容，并 用getframe函数把它们写入一个矩阵ax = axes('u nits','pixels','positio n',[1 1 500 400],'color','k'); text('u nits', 'pixels', 'positi on', [130,255,0],...'stri ng','MCM','color','w','fo ntn ame','helvetica','fo ntsize',100)text(' un its', 'pixels', 'positi on', [10,120,0],...'stri ng','Cellular Automata','color','w','fo ntn ame','helvetica','fo ntsize',50)in itial = getframe(gca);[a,b,c]=size(i nitial.cdata); z=zeros(a,b);cells = double(i nitial.cdata(:,:,1)==255); visit = z ;sum = z;经过几十个时间间隔（从MCM Cellular Automata这个图像开始），我们可以得 到以下的图像。

50100 150 200 250 300 350400 450500 5504 .激发介质(BZ reaction or heart ) 规则：元胞有10个不同的状态状态0是体眠1-5为活跃状态，、6-9为是极活跃 状态计算每一个处于活跃的状态的元胞近邻的 8个元胞如果和大于或等于3 （至少有三个活跃的邻居），则下一时刻该元胞=1 < 不需要其它输入，1至9种状态依次出现如果该时刻状态=1那么下一时刻 状态=2如果该时刻状态=2，然后下一时刻状态=3，对于其它的状态 依次类推，直到第9种状态如果状态=9，然后下一状态=0并且元胞回 到休息状态核心代码：x = [2: n-1]; y = [2： n-1];sum(x,y) = ((cells(x,y-1)> 0)&(cells(x,y-1)< t)) + ((cells(x,y+1)> 0)&(cells(x,y+1)< t)) + ... ((cells(x-1, y)> 0)&(cells(x-1, y)< t)) + ((cells(x+1,y)> 0)&(cells(x+1,y)< t)) + ... ((cells(x-1,y-1)> 0)&(cells(x-1,y-1)< t)) + ((cells(x-1,y+1)> 0)&(cells(x-1,y+1)< t)) + ... ((cells(x+1,y-1)> 0)&(cells(x+1,y-1)< t)) + ((cells(x+1,y+1)> 0)&(cells(x+1,y+1)< t));cells = ((cells==0) & (sum>=t1)) + ...2*(cells=1) + .。