[Python练习向] 简易元胞自动机框架

　　关于元胞自动机（Cellular Automata, CA）的原理这里就不再叙述了。

　　以前的学校地理系有用CA做城市发展研究的大牛，无论在课堂上还是在项目中都会用得到CA模型。虽然CA模型本身并不复杂，但是每次从新写起也是十分麻烦，因此一个通用的CA模型框架能够减少很多工作量。幸好对于大多数CA模型的应用来说，其邻域定义大同小异，最主要的不同在于元胞的转换规则会根据问题和方法的不同而有所变化。因此这个通用框架的核心在于，用户可以自行定义元胞的转换规则。

　　PyCA构建CA模型的流程如下：

　　创建转换规则类需要两类参数：

　　1）邻域类型和大小

　　在PyCA中定义了三种类型的邻域： von、moore和custom，分别对应Von Neumann邻域（大小为1时即四邻域）、Moore领域（大小为1时即八邻域）和自定义邻域。

　　其中2代表中心元胞。如果需要自定义领域，则需要如上图所示自行定义邻域矩阵，使用1代表有效邻域，2代表中心元胞，其他位置则用0表示。需要注意的是，PyCA没有其他依赖numpy，定义矩阵时必须使用原生的list。

　　2）自定义转换规则

　　在定义转换规则时，用户不需要了解具体的模型运行流程参数。需要处理的参数有三个，包括当前中心元胞的值、有效邻域元胞的值和环境因子，并根据上述参数给出下一次循环时中心元胞的值。因此转换规则的定义如下：

newvlaue = rule(oldvalue, neighborhood, envi)

　　neighborhood和envi为邻域元胞值和环境参数。这个两个参数会传入和之前定义邻域大小相一致的矩阵，有效邻域中保留对应的元胞或环境参数，无效邻域将赋值为None。由于该框架支持多层环境参数，当输入的环境参数有多个图层时，neighborhood是一个包含每一层邻域的list。

　　以最简单的生命游戏为例，其元胞转换规则应该如下：

def rule_lifegame(pix, neighborhood, envi):
    # Replace all the None with 0
    for r in range(3):
        for c in range(3):
            if neighborhood[r][c] == None:
                neighborhood[r][c] = 0
    
    # 1 means living, and 0 means dead
    if pix == 1:
        # The central cell can not be counted.
        live = numpy.sum(neighborhood) - 1
        if live < 2:
            return 0
        elif live <= 3:
            return 1
        else:
            return 0
    else:
        live = numpy.sum(neighborhood)
        if agent == 3:
            return 255
        else:
            return 0;

　　当定义好转换规则时，就可以用rulebuilder类来创建模型规则。rulebuilder类有两种初始化方法，分别对应预定义邻域（von、moore）和自定义领域:

# Predefined neighborhood: 'von' and 'moore'
rulebuilder = PyCA.rulebuilder(rule, ntype, nsize)
# Custome neighborhood:
rulebuilder = PyCA.rulebuilder(rule, neighboorhoodregion)

　　当新建好rulebuilder类后就可以调用build()方法类创建模型规则，以上面的生命游戏为例：

rulebuilder = PyCA.rulebuilder(rule_lifegame, 'von', 1)
carule = rulebuilder.build()

　　最后一步是利用CA类来创建CA模型。新建CA类需要三个参数：由rulebuilder.build()生成的模型规则、初始元胞分布、环境参数：

ca = PyCA.CA(envi, initialcells, carule)

　　然后调用CA类中的evolve()方法来进行迭代：

ca.evolve(time, returnresult, storehistory)

　　其中time指定迭代次数，默认值为1。returnresult指定在完成所有迭代后是否返回当前元胞分布，默认为False。storehistory指定每次迭代后是否将元胞分布保存到ca.history中，默认值为False。ca.history是一个dict，定义为 ca.history = {time: result of that time}。用户可以在完成迭代后调用ca.history得到每次迭代结果。

　　至此，一个完整的CA模型创建和运行流程完成。PyCA源码地址：http://pan.baidu.com/s/1mgwjkLE

　　最后，给出一个完整生命游戏Demo:

import PyCA, numpy, random
from PIL import Image

# rule for life game
def rule_lifegame(pix, neighborhood, envi):
    for r in range(3):
        for c in range(3):
            if neighborhood[r][c] == None:
                neighborhood[r][c] = 0
    
    if pix == 255:
        live = numpy.sum(neighborhood) - 255
        if live < 2 * 255:
            return 0
        elif live <= 3 * 255:
            return 255
        else:
            return 0
    else:
        live = numpy.sum(neighborhoodn)
        if agent == 3 * 255:
            return 255
        else:
            return 0;
            
# build initial cell distributioin, randomly pick up 3000 cells as living cells
(row, col) = (100, 100)
cells= numpy.zeros((row, col)).tolist()
for i in range(3000):
    cells[random.randint(0, 99)][random.randint(0, 99)] = 255

# build trasformation rules
rulebuilder = PyCA.rulebuilder(rule_lifegame, 'von', 1)
carule = rulebuilder.build()

# create CA model
ca = PyCA.CA(cells, cells, carule)
ca.evolve(50, False, True)    

# save the array of result as bmp file
for (time, result) in ca.history.items():
    arr = numpy.array(result, dtype = numpy.uint8)
    im = Image.fromarray(arr)
    im.save('time_{0}.bmp'.format(time))

P.S: 为避免严重的边界效应，该模型会忽略处于边界的元胞。

如有什么问题欢迎私信讨论。

[Python练习向] 简易元胞自动机框架

推荐阅读