元胞自动机在交通流建模中的应用（含 Python 完整代码）
数学建模国赛有一个经久不衰的题型——交通流。
2023 年国赛 B 题是高铁站客流，2024 年 A 题是城市路网潮汐车道，2025 年 C 题又考了共享单车调度。每年都有人问：交通流问题到底用什么模型？
答案之一是：元胞自动机（Cellular Automaton, CA）。
它不需要 Navier-Stokes 那种复杂的偏微分方程，几行代码就能模拟出真实的堵车、走走停停、幽灵堵车（ghost traffic jam）现象。这篇文章带你从原理到代码，30 分钟掌握这套工具。
一、什么是元胞自动机？
把一条公路想象成一排格子，每个格子代表一辆车或一段空地。元胞自动机用简单的规则决定每一秒每辆车的行为：
- 加速
- 减速
- 随机慢化
- 位置更新
就这么四步。上百辆车同时遵守这套规则，就会涌现出复杂的宏观现象——堵车、消散、绿波带。这种"局部简单 → 全局复杂"的特性，正是元胞自动机的魅力。
二、Nagel-Schreckenberg 模型（NaSch）
1992 年两位德国物理学家提出的 NaSch 模型是交通流 CA 的标准。下面是完整实现：
import numpy as np import matplotlib.pyplot as plt from matplotlib.animation import FuncAnimation # 参数 ROAD_LENGTH = 200 # 道路长度（格子数） NUM_CARS = 40 # 车辆数 MAX_SPEED = 5 # 最大速度 P_SLOW = 0.3 # 随机慢化概率 STEPS = 200 # 模拟步数 class NaSchSimulation: def __init__(self): # 初始化：随机放置车辆，速度随机为 0~MAX_SPEED self.positions = np.zeros(ROAD_LENGTH, dtype=int) car_idx = np.random.choice(ROAD_LENGTH, NUM_CARS, replace=False) self.positions[car_idx] = 1 self.speeds = np.zeros(ROAD_LENGTH, dtype=int) self.speeds[car_idx] = np.random.randint(0, MAX_SPEED+1, NUM_CARS) def step(self): new_positions = np.zeros(ROAD_LENGTH, dtype=int) new_speeds = np.zeros(ROAD_LENGTH, dtype=int) for i in range(ROAD_LENGTH): if self.positions[i] == 0: continue # 1. 计算与前车的距离 gap gap = 1 while self.positions[(i + gap) % ROAD_LENGTH] == 0: gap += 1 gap -= 1 # 实际空格数 # 2. 加速 v = min(self.speeds[i] + 1, MAX_SPEED) # 3. 减速 v = min(v, gap) # 4. 随机慢化 if np.random.rand() < P_SLOW: v = max(v - 1, 0) # 5. 位置更新 new_pos = (i + v) % ROAD_LENGTH new_positions[new_pos] = 1 new_speeds[new_pos] = v self.positions = new_positions self.speeds = new_speeds def get_density(self): return NUM_CARS / ROAD_LENGTH def get_mean_speed(self): return np.mean(self.speeds[self.positions == 1]) # 运行模拟 sim = NaSchSimulation() history = [] for t in range(STEPS): history.append(sim.positions.copy()) sim.step() # 可视化时空图（车流演化） plt.figure(figsize=(12, 6)) plt.imshow(np.array(history).T, cmap='Greys', aspect='auto', origin='lower') plt.xlabel('Time Step') plt.ylabel('Position') plt.title('NaSch Traffic Simulation (Spatiotemporal Diagram)') plt.colorbar(label='Vehicle Presence') plt.tight_layout() plt.savefig('nasch_spacetime.png', dpi=200) plt.show()
跑完你会看到一张经典的"时空图"——横轴是时间，纵轴是位置，黑点代表车。你会看到斜线代表车辆移动，斜线汇聚的地方就是堵车点。
三、复现"幽灵堵车"
什么是幽灵堵车？高速路上没有任何事故，突然某一段开始堵，几分钟后堵车消散，再过一会儿又堵起来——这是元胞自动机最经典的涌现现象。
复现方法很简单：把上面代码里的 P_SLOW = 0.3 调高到 0.5，让更多司机随机刹车。跑 500 步后你会看到堵车波自发形成。
下面这段代码可以画出"流量-密度"基本图，这是交通工程领域的标准分析：
densities = np.arange(0.05, 0.95, 0.05) flows = [] for rho in densities: cars = int(ROAD_LENGTH * rho) sim = NaSchSimulation.__new__(NaSchSimulation) sim.positions = np.zeros(ROAD_LENGTH, dtype=int) car_idx = np.random.choice(ROAD_LENGTH, cars, replace=False) sim.positions[car_idx] = 1 sim.speeds = np.zeros(ROAD_LENGTH, dtype=int) sim.speeds[car_idx] = np.random.randint(0, MAX_SPEED+1, cars) speeds = [] for _ in range(100): sim.step() speeds.append(sim.get_mean_speed()) flows.append(np.mean(speeds) * rho) plt.figure(figsize=(8, 5)) plt.plot(densities, flows, 'o-', linewidth=2) plt.xlabel('Density (cars/cell)') plt.ylabel('Flow (cars/step)') plt.title('Fundamental Diagram of Traffic Flow') plt.grid(True, alpha=0.3) plt.tight_layout() plt.savefig('fundamental_diagram.png', dpi=200)
你会看到一条倒 U 形曲线——密度低时流量随密度上升，密度高时流量反而下降。曲线最高点就是道路的临界密度，对应的流量是最大通行能力。
四、模型扩展：红绿灯与多车道
4.1 加入红绿灯
class NaSchWithTrafficLight(NaSchSimulation): def __init__(self, light_position=100, light_cycle=20): super().__init__() self.light_pos = light_position self.cycle = light_cycle self.light_state = 'red' # 初始红灯 def step(self): # 红绿灯状态切换 if self.step_count % self.cycle < self.cycle / 2: self.light_state = 'red' else: self.light_state = 'green' self.step_count += 1 # 改写减速规则：红灯时靠近路口要停 for i in range(ROAD_LENGTH): if self.positions[i] == 0: continue gap = self._compute_gap(i) v = min(self.speeds[i] + 1, MAX_SPEED) v = min(v, gap) # 红灯逻辑：在路口前停下 if self.light_state == 'red': dist_to_light = (self.light_pos - i) % ROAD_LENGTH if 0 < dist_to_light <= v + 1: v = dist_to_light - 1 if np.random.rand() < P_SLOW: v = max(v - 1, 0) new_pos = (i + v) % ROAD_LENGTH # ... 后续同父类
4.2 加入多车道
多车道 CA 还需要加入换道规则——车辆在什么条件下变道：
def lane_change_decision(car): # 换道动机：当前车道太慢 if car.speed < MAX_SPEED: # 换道条件：本车道有前车阻挡 + 目标车道有空间 if front_gap < safe_distance and side_gap > safe_distance: return True return False
加入多车道后，你会发现道路通行能力显著提升——这正是现实高速的写照。
五、参数敏感性分析（O 奖必做）
光跑通模型还不够，O 奖论文必须有 sensitivity analysis：
# 扫描随机慢化概率 P_SLOW p_range = [0.0, 0.1, 0.2, 0.3, 0.4, 0.5] mean_flows = [] for p in p_range: P_SLOW = p sim = NaSchSimulation() flows_at_p = [] for _ in range(200): sim.step() flows_at_p.append(sim.get_mean_speed() * NUM_CARS / ROAD_LENGTH) mean_flows.append(np.mean(flows_at_p[100:])) # 去掉瞬态 plt.plot(p_range, mean_flows, 's-', linewidth=2) plt.xlabel('Random Slowdown Probability P') plt.ylabel('Mean Flow') plt.title('Sensitivity: How Driver Distraction Affects Traffic') plt.grid(alpha=0.3) plt.tight_layout() plt.savefig('sensitivity_p_slow.png', dpi=200)
你会在论文里写："P_slow 从 0.0 升到 0.5，平均流量下降 67%，临界密度从 0.31 降到 0.18。这意味着减少 30% 的急刹行为就能让道路通行能力翻倍。"
这种"模型 → 政策建议"的链条，是 O 奖和 M 奖的分水岭。
六、CA 模型在国赛中的实战套路
如果你的赛题是交通流，下面这套模板可以直接套：
- 问题 1（数据预处理）：清洗轨迹数据 / 收费站数据 → 得到 density、flow、speed
- 问题 2（基础模型）：用 NaSch 模型拟合现状，画 fundamental diagram
- 问题 3（场景模拟）
- 问题 4（策略优化）
每个问题配一张时空图 + 一张流量密度图，论文图表就达标了。
七、局限性与改进方向
CA 模型的局限性也写在论文里能加分：
- 离散化误差：格子长度通常取 7.5m（一个车身+安全距离），但实际车长差异大
- 同质化假设
- 忽略横向动力学
进一步改进： - 换成 连续模型：最优速度模型（OVM）、智能驾驶模型（IDM） - 加入 强化学习：用 RL 训练自适应巡航 - 引入 多智能体：用 SUMO 软件做精细仿真
关注「数学建模小秘籍」，获取更多国赛/美赛备赛干货。 下一篇会讲《多目标优化算法实战对比：NSGA-III / MOEA/D / SPEA2 怎么选》，把数学建模里另一大类"既要又要"的问题彻底讲清楚。
10个月宝宝每天需要喝多少奶粉？
