简单二维元胞自动机 MATLAB实现

20世纪50年代，乌尔姆和冯·诺依曼（对此人真是崇拜的五体投地）为了研究机器人自我复制的可能性，提出了一种叫做元胞自动机（Cellular Automaton,CA）的算法。该算法采用局相互作用规则，最终产生整体的自复制构型（和蒙特卡罗法简直就是完全两种不同的哲学思想，竟然和同一个人有关）。元胞自动机已经成为研究复杂系统行为的一个理论框架（蒙特卡洛法也是！！），也是人工智能的雏形（alpha GO 也用到了一种蒙特卡洛树搜索的技术！！！）。

原本打算每天晚上抽出点时间研究一个分形图形的，可惜时间管理太差，没能坚持，但学习分形学确实让人收获颇丰！而且一维元胞自动机大一时C语言简单实现过。所以，实现二维的也没有太大难度。但就实现速度和代码质量而言，已经深刻认识到了自己内力（算法）太浅薄，写的就是白菜代码。

规则

设有axa的网格（a为奇数），中间点设为1，其他点设为0；则此网格共有q/2层，设k为层数k=1~a/2;进行如下循环：

1.设k层中任一点（Xi, Yi）,则判断其八个邻居的现有状态；其和为奇数则（Xi，Yi）为1，否则为0；

2.对第k层中（Xi, Yi）= 1的点着色；

3.k = k + 1;返回步骤1.

代码

 1 % 二维元胞自动机
 2 % 2016-4-8
 3 % zhaoyu
 4 % 规则，先把中间点置为1，每一时间步对每一点，如果周围
 5 % 八个点和为偶数，则变为0，为奇数则变为 1
 6
 7 % 颜色控制
 8 Map = [1 1 1; 0 0 0];
 9 colormap(Map);
10 % 设置网格大小
11 S = 121;
12 L = zeros(S);
13 % 把中间一个数设置为 1 作为元胞种子
14 M = (S+1)/2;
15 L(M, M) = 1;
16 Temp = L;
17 imagesc(L);
18
19 % 计算层数
20 Layer = (S-1)/2 + 1;
21
22 for t=2:Layer
23     for x=M-t+1:M+t-1
24        if x==M-t+1 || x==M+t-1
25
26           for y=M-t+1:M+t-1
27             SUM = 0;
28             for m=-1:1
29                for n=-1:1
30                   if x+m>0 && x+m<=S && y+n>0 && y+n<=S
31                      SUM = SUM + L(x+m, y+n);
32                   end
33                end
34             end
35             SUM = SUM - L(x, y);
36             Temp(x, y) = mod(SUM, 2);
37           end
38
39        else
40             y = M-t+1;
41             SUM = 0;
42             for m=-1:1
43                for n=-1:1
44                   if x+m>0 && x+m<=S && y+n>0 && y+n<=S
45                      SUM = SUM + L(x+m, y+n);
46                   end
47                end
48             end
49             SUM = SUM - L(x, y);
50             Temp(x, y) = mod(SUM, 2);
51
52             y = M+t-1;
53             SUM = 0;
54             for m=-1:1
55                for n=-1:1
56                   if x+m>0 && x+m<=S && y+n>0 && y+n<=S
57                      SUM = SUM + L(x+m, y+n);
58                   end
59                end
60             end
61             SUM = SUM - L(x, y);
62             Temp(x, y) = mod(SUM, 2);
63        end
64     end
65     L = Temp;
66     imagesc(L);
67     % 速度控制
68     pause(0.2);
69 end

生成图形