模拟疫情感染情况的并行化优化的元胞自动机及其代码

来源：https://maiya.fan/blog?id=6566173a6aa58e39d930190b

引言

元胞自动机(cellular automata，CA) 是一种时间、空间、状态都离散，空间相互作用和时间因果关系为局部的网格动力学模型，具有模拟复杂系统时空演化过程的能力。这里用元胞自动机来模拟疫情病毒感染情况。
相比于其他元胞自动机的实现，该项目使用pytorch实现并行化处理元胞自动机的状态更新。并且适配了CPU和GPU，极大提高了运行速度。修改后可以适用于其他元胞自动机的情况。

基于元胞自动机的病毒感染模型

在元胞自动机内，每个点代表一个人。

状态分类

每个人具有不同的状态，具体类别如下：

死亡人群0
痊愈人群1
易感人群2
潜伏人群3
无症状感染者4
轻症感染者5
重症感染者6
危重症感染者7
隔离者8

状态变化关系

状态变化基于某点周围的点的感染情况进行改变。这里以易感人群为例，具体的变化代码如下，完整代码见文章结尾：

        x_clone[((((1 - p_isolate2) * p_infection2 + p_isolate2) > rand) * (rand >= p_isolate2)) * (x == 2)] = 3
        if day >= self.isolate_start:
            need_isolate = (rand < p_isolate2) * (x == 2)
            num_need_isolate = torch.sum(need_isolate)
            if res_num_isolate >= num_need_isolate:
                x_clone[need_isolate] = 8
                # res_num_isolate = res_num_isolate - num_need_isolate
            else:
                need_isolate_index = torch.where(need_isolate)
                num_cannot_isolate = num_need_isolate - res_num_isolate
                rand_isolate = torch.tensor(random.sample(range(num_need_isolate), num_cannot_isolate)).to(self.device)
                need_isolate[0, need_isolate_index[1][rand_isolate], need_isolate_index[2][rand_isolate]] = False
                x_clone[need_isolate] = 8

项目输入

该项目可以调节元胞自动机的各种参数。其中，节点数量分为横轴和竖轴，两者相乘即为总的人数。num_person3为初始感染人数。p_infection2、p_infection3、…、p_infection8为对应类型的人的状态转移的概率。如p_infection6代表重症感染者的转移概率，其由两个参数构成一个列表[a,b]，a代表重症患者转化为危重症患者的概率，b代表重症患者痊愈的概率。
输入部分的代码如下，完整代码见文末：

    hight = 20 * 10
    width = 30 * 10
    num_person3 = 100
    num_classes = 9
    color_mul = [53, 19, 41]

    env = init_env(hight, width, num_person3=num_person3)
    figure, plot1, trans, ax1, data_list, lines_list3, ax3, figure2, ax4, lines_list4 = init_plot(env, num_classes,
                                                                                                  color_mul)
    op = Operate(num_classes=num_classes, p_infection2=[0.0, 0.3, 0.3, 0.2, 0.2, 0.9], isolate_start=152,
                 p_infection3=[0.1, 0.1],
                 p_infection4=[0.12, 0.8], p_infection5=[0.075, 0.1667, 0.750717], p_infection6=[0.0417, 0.0917],
                 p_infection7=[0.05, 0.075], p_infection8=[0.3], kernel_size=(3, 3), neighbor_mode=0,
                 num_isolate=2579, device=device).to(device)

结果展示

ticks代表该状态的轮数。

总人数为60000，初始感染数为100

总人数为6000000，初始感染数为100

完整代码：