跳转到内容
Jason's Study

Feautrier 多面体调度 — 把循环并行化变成解几何方程

是什么

Feautrier 多面体调度是一种让编译器自己看一段嵌套循环，然后决定哪几次循环可以同时跑、哪几次必须等前面跑完的方法。日常类比：像帮一个工厂规划装配流水线——同一道工序的零件能并排做就并排做，下一道工序得等上一道完才能开始。编译器要在不破坏依赖的前提下，让”同时做”的部分越多越好。

更具体地说，它把循环里的每一次迭代当成几何空间里的一个点，把”这次迭代要等那次迭代算完”当成空间里的一条箭头，然后求一个线性函数把每个点映射到时间步——这就是它名字里”仿射”和”多面体”的来源。

为什么重要

这篇 1992 年的论文是自动并行化和循环优化这两件事的数学根。不理解它，下面的事都讲不清：

  • 为什么 LLVM 里有个叫 Polly 的优化器能把朴素的三重循环矩阵乘法自动 tiling、变快几十倍
  • 为什么 MLIR 现在专门搞了一套叫 affine 的方言，作为 IR 一等公民
  • 为什么 AI 编译器（TVM / XLA / IREE）做算子融合、tile size 选择时还在引用 30 多年前的论文
  • 为什么”循环交换、循环倾斜、循环分块”这些经典优化能被一个统一的数学框架包起来，不再是一堆零散启发式

简单说：你今天写 for i: for j: C[i][j] += A[i][k]*B[k][j]，编译器能自动让它在 8 核 CPU 上跑满，背后就是这套理论。

核心要点

整个想法可以拆成 三层抽象：

  1. 把循环画成几何图形：嵌套循环 for i in 1..N: for j in 1..i: ... 的所有 (i,j) 取值组合，构成一个三角形区域。这个区域叫迭代域（iteration domain），它是一个多面体——边界都是直线（所以叫”多面体编译”）。

  2. 把依赖画成箭头：如果第 (2,1) 次迭代要用到第 (1,1) 次迭代写入的值，就在两点之间画一条箭头。所有这种箭头的集合叫依赖关系。约束：箭头方向必须从”早”指向”晚”。

  3. 找一个把点映射到时间的线性函数 θ：θ(i,j) = c·i + d·j + e，要满足”每条箭头的源点 θ 值 < 汇点 θ 值”，并且让所有点的 θ 最大值尽量小（θ 越小代表用的时间步越少，并行度越高）。求这个 θ 是一个线性规划问题，有限维度，能在多项式时间求解。

关键数学技巧叫 Farkas 引理——把”对多面体里无穷多个整数点都成立”这种看似无限的约束，等价转换成”对有限个系数成立”的线性不等式，于是可以丢给 ILP 求解器。

一句话总结：循环里有什么值要算 → 几何点；谁等谁 → 箭头；怎么排时间 → 线性函数；怎么找最优 → 解 ILP。整套思路的精妙在于”工程问题完整翻译成数学问题”。

一个对应表：

循环代码概念几何对象数学工具
for 嵌套多面体（迭代域）线性不等式系统
数据依赖多面体之间的有向连接仿射映射
调度（先后顺序）时间轴上的投影线性函数 θ
”找最快方案”投影最紧的那个线性规划求解

实践案例

案例 1：矩阵乘法的循环交换

for (i = 0; i < N; i++)
  for (j = 0; j < N; j++)
    for (k = 0; k < N; k++)
      C[i][j] += A[i][k] * B[k][j];

读不读得出来：i 和 j 两层完全独立（每个 C[i][j] 互不打架），只有 k 必须串行累加。Feautrier 算法读完依赖：

  • 把迭代点 (i,j,k) 看成 3D 空间里的立方体
  • 依赖箭头：(i,j,k) → (i,j,k+1)（沿 k 轴累加）
  • 求出 θ(i,j,k) = k 是合法调度——意味着”按 k 分层执行，每一层内 (i,j) 全部并行”

这就自动等价于”让 k 跑外循环、(i,j) 双层 omp parallel”。手写时这一步要靠经验，多面体编译能从依赖图机械推出来。换句话说，“凭经验”被替换成”解一个 LP”。

更进一步：编译器还能根据缓存大小自动决定 tile size——比如把 N×N 矩阵切成 32×32 的小块来跑，让每个块都能放进 L1 缓存。这一切都是从同一套调度推导框架里自然延伸出来的。

案例 2：循环倾斜（loop skewing）解锁波前并行

for (i = 1; i < N; i++)
  for (j = 1; j < N; j++)
    A[i][j] = A[i-1][j] + A[i][j-1];

依赖箭头：(i-1,j) → (i,j)、(i,j-1) → (i,j)。两层都不能直接并行，因为下一行依赖上一行。

但 Feautrier 的算法会自动解出 θ(i,j) = i + j——对角线方向才是真正的时间轴。同一条对角线上所有点（i+j = const）互不依赖，可以并排算。这就是大名鼎鼎的”波前并行”（wavefront）。手写时要”先想到沿对角线切”，多面体编译能从依赖图直接推出 i+j。这种”找到隐藏的并行轴”在手工优化里需要灵感，多面体框架直接当成线性规划的解。

案例 3：今天写代码时它在哪里

# JAX / XLA 看到这种代码
def matmul(A, B):
    return jnp.einsum("ij,jk->ik", A, B)

XLA HLO 里有一组 pass 在做 tiling、fusion、layout 选择，部分基于多面体表示。再下一层，MLIR affine dialect 直接把”循环 + 仿射访问”做成 IR 类型，可以调 ISL 求解依赖。这都是 Feautrier 1992 的工程后代。

案例 4：“调度”和”代码”是两个东西

在 Halide / TVM 这类 AI 编译器里，“算法”（要算什么）和”调度”（怎么算）是分开写的。这个把”算法 vs 调度”分离的设计哲学，思想根源就是 Feautrier——仿射调度 θ 本身就是一个独立于循环代码的数学对象，可以单独求解、单独优化。

踩过的坑

  1. “仿射”是硬约束：循环边界、数组下标都必须是循环变量的线性函数。A[i*i]、while (cond)、间接访问 A[idx[i]] 都直接不适用。这是这套理论用不到通用代码、只能用在科学计算/张量计算的根本原因。

  2. ILP 复杂度爆炸：循环嵌套深 / 数组多时，求解器可能跑几分钟。工业实现（PLuTo、ISL）做了大量启发式和近似来控制时间。

  3. 找到调度 ≠ 写出代码：θ 是一个数学函数，但要把它变回一段真正能跑的 for 循环（这个过程叫 code generation），需要 CLooG / ISL codegen 这样的专门工具，本身也很难。

  4. GPU / SIMD 还要再加一层：1992 年 Feautrier 想的是 SMP 多核。今天要把循环映射到 GPU block / thread / shared memory，得在多面体之上再叠一层 tiling 和 mapping 决策（PLuTo 的贡献）。

  5. 多维调度 vs 一维调度：1992 论文 Part I 处理一维调度（一个 θ），Part II 才处理多维调度。多维情况指数级复杂，工业实现几乎都用启发式裁剪而非真正的最优解。

适用 vs 不适用场景

适用：

  • 嵌套循环 + 静态边界 + 仿射数组下标的科学计算（线性代数 / 卷积 / 模板计算 / 偏微分方程求解）
  • AI 编译器优化张量算子（matmul、conv、reduce）
  • 自动并行化、自动 tiling、循环融合 / 分裂 / 倾斜
  • 编译器自动向量化（找出可以打包成 SIMD 指令的迭代组）

不适用：

  • 不规则数据结构（链表、图、稀疏矩阵）→ 用其他依赖分析（如 steensgaard-pointer）
  • 控制流复杂、循环边界依赖运行时变量 → 退回 kildall-dataflow 这种格论框架
  • 间接索引、while 循环、递归 → 需要不同的并行化思路

历史小故事（可跳过）

  • 1967 年：Lamport 写出”超平面法”，第一次用线性代数找循环并行轴。当时的目标是让 ILLIAC IV 这种向量机吃满。
  • 1988 年：Feautrier 自己先做了 Parametric Integer Programming 工具——求带参数的整数线性规划，为后面铺路。
  • 1991 年：他发表 Dataflow Analysis of Array and Scalar References，把数组的精确依赖分析建好。
  • 1992 年：本论文出场，分两部分发在 IJPP（International Journal of Parallel Programming），把”找最佳仿射调度”完整解决。
  • 2000s：Bondhugula 的 PLuTo 和 Verdoolaege 的 ISL 把它工业化。
  • 2010s：LLVM Polly、MLIR affine dialect 把这套搬进主流编译器生态。
  • 2020s：AI 编译器（XLA、TVM、IREE、Halide）继续在多面体之上加 schedule auto-tuning，把”找最优 θ”变成机器学习调参。

40 年下来，这套理论从”几个法国数学家的纸上推导”变成了 AI 编译器栈的底层依赖。Feautrier 本人在巴黎的 ENS Mines 做了一辈子，他自己的工具 PIP 至今仍是这条线的源头。

学到什么

  1. 把工程问题翻译成几何 + 线性规划——是计算机科学最强的招式之一，编译优化、机器学习、调度问题反复用
  2. Farkas 引理这种”无限到有限”的等价转换，能把看起来不可解的问题压成 ILP
  3. 限制换效率：仿射限制看起来吓人，但正因为限制，问题才变得可解、才能机械化
  4. 理论 → 算法 → 工程，1967 → 1992 → 2010s，每一代隔 10-20 年。今天 AI 编译器仍在吃这套老底
  5. 抽象层次跃迁：从一段 C 循环 → 几何对象 → 线性约束 → ILP，这种”层层换皮”是把复杂问题搬到合适数学工具上的通用招式

延伸阅读

关联

  • kildall-dataflow —— 数据流分析框架，是循环优化的”另一种几何”，基于格论而非多面体
  • ssa —— Polly 等多面体后端通常在 SSA 形式上工作
  • llvm —— Polly 是 LLVM 的多面体优化插件
  • lamport-time —— Lamport 把时间作为偏序，思路上和 Feautrier 把迭代映射到时间步同源
  • hughes-fp-matters —— 函数式”算法和调度分离”的思想在多面体编译里再次出现：θ 是独立于循环的数学对象

反向链接

  • cousot-halbwachs-polyhedra-1978 —— Cousot-Halbwachs 凸多面体域 — 让分析器自己发现变量间的线性关系
  • fpga-hls-2011 —— FPGA HLS 2011 — 把 C 代码自动翻译成芯片电路的范式
  • halide —— Halide — 把”算什么”和”怎么算”分开写
  • hughes-fp-matters —— Why FP Matters — 函数式真正赢在能拆能粘
  • kildall-dataflow —— Kildall 数据流框架 — 用一套格论统一所有全局编译优化
  • llvm —— LLVM — 模块化编译器框架
  • mlir —— MLIR — 给编译器一套乐高，每层抽象都能搭自己的方言
  • ssa —— SSA — 静态单赋值形式
  • steensgaard-pointer —— Steensgaard 指针分析 — 用等价合并把指针分析压到几乎线性
  • tvm —— TVM — 让一份模型能在所有硬件上跑得快
  • tvm-2018 —— TVM OSDI 2018 — 把 Halide 思想搬到深度学习