第 27 章 多线程算法
本章对应原书第七部分导言与第 27 章。前面的算法都假设一次只执行一条指令。今天每台笔记本都是多核共享内存机器,量化研究里的参数扫描、蒙特卡洛模拟、因子批量计算又天然可以并行——问题是:并行之后能快多少?加机器还有没有用?原书给出了一套干净的回答:把计算看成一张有向无环图,用工作量(总的计算量)和跨度(关键路径长度)两个数刻画它,再用两条定律和一个调度定理估计任意处理器数下的运行时间。本章把这套工具讲清楚,然后用它分析并行回测、竞争条件和并行前缀计算。
学习目标
读完本章,你应当能够:
- 用 spawn、sync、parallel 三个关键字描述动态多线程算法,画出计算 DAG,区分逻辑并行与逻辑串行。
- 定义并计算工作量 \(T_1\)、跨度 \(T_\infty\)、加速比、并行度、松弛度,运用工作量定律、跨度定律和贪心调度界 \(T_P\le T_1/P+T_\infty\) 估算并行性能。
- 按"串联相加、并联取最大"的规则写出分治多线程算法的工作量与跨度递归式,分析并行循环、矩阵乘法和归并排序。
- 识别确定性竞争(determinacy race),理解它为什么难以通过测试发现,以及在交易系统中的典型表现。
- 用进程池并行化参数扫描回测,用模型解释实测加速比;理解并行前缀(scan)如何让 EWMA 这类线性递推也能并行。
读前导读
这一章在解决什么问题。 电脑有 8 个核,参数扫描能快 8 倍吗?买 64 核的服务器能快 64 倍吗?本章给出一套只用两个数就能回答的估算方法。工作量 \(T_1\) 是把所有活一个人干完要多久;跨度 \(T_\infty\) 是"即使人手无限,也绕不过去的那条最长先后依赖链"要多久。有 \(P\) 个处理器时,运行时间大约是
第一项是"人多力量大",第二项是"有些事只能一件接一件做"。
作为 CPA,你对这件事有切身体会:月末关账。几十个子公司的银行对账、应收应付核对可以分给很多会计同时做(工作量大,但能分摊);但"子公司出数 → 抵销内部往来 → 合并报表 → 审阅 → 签字"这一条链必须按先后顺序走(这就是跨度)。加人只能压缩第一部分;如果合并这一步本身要两天,关账就不可能少于两天。本章把这个直觉精确化,并用它分析矩阵乘法、排序等算法"值不值得并行",以及并行程序最阴险的 bug——两个线程同时改同一个数(竞争条件),在交易系统里表现为"成交回报丢了"。
需要先想起来的数学。
- 递归式及其求解。 分治算法的时间常写成 \(T(n)=T(n/2)+\Theta(\lg n)\) 这类"自己调用自己"的式子,意思是"规模 \(n\) 的时间 = 规模减半的时间 + 本层的额外开销"。解法是一层层展开再求和,或者套主定理。例:\(T(n)=T(n/2)+1\),每层加 1,减半 \(\lg n\) 次到底,所以 \(T(n)=\Theta(\lg n)\)。详见本册第 04b 章。
- 对数与 \(\lg\)。 \(\lg n\) 是以 2 为底的对数,含义是"把 \(n\) 反复对半分,几次分到 1";\(\lg1024=10\)。\(\lg^2n\) 是 \((\lg n)^2\),不是 \(\lg\lg n\)。见 第 00 册第 04 章 级数与收敛。
- 渐近记号 \(\Theta\)、\(O\)。 \(\Theta(n^2)\) 表示"随 \(n\) 增长的速度与 \(n^2\) 相同,常数不计"。比较两个算法时,先比阶数再比常数。详见本册第 03 章。
- 结合律。 运算 \(\otimes\) 满足 \((a\otimes b)\otimes c=a\otimes(b\otimes c)\),就可以随意加括号,从而把一长串运算分给多个人分别算一段、最后拼起来。加法、乘法、取最大、矩阵乘法都满足;减法不满足(\((5-3)-1\ne5-(3-1)\))。见 第 00 册第 06 章 线性代数速成。
- 取整与最大值记号。 \(\lfloor x\rfloor\) 是向下取整(\(\lfloor 3.7\rfloor=3\));\(\max(a,b)\) 取较大者;\(P\ll Q\) 读作"\(P\) 远小于 \(Q\)"。期望 \(E[\cdot]\) 在思考题 27-6 中出现,见 第 00 册第 07 章 概率中的分析工具。
怎么读这一章。 核心是 27.1.4–27.1.6:工作量、跨度、两条定律、贪心调度界、"串联相加、并联取最大",务必读懂;27.1.8 竞争条件和 27.1.9 国际象棋的例子也很重要且好读。27.1.1–27.1.3 的术语(spawn、sync、链、四类边)读一遍有印象即可。27.2 矩阵乘法和 27.3 归并排序是练习分析方法的例子,第一次可以只看每个算法的工作量、跨度和并行度结论,跳过 P-MERGE 的细节。27.5 量化实战建议精读"读结果"部分,尤其是参数扫描的过拟合提醒和"EWMA 也能并行"。
27.0 第七部分导言(简述)
原书第七部分"算法问题选编"收录了若干扩展专题:第 27 章多线程算法;第 28 章矩阵运算;第 29 章线性规划;第 30 章多项式与 FFT;第 31 章数论算法;第 32 章字符串匹配;第 33 章计算几何;第 34 章 NP 完全性;第 35 章近似算法。本册第 27、28、29 章依次对应原书第 27–29 章;第 30 章起见本册后续各章,其中原书第 31–33 章合并为本册第 31 章。
27.1 动态多线程模型
27.1.1 背景
并行计算机有多种形态:单芯片多核(多个核共享内存)、由普通机器经网络互连的集群、定制架构的超级计算机。按内存模型分为共享内存(shared memory,每个处理器可直接访问任意内存位置)和分布式内存(distributed memory,各自私有内存,靠消息通信)。多核普及之后,共享内存成为主流,本章采用这一模型。
直接用静态线程(static threading)编程——自己创建固定数量的线程、自己分配任务——困难且易错,最难的是负载均衡。于是出现了并发平台(concurrency platform):一层负责调度和管理并行资源的软件。动态多线程(dynamic multithreading)让程序员只描述"哪些部分可以并行",由平台的调度器自动做负载均衡。Cilk、OpenMP、Intel TBB、.NET Task Parallel Library 都支持这种模型;Python 的 concurrent.futures、Dask、Ray 也是同一思路。
27.1.2 三个关键字
伪代码只加三个关键字:
- spawn:放在过程调用前,表示子过程可以与调用者并行执行(嵌套并行,nested parallelism);
- sync:等待本过程派生的所有子过程完成;每个过程返回前隐式执行 sync;
- parallel:放在 for 前,表示各次迭代可以并发(并行循环)。
删去这三个关键字,就得到同一问题的串行算法,称为多线程算法的串行化(serialization)。
引例:斐波那契数。串行递归 FIB(n) 的运行时间 \(T(n)=T(n-1)+T(n-2)+\Theta(1)=\Theta(\phi^n)\),\(\phi=(1+\sqrt5)/2\),是一个很差的算法,但两次递归调用相互独立,适合用来说明概念:
P-FIB(n):
if n <= 1: return n
x = spawn P-FIB(n-1) # 子任务可与父任务并行
y = P-FIB(n-2) # 父任务继续执行
sync # 等待所有派生的子任务完成
return x + y
spawn 表达的是逻辑并行性,不是"必须"并行——哪些子计算真正同时执行,由运行时调度器决定。sync 必不可少:没有它,可能在 \(x\) 算出之前就做加法。
27.1.3 计算 DAG
把一次多线程执行建模为计算 DAG(computation dag)\(G=(V,E)\):顶点是指令,边 \((u,v)\) 表示 \(u\) 必须先于 \(v\) 执行。把不含 spawn、sync 和返回的连续指令合并成一个链(strand)。DAG 中有 \(u\) 到 \(v\) 的路径,则二者逻辑串行,否则逻辑并行。
边分四类(原书图 27.2,P-FIB(4)):延续边(continuation edge)连接同一过程实例内相继的链;派生边(spawn edge)指向被派生的子过程;调用边(call edge)指向普通调用的子过程;返回边(return edge)从子过程指回调用者 sync 之后的链。派生和调用的区别是:派生还产生一条延续边,表示其后继可以与子过程同时执行。
理想并行计算机:若干个计算能力相同的处理器,加上顺序一致(sequentially consistent)的共享内存——无论实际有多少并发读写,结果都像是把所有指令交错成一个与 DAG 偏序一致的线性顺序逐条执行。忽略调度开销。
27.1.4 工作量与跨度
- 工作量(work)\(T_1\):在单处理器上执行全部计算的时间,即所有链的时间之和。
- 跨度(span)\(T_\infty\):DAG 中最长路径上链的时间之和,即关键路径(critical path)的长度——无限多处理器时的运行时间。它可以用第 24 章的 DAG 最长路径算法在 \(\Theta(V+E)\) 内求出。
- \(T_P\):\(P\) 个处理器上的运行时间。
原书图 27.2 中 P-FIB(4) 的工作量为 17、跨度为 8(每条链计为单位时间)。
金融直觉:用一个日终批处理算一遍。任务和耗时(分钟):行情清洗 10 → 复权 5 → 然后三件互不依赖的事并行:动量因子 20、波动率因子 30、流动性因子 15 → 三者都完成后组合优化 25。工作量 \(T_1=10+5+20+30+15+25=105\) 分钟,这是一台机器从头做到尾的时间。跨度 \(T_\infty=10+5+30+25=70\) 分钟,即最长依赖链"清洗→复权→波动率→优化",这是机器再多也省不掉的时间。并行度 \(105/70=1.5\):最多快 1.5 倍;两台机器时 \(\max(105/2,70)=70\) 已经顶到跨度,再加第三台机器没有意义。
两条下界:
前者因为 \(P\) 个处理器每步至多做 \(P\) 单位工作;后者因为无限处理器的机器总能模拟 \(P\) 处理器的机器。
白话解释:工作量定律是"总活量 ÷ 人数":105 分钟的活交给 3 台机器,即使分得再均匀,也至少要 35 分钟。跨度定律是"最长的先后依赖链":链上每一步必须等上一步完成,人再多也只能一步一步走,至少 70 分钟。两个下界同时成立,所以 \(T_P\ge\max(T_1/P,\ T_\infty)\);上例中 3 台机器至少要 70 分钟,跨度是瓶颈。
加速比(speedup)\(T_1/T_P\le P\)。\(T_1/T_P=\Theta(P)\) 称线性加速,\(=P\) 称完美线性加速。
并行度(parallelism)\(T_1/T_\infty\) 有三种理解:关键路径上每一步平均可并行的工作量;无论多少处理器都不可能超过的加速比上限;能否达到完美线性加速的分界——一旦 \(P>T_1/T_\infty\),由跨度定律 \(T_1/T_P\le T_1/T_\infty<P\)。P-FIB(4) 的并行度只有 \(17/8=2.125\),处理器再多也难超过 2 倍加速。
松弛度(slackness)\((T_1/T_\infty)/P\):并行度超过处理器数的倍数。小于 1 不可能完美线性加速;越大,好的调度器越能接近完美线性加速。
27.1.5 贪心调度
调度器必须在线工作(事先不知道何时派生、何时结束)。原书分析一个简单的集中式贪心调度器(greedy scheduler):每个时间步尽可能多地分配就绪的链。若就绪链 \(\ge P\) 个,是完全步,任选 \(P\) 个执行;否则是不完全步,全部执行。
定理 27.1 在 \(P\) 处理器的理想并行计算机上,贪心调度器执行工作量 \(T_1\)、跨度 \(T_\infty\) 的计算,用时
证明:完全步每步做 \(P\) 单位工作,若完全步多于 \(\lfloor T_1/P\rfloor\) 个,总工作量将超过 \(T_1\),矛盾。对不完全步:它执行了剩余 DAG 中所有入度为 0 的链,而剩余 DAG 的最长路径必从某个入度为 0 的链开始,所以每个不完全步使剩余最长路径长度减 1,不完全步至多 \(T_\infty\) 个。\(\square\)
推导拆解:证明的思路是把每个时间步分成两类分别数。 完全步(人手全部占满):每步完成 \(P\) 单位工作,而总工作只有 \(T_1\),所以这类步最多 \(T_1/P\) 个。 不完全步(有人闲着):说明此刻所有"可以开工"的任务都在做了。剩下没做的任务里,最长依赖链的第一个任务必然是"可以开工"的(它前面没有未完成的任务),所以这一步一定推进了最长链,最长链缩短 1。最长链一开始是 \(T_\infty\),所以这类步最多 \(T_\infty\) 个。 两类相加,\(T_P\le T_1/P+T_\infty\)。对比下界 \(\max(T_1/P,T_\infty)\):两个数之和不会超过较大者的 2 倍,这就是推论 27.2"贪心调度不差于最优的 2 倍"。直观地说,只要调度器"不让任何人在有活可干时闲着",就已经足够好。
推论 27.2 贪心调度不差于最优调度的 2 倍:最优时间 \(T_P^*\ge\max(T_1/P,T_\infty)\),而 \(T_P\le T_1/P+T_\infty\le2\max(T_1/P,T_\infty)\)。
推论 27.3 若 \(P\ll T_1/T_\infty\),则 \(T_P\approx T_1/P\),加速比约为 \(P\)。经验法则:松弛度至少 10 通常就够了,此时跨度项不到每处理器工作量项的 10%。
习题 27.1-3 给出一个略紧的界 \(T_P\le(T_1-T_\infty)/P+T_\infty\)。实际系统多用随机工作窃取(work-stealing)调度(Blumofe 与 Leiserson),它是分布式的,期望时间 \(E[T_P]\le T_1/P+O(T_\infty)\)。
27.1.6 分析多线程算法
工作量就是串行化的运行时间。跨度的组合规则(原书图 27.3):
- 两段计算串联:工作量相加,跨度相加;
- 两段计算并联:工作量相加,跨度取最大值。
P-FIB:\(T_1(n)=\Theta(\phi^n)\);跨度 \(T_\infty(n)=\max(T_\infty(n-1),T_\infty(n-2))+\Theta(1)=\Theta(n)\);并行度 \(\Theta(\phi^n/n)\),随 \(n\) 急剧增长。
推导拆解:P-FIB 的跨度递归为什么取 \(\max\)?FIB\((n-1)\) 和 FIB\((n-2)\) 是并联的两段,墙上时间取决于较慢的那段,即 FIB\((n-1)\);加上 spawn 前和 sync 后的常数时间,得 \(T_\infty(n)=T_\infty(n-1)+\Theta(1)\)。每层只加常数、一共 \(n\) 层,所以 \(T_\infty(n)=\Theta(n)\)。工作量则是两段相加,与串行算法一样按 \(\phi^n\) 增长(\(\phi\approx1.618\) 是黄金分割比)。并行度 = 工作量 ÷ 跨度,指数除以线性,增长极快:\(n=20\) 时约 1095(见 27.5.1 节)。
27.1.7 并行循环
例:矩阵-向量乘法 \(y_i=\sum_ja_{ij}x_j\)。
MAT-VEC(A, x):
n = A.rows; y = new vector(n)
parallel for i = 1 to n: y[i] = 0
parallel for i = 1 to n:
for j = 1 to n:
y[i] = y[i] + a[i][j] * x[j]
return y
编译器用分治把 parallel for 实现为一棵二叉派生树(原书 MAT-VEC-MAIN-LOOP、图 27.4):区间一分为二,左半 spawn、右半直接递归,然后 sync。叶子是单次迭代。
- 工作量:串行化 \(\Theta(n^2)\)。递归派生树的内部结点比叶子少一个,每个内部结点常数工作,只增加常数因子。实践中常把若干次迭代合到一个叶子里,叫粗化(coarsening),降低开销但也降低并行度。
- 跨度:\(n\) 次迭代、第 \(i\) 次迭代跨度为 \(iter_\infty(i)\) 的并行循环,
MAT-VEC 每次外层迭代内含 \(n\) 次串行迭代,跨度 \(\Theta(n)\),并行度 \(\Theta(n^2)/\Theta(n)=\Theta(n)\)。习题 27.1-6 用并行归约把内层也并行化,并行度达到 \(\Theta(n^2/\lg n)\)。
白话解释:并行循环的跨度里为什么有 \(\Theta(\lg n)\)?因为"把 \(n\) 个迭代分发出去"本身也要时间。调度方式像淘汰赛的反向:先把 \(n\) 个任务一分为二交给两个人,每人再一分为二……分到每人一个为止,需要 \(\lg n\) 轮。1024 个迭代只需 10 轮分发,而不是逐个派发 1024 次。所以并行循环的跨度 = 分发的 \(\lg n\) + 最慢那次迭代本身的跨度。MAT-VEC 中最慢的一次迭代是内层 \(n\) 步串行求和,所以跨度是 \(\Theta(n)\)。
27.1.8 竞争条件
多线程算法若对同一输入、无论如何调度都行为相同,称为确定性的(deterministic)。本应确定的程序常因确定性竞争(determinacy race)而出错:两条逻辑并行的指令访问同一内存位置,且至少一条是写。著名的竞争 bug 包括 Therac-25 放疗机事故(3 人死亡)和 2003 年北美大停电(逾 5000 万人断电)。
RACE-EXAMPLE():
x = 0
parallel for i = 1 to 2:
x = x + 1
print x
串行化总打印 2,但并行时可能打印 1。原因是 x = x + 1 不是原子操作,而是"读入寄存器 → 加 1 → 写回"三步。若处理器 1 读到 \(x=0\),处理器 2 也读到 \(x=0\),两者各自加 1 后写回,就丢失了一次更新(原书图 27.5)。大多数交错顺序结果正确,只有少数出错,所以实验室里测几天也可能测不出来,上线后却偶发崩溃。
金融直觉:这就是会计里的"并发记账丢更新"。两个记账员同时处理同一个科目:甲看到余额 100,准备加 30;乙也看到余额 100,准备加 50。甲写回 130,乙随后写回 150,甲的 30 就凭空消失了,正确余额应是 180。手工账靠"一个科目同一时间只有一个人能改"(相当于加锁)来避免;交易系统里对持仓、可用资金、风控额度的更新同样需要这种保护,见 27.5.3 节。这类错误难查,是因为只有两人的读写恰好交错时才出错,绝大多数时候账是对的。
原书的约定是让并行的链相互独立:parallel for 的各次迭代互不干扰;spawn 与对应 sync 之间,子任务与父任务(及其他子任务)的代码互不干扰。反例 MAT-VEC-WRONG 把内层循环也改为 parallel for 以求 \(\Theta(\lg n)\) 跨度,结果所有 \(j\) 同时更新 \(y_i\),产生竞争。
27.1.9 国际象棋程序的教训
★Socrates 国际象棋程序在 32 处理器的机器上开发,最终要跑在 512 处理器的超级计算机上。某项"优化"让 32 处理器上的基准从 65 秒降到 40 秒。但原版 \(T_1=2048\)、\(T_\infty=1\),优化版 \(T'_1=1024\)、\(T'_\infty=8\)。用 \(T_P\approx T_1/P+T_\infty\):
| 原版 | "优化"版 | |
|---|---|---|
| \(P=32\) | \(2048/32+1=65\) | \(1024/32+8=40\) |
| \(P=512\) | \(2048/512+1=5\) | \(1024/512+8=10\) |
在 512 个处理器上,"优化"版反而慢一倍——跨度 8 在大处理器数时成为主导项。两版在 \(P=1024/7\approx146\) 时一样快(习题 27.1-9)。开发者据此放弃了这项优化。教训:用工作量和跨度外推性能,比直接用某个处理器数上的实测时间更可靠。
27.2 多线程矩阵乘法
27.2.1 并行化三重循环
P-SQUARE-MATRIX-MULTIPLY(A, B):
n = A.rows; C = new n×n matrix
parallel for i = 1 to n:
parallel for j = 1 to n:
c[i][j] = 0
for k = 1 to n:
c[i][j] += a[i][k] * b[k][j]
return C
工作量 \(\Theta(n^3)\);跨度 \(\Theta(\lg n)+\Theta(\lg n)+\Theta(n)=\Theta(n)\);并行度 \(\Theta(n^2)\)。最内层 \(k\) 循环不能直接改成 parallel for——所有 \(k\) 都写 \(c_{ij}\),会产生竞争。
27.2.2 分治版本
把 \(A,B,C\) 各分成四个 \(n/2\times n/2\) 的块:
P-MATRIX-MULTIPLY-RECURSIVE 并行派生 8 个 \(n/2\) 阶乘法(前 4 个写入 \(C\),后 4 个写入临时矩阵 \(T\)),sync 后再用双重 parallel for 把 \(T\) 加到 \(C\)。
- 工作量:\(M_1(n)=8M_1(n/2)+\Theta(n^2)=\Theta(n^3)\)。
- 跨度:8 个递归并行,取其一;加法循环跨度 \(\Theta(\lg n)\):
推导拆解:把递归式逐层展开。第 1 层额外开销 \(\lg n\),第 2 层规模减半,开销 \(\lg(n/2)=\lg n-1\),第 3 层 \(\lg n-2\)……一直到规模为 1,共 \(\lg n\) 层。总和是 \(\lg n+(\lg n-1)+\cdots+1=\frac{\lg n(\lg n+1)}{2}=\Theta(\lg^2n)\)。这和年金现值里"逐期加总再用公式求和"是同一类操作,只是这里每期的金额递减 1。以 \(n=1024\) 为例,\(\lg n=10\),跨度约 \(10\times11/2=55\) 个单位,而工作量约 \(10^9\),并行度极高。
- 并行度 \(\Theta(n^3/\lg^2n)\),非常高。代价是临时矩阵 \(T\)。思考题 27-2 去掉 \(T\)(先并行算 4 个乘积累加到 \(C\),sync,再算另外 4 个),跨度变成 \(\Theta(n)\);对 \(1000\times1000\) 的矩阵,并行度仍约 \(10^6\),远超实际处理器数。
多线程 Strassen:分块、构造 10 个和/差矩阵 \(S_1..S_{10}\)(\(\Theta(n^2)\) 工作量、\(\Theta(\lg n)\) 跨度)、并行递归计算 7 个乘积、再组合。工作量 \(\Theta(n^{\lg7})\),跨度仍为 \(\Theta(\lg^2n)\),并行度 \(\Theta(n^{\lg7}/\lg^2n)\)。
并行 BLAS(如 OpenBLAS、MKL)里的分块矩阵乘法正是这个结构,这也是协方差矩阵估计、因子暴露矩阵乘法等运算能在多核上接近线性加速的原因。
27.3 多线程归并排序
27.3.1 只并行递归不够
MERGE-SORT'(A, p, r):
if p < r:
q = floor((p + r) / 2)
spawn MERGE-SORT'(A, p, q)
MERGE-SORT'(A, q+1, r)
sync
MERGE(A, p, q, r) # 串行合并,Θ(n) 工作量与跨度
工作量 \(\Theta(n\lg n)\);跨度 \(MS'_\infty(n)=MS'_\infty(n/2)+\Theta(n)=\Theta(n)\);并行度只有 \(\Theta(\lg n)\)——排 1000 万个元素,也只能在几个处理器上线性加速。瓶颈是串行的 MERGE。
27.3.2 并行合并 P-MERGE
把 \(T[p_1..r_1]\)(长 \(n_1\))与 \(T[p_2..r_2]\)(长 \(n_2\))合并到 \(A[p_3..]\),不妨设 \(n_1\ge n_2\)(原书图 27.6):
- 取较长子数组的中位元素 \(x=T[q_1]\),\(q_1=\lfloor(p_1+r_1)/2\rfloor\);
- 在另一个子数组中二分查找 \(q_2\),使 \(x\) 插入 \(T[q_2-1]\) 与 \(T[q_2]\) 之间后仍有序;
- \(q_3=p_3+(q_1-p_1)+(q_2-p_2)\),令 \(A[q_3]=x\);
- 并行递归合并 \(x\) 左边的两段与 \(x\) 右边的两段。
P-MERGE(T, p1, r1, p2, r2, A, p3):
n1 = r1 - p1 + 1; n2 = r2 - p2 + 1
if n1 < n2: 交换 (p1,r1,n1) 与 (p2,r2,n2) # 保证 n1 >= n2
if n1 == 0: return # 两个都为空
q1 = floor((p1 + r1) / 2)
q2 = BINARY-SEARCH(T[q1], T, p2, r2)
q3 = p3 + (q1 - p1) + (q2 - p2)
A[q3] = T[q1]
spawn P-MERGE(T, p1, q1-1, p2, q2-1, A, p3)
P-MERGE(T, q1+1, r1, q2, r2, A, q3+1)
sync
分析(\(n=n_1+n_2\)):
- 跨度:\(n_2\le n/2\),最坏时一个递归合并 \(\lfloor n_1/2\rfloor\) 个与全部 \(n_2\) 个元素,\(\lfloor n_1/2\rfloor+n_2\le n/2+n_2/2\le3n/4\)。加上二分查找:
- 工作量:每个元素都要复制,\(\Omega(n)\);上界由递归式 \(PM_1(n)=PM_1(\alpha n)+PM_1((1-\alpha)n)+O(\lg n)\)(\(1/4\le\alpha\le3/4\))用代入法 \(PM_1(n)\le c_1n-c_2\lg n\) 证得 \(\Theta(n)\)。
- 并行度 \(\Theta(n/\lg^2n)\)。
27.3.3 多线程归并排序 P-MERGE-SORT
两个子数组递归并行排序到临时数组,再用 P-MERGE 合并。工作量 \(\Theta(n\lg n)\),跨度 \(PMS_\infty(n)=PMS_\infty(n/2)+\Theta(\lg^2n)=\Theta(\lg^3n)\),并行度 \(\Theta(n/\lg^2n)\),远好于 \(\Theta(\lg n)\)。实践中应粗化基本情形:规模足够小时改用串行排序,牺牲一点并行度换取更小的常数。
27.4 思考题选讲
- 27-1 粒度:把 \(n\) 次迭代按粒度 \(g\) 切成 \(n/g\) 段,用一个串行 for 循环依次 spawn 各段。跨度 \(\Theta(n/g+g)\),取 \(g=\sqrt n\) 时并行度最大为 \(\Theta(\sqrt n)\)。若改用分治派生(parallel for 的标准实现),跨度可降到 \(\Theta(\lg n+g)\)。量化实践中的
chunksize、批大小就是这个粒度参数。 - 27-2 节省临时空间:见 27.2.2 节。
- 27-3 多线程矩阵算法:多线程化第 28 章的 LU、LUP 分解和 LUP-SOLVE,以及对称正定矩阵的分块求逆。
- 27-4 并行归约与前缀计算:\(\otimes\) 为满足结合律的运算。归约 \(x_1\otimes\cdots\otimes x_n\) 可用分治在 \(\Theta(n)\) 工作量、\(\Theta(\lg n)\) 跨度内完成。前缀计算(scan)\(y_i=x_1\otimes\cdots\otimes x_i\) 的串行循环有依赖,不能直接 parallel for。三种并行方案:每个位置独立归约(工作量 \(\Theta(n^2)\),跨度 \(\Theta(\lg n)\));递归求左右两半再合并(\(\Theta(n\lg n)\)、\(\Theta(\lg^2n)\));两遍扫描(上行在内部结点记录左半区间的归约值,下行把"左边所有元素的归约"传给右半),工作量 \(\Theta(n)\)、跨度 \(\Theta(\lg n)\)。GPU 和向量化库的
cumsum正是用两遍法实现的。 - 27-5 模板计算:\(A[i,j]\) 只依赖左上方已算好的元素(如最长公共子序列的 DP 表)。2×2 分块(先 \(A_{11}\),再并行 \(A_{12}\)、\(A_{21}\),最后 \(A_{22}\))跨度 \(\Theta(n^{\lg3})\);按反对角线"波前"推进可以达到 \(\Theta(n)\) 跨度、\(\Theta(n)\) 并行度。
- 27-6 随机化多线程算法:工作量、跨度定律改为期望形式。加速比应定义为 \(E[T_1]/E[T_P]\),而不是 \(E[T_1/T_P]\)——原书用一个 1%/99% 的例子说明后者会严重误导。
27.5 量化实战
27.5.1 工作量、跨度与贪心调度的数值验证
下面把 P-FIB(n) 的计算 DAG 显式建出来:每个非基本实例有三条链(spawn 之前、spawn 与 sync 之间、sync 之后),基本实例一条链。用第 24 章的 DAG 最长路径求跨度,再模拟贪心调度器,检验定理 27.1。然后复算国际象棋程序的例子,并实现思考题 27-4 的两遍扫描前缀计算。
import math
from collections import defaultdict
# ---------- 1. 把 P-FIB(n) 的计算 DAG 显式建出来,算工作量、跨度,并模拟贪心调度 ----------
def build_pfib_dag(n):
succ, indeg, cnt = defaultdict(list), defaultdict(int), [0]
def new():
cnt[0] += 1; return cnt[0]
def edge(a, b):
succ[a].append(b); indeg[b] += 1
def pfib(n, entry):
"""entry: 进入本实例前的链;返回本实例最后一条链。"""
if n <= 1:
a = new(); edge(entry, a) if entry else None; return a
a = new(); edge(entry, a) if entry else None # spawn 之前的链
x_end = pfib(n - 1, a) # 派生边 a -> 子实例
b = new(); edge(a, b) # 延续边 a -> b(与子实例并行)
y_end = pfib(n - 2, b) # 调用边 b -> 子实例
c = new(); edge(x_end, c); edge(y_end, c) # sync 之后的链(返回边汇合)
return c
pfib(n, None)
return cnt[0], succ, indeg
def span(V, succ, indeg):
indeg = dict(indeg); depth = {v: 1 for v in range(1, V + 1)}
ready = [v for v in range(1, V + 1) if indeg.get(v, 0) == 0]
while ready: # 拓扑序上求最长路径(24.2 节)
u = ready.pop()
for v in succ[u]:
depth[v] = max(depth[v], depth[u] + 1); indeg[v] -= 1
if indeg[v] == 0: ready.append(v)
return max(depth.values())
def greedy_schedule(V, succ, indeg, P):
indeg = dict(indeg); ready = [v for v in range(1, V + 1) if indeg.get(v, 0) == 0]; t = 0
while ready:
run, ready = ready[:P], ready[P:] # 完全步取 P 个,不完全步全取
t += 1
for u in run:
for v in succ[u]:
indeg[v] -= 1
if indeg[v] == 0: ready.append(v)
return t
V, succ, indeg = build_pfib_dag(4)
print(f"P-FIB(4): 工作量 T1 = {V},跨度 T_inf = {span(V, succ, indeg)}")
V, succ, indeg = build_pfib_dag(20)
T1, Tinf = V, span(V, succ, indeg)
print(f"P-FIB(20): T1 = {T1},T_inf = {Tinf},并行度 = {T1 / Tinf:.0f}")
for P in (1, 4, 16, 64, 256, 1024):
TP = greedy_schedule(V, succ, indeg, P)
print(f" P={P:5d}: 贪心调度 T_P = {TP:6d},下界 max(T1/P, T_inf) = {max(math.ceil(T1/P), Tinf):6d},"
f"上界 T1/P + T_inf = {T1 / P + Tinf:8.1f},加速比 {T1 / TP:6.1f}")
# ---------- 2. 国际象棋程序的教训:用 T_P ≈ T1/P + T_inf 外推 ----------
for P in (32, 512):
print(f"P={P}: 原版 {2048 / P + 1:.0f} 秒,'优化'版 {1024 / P + 8:.0f} 秒")
print("两版同样快的处理器数 P =", 1024 / 7)
# ---------- 3. 思考题 27-4:两遍扫描的并行前缀(P-SCAN-3),并统计工作量与跨度 ----------
def p_scan(x, op, identity):
n = len(x); t = [None] * n; y = [None] * n; work = [0]
def up(i, j): # 返回 (区间归约值, 跨度);t[k] 存左半区间的归约
if i == j: return x[i], 1
k = (i + j) // 2
(l, sl), (r, sr) = up(i, k), up(k + 1, j) # 两个递归可并行:跨度取 max
t[k] = l; work[0] += 1
return op(l, r), max(sl, sr) + 1
def down(v, i, j): # v = x[0] ⊗ ... ⊗ x[i-1]
if i == j:
y[i] = op(v, x[i]); work[0] += 1; return 1
k = (i + j) // 2
work[0] += 1
return max(down(v, i, k), down(op(v, t[k]), k + 1, j)) + 1
_, s_up = up(0, n - 1)
s_down = down(identity, 0, n - 1)
return y, work[0], s_up + s_down
import numpy as np
rng = np.random.default_rng(0)
for n in (1024, 16384):
x = list(rng.normal(size=n))
y, W, S = p_scan(x, lambda a, b: a + b, 0.0)
print(f"n={n}: 前缀和正确 {np.allclose(y, np.cumsum(x))},工作量 {W}(≈{W / n:.1f}n),"
f"跨度 {S}(≈{S / math.log2(n):.1f} lg n)")
# EWMA 是线性递推 s_t = a*s_{t-1} + (1-a)*x_t。仿射映射 s -> A*s + B 的复合满足结合律,
# 所以 EWMA 也能用并行前缀计算:先 f 后 g 的复合为 (A_g*A_f, A_g*B_f + B_g)
a = 0.94; r = rng.normal(0, 0.01, 5000); x2 = r ** 2 # RiskMetrics 式方差递推
maps = [(a, (1 - a) * v) for v in x2]
compose = lambda f, g: (g[0] * f[0], g[0] * f[1] + g[1])
y, W, S = p_scan(maps, compose, (1.0, 0.0))
s0 = x2[:20].mean()
ewma_scan = np.array([A * s0 + B for A, B in y])
ewma_loop = np.empty(len(x2)); s = s0
for i, v in enumerate(x2):
s = a * s + (1 - a) * v; ewma_loop[i] = s
print(f"EWMA 方差:并行前缀与串行递推一致 {np.allclose(ewma_scan, ewma_loop)},跨度 {S}(串行需 {len(x2)} 步)")
输出:
P-FIB(4): 工作量 T1 = 17,跨度 T_inf = 8
P-FIB(20): T1 = 43781,T_inf = 40,并行度 = 1095
P= 1: 贪心调度 T_P = 43781,下界 max(T1/P, T_inf) = 43781,上界 T1/P + T_inf = 43821.0,加速比 1.0
P= 4: 贪心调度 T_P = 10949,下界 max(T1/P, T_inf) = 10946,上界 T1/P + T_inf = 10985.2,加速比 4.0
P= 16: 贪心调度 T_P = 2745,下界 max(T1/P, T_inf) = 2737,上界 T1/P + T_inf = 2776.3,加速比 15.9
P= 64: 贪心调度 T_P = 698,下界 max(T1/P, T_inf) = 685,上界 T1/P + T_inf = 724.1,加速比 62.7
P= 256: 贪心调度 T_P = 191,下界 max(T1/P, T_inf) = 172,上界 T1/P + T_inf = 211.0,加速比 229.2
P= 1024: 贪心调度 T_P = 68,下界 max(T1/P, T_inf) = 43,上界 T1/P + T_inf = 82.8,加速比 643.8
P=32: 原版 65 秒,'优化'版 40 秒
P=512: 原版 5 秒,'优化'版 10 秒
两版同样快的处理器数 P = 146.28571428571428
n=1024: 前缀和正确 True,工作量 3070(≈3.0n),跨度 22(≈2.2 lg n)
n=16384: 前缀和正确 True,工作量 49150(≈3.0n),跨度 30(≈2.1 lg n)
EWMA 方差:并行前缀与串行递推一致 True,跨度 28(串行需 5000 步)
读结果:
- P-FIB(4) 的工作量 17、跨度 8 与原书图 27.2 一致。P-FIB(20) 的跨度为 40(\(T_\infty(n)=2n\)),并行度约 1095。
- 贪心调度的实测时间始终夹在下界 \(\max(T_1/P,T_\infty)\) 与上界 \(T_1/P+T_\infty\) 之间。\(P\le64\) 时松弛度 \(\ge17\),加速比几乎等于 \(P\);\(P=1024\) 时松弛度约 1.07,加速比只有 644——这就是推论 27.3"松弛度至少 10"的含义。
- 两遍扫描的工作量约 \(3n\)、跨度约 \(2\lg n\),与 \(\Theta(n)\)、\(\Theta(\lg n)\) 相符。
- EWMA 也能并行。\(s_t=a\,s_{t-1}+(1-a)x_t\) 看上去是典型的串行递推,但每一步都是一个仿射映射 \(s\mapsto A s+B\),而仿射映射的复合满足结合律。把映射序列做前缀"复合",就得到 \(s_t=A_{1..t}s_0+B_{1..t}\)。5000 步的递推在并行模型下跨度只有 28。同样的技巧适用于任何一阶线性递推:指数加权均值与方差、AR(1) 滤波、复利累积、卡尔曼滤波中的某些线性部分(第 06 册)。
推导拆解:为什么一个"每步都依赖上一步"的递推能并行?用两步把仿射映射的复合算出来。第 \(t\) 步是 \(s\mapsto a s+(1-a)x_t\),记作 \((A_t,B_t)=(a,(1-a)x_t)\)。先做第 1 步再做第 2 步:\(s_2=a\big(a s_0+(1-a)x_1\big)+(1-a)x_2=a^2s_0+\big[a(1-a)x_1+(1-a)x_2\big]\),即复合后的映射是 \((a^2,\ a(1-a)x_1+(1-a)x_2)\),形状不变,仍是"乘一个数再加一个数"。 关键在于复合满足结合律:先把第 1–2 步合成一个映射、第 3–4 步合成一个映射,再把这两个合起来,与按顺序一步步复合结果相同。于是可以像淘汰赛一样两两合并:第一轮 2500 对同时合并,第二轮 1250 对……约 \(\lg5000\approx13\) 轮就得到全程映射,再用"下行"一遍把每个中间时点的结果分发出去,所以跨度只有几十步。这与复利计算相同:\((1+r_1)(1+r_2)(1+r_3)(1+r_4)\) 可以先算前两期和后两期再相乘,结果不变。
27.5.2 并行回测:参数扫描
参数网格搜索是最典型的"各次迭代相互独立"的 parallel for。下面在 500 只模拟股票、10 年日线上扫描 108 组双均线参数,每组还做 2000 次自助法估计夏普比率的标准误。用进程池(而不是线程)并行,并与 \(T_P\approx T_1/P+T_\infty\) 比较。
import os
os.environ.setdefault("OMP_NUM_THREADS", "1") # 每个进程只用一个 BLAS 线程,避免超额订阅
import time, threading
import numpy as np
from concurrent.futures import ProcessPoolExecutor
def make_prices(n_assets=500, n_days=2520, seed=1):
rng = np.random.default_rng(seed)
r = rng.normal(0.0002, 0.015, (n_days, n_assets))
return np.cumprod(1 + r, axis=0)
def backtest(args):
"""一组均线参数 (fast, slow) 在全部资产上的回测;各组参数之间完全独立。"""
fast, slow = args
P = PRICES
c = np.cumsum(P, axis=0); c = np.vstack([np.zeros((1, P.shape[1])), c])
ma_f = (c[slow:] - c[slow - fast:-fast]) / fast # 与 ma_s 对齐到同一天
ma_s = (c[slow:] - c[:-slow]) / slow
pos = np.sign(ma_f - ma_s)[:-1] # 今天的信号,明天持有
ret = P[slow:] / P[slow - 1:-1] - 1
pnl = (pos * ret).mean(axis=1) # 等权组合的日收益
sr = pnl.mean() / pnl.std() * np.sqrt(252)
rng = np.random.default_rng(fast * 1000 + slow) # 每组参数固定种子:串行、并行结果相同
idx = rng.integers(0, len(pnl), (2000, len(pnl))) # 自助法估计夏普比率的标准误
boot = pnl[idx]; sr_boot = boot.mean(1) / boot.std(1) * np.sqrt(252)
return fast, slow, sr, sr_boot.std()
def init(seed):
global PRICES
PRICES = make_prices(seed=seed)
GRID = [(f, s) for f in range(5, 65, 5) for s in range(80, 260, 20) if f < s]
if __name__ == "__main__":
t0 = time.perf_counter(); init(1); t_load = time.perf_counter() - t0 # 串行部分:数据准备
t0 = time.perf_counter(); res1 = [backtest(g) for g in GRID]; T1 = time.perf_counter() - t0
per_task = T1 / len(GRID)
print(f"{len(GRID)} 组参数,串行 T1 = {T1:.2f} s(每组 {1e3 * per_task:.0f} ms),数据准备 {t_load:.2f} s")
for P in (2, 4, 8):
with ProcessPoolExecutor(P, initializer=init, initargs=(1,)) as ex:
list(ex.map(backtest, GRID[:P])) # 预热:进程启动与数据生成
t0 = time.perf_counter(); resP = list(ex.map(backtest, GRID)); TP = time.perf_counter() - t0
model = T1 / P + per_task # 跨度 ≈ 单个任务时长
same = np.allclose([x[2] for x in res1], [x[2] for x in resP])
print(f"P={P}: T_P = {TP:.2f} s,加速比 {T1 / TP:.2f},模型 T1/P+T_inf = {model:.2f} s,结果一致 {same}")
best = max(res1, key=lambda x: x[2])
print(f"最优参数 fast={best[0]}, slow={best[1]},年化夏普 {best[2]:.2f} ± {best[3]:.2f}(随机游走数据,真实夏普为 0)")
# 竞争条件:两个线程同时给同一持仓"读-改-写"
def add_fills(book, n, lock=None):
for _ in range(n):
if lock: lock.acquire()
q = book["AAPL"] # 读
time.sleep(0) # 让出 CPU,放大交错的概率
book["AAPL"] = q + 1 # 写
if lock: lock.release()
for lock in (None, threading.Lock()):
book = {"AAPL": 0}
ts = [threading.Thread(target=add_fills, args=(book, 2000, lock)) for _ in range(2)]
[t.start() for t in ts]; [t.join() for t in ts]
print(("有锁" if lock else "无锁"), ":两个线程各成交 2000 股,记录的持仓 =", book["AAPL"])
输出(计时随机器和负载而异):
108 组参数,串行 T1 = 2.22 s(每组 21 ms),数据准备 0.02 s
P=2: T_P = 1.16 s,加速比 1.91,模型 T1/P+T_inf = 1.13 s,结果一致 True
P=4: T_P = 0.65 s,加速比 3.43,模型 T1/P+T_inf = 0.57 s,结果一致 True
P=8: T_P = 0.39 s,加速比 5.63,模型 T1/P+T_inf = 0.30 s,结果一致 True
最优参数 fast=20, slow=80,年化夏普 0.47 ± 0.33(随机游走数据,真实夏普为 0)
无锁 :两个线程各成交 2000 股,记录的持仓 = 2000
有锁 :两个线程各成交 2000 股,记录的持仓 = 4000
读结果:
- 模型与实测。这个计算的"DAG"很简单:108 个独立任务并联,工作量 \(T_1\) 是所有任务之和,跨度约等于一个任务的时长(21 ms),并行度约 108。\(P=2,4\) 时实测与 \(T_1/P+T_\infty\) 很接近;\(P=8\) 时加速比 5.6,低于模型的 7.4。差距来自理想模型忽略的东西:自助法的随机下标访问受内存带宽限制,8 个进程同时访存时互相争抢;进程间传参和返回结果的开销;以及现代 CPU 的大小核、睿频差异。理想模型给的是上界,越靠近硬件极限,偏差越大。
- 结果可复现。每组参数用自己的固定种子,所以串行和并行得到完全相同的结果。并行回测中,随机数种子要跟着任务走,而不是跟着进程走,否则结果依赖于调度顺序,变成"非确定性"的。
- 数据挖掘偏差。数据是无漂移差异的随机游走,真实夏普为 0,但 108 组参数里最好的那组年化夏普 0.47,约 1.4 倍标准误。并行让参数扫描变得便宜,也让过拟合变得更便宜——多重检验校正(第 03 册)在这里必不可少。
- Python 的特殊情况。CPython 默认构建有全局解释器锁(GIL),纯 Python 的 CPU 密集代码用线程不能并行,所以这里用进程池。numpy 的大部分运算会释放 GIL,因此"线程 + 大块 numpy 运算"也能获得一些并行;但每个进程内若 BLAS 也开多线程,就会出现"进程数 × BLAS 线程数"远超核数的超额订阅,代码第一行把
OMP_NUM_THREADS设为 1 就是为了避免它。 - 粒度。每个任务 21 ms,远大于进程间通信开销,粒度合适。若把任务切成每个 0.1 ms,调度和通信开销就会吃掉收益——这是思考题 27-1 粒度问题的实际版本,
ex.map的chunksize参数就是用来调粒度的。
27.5.3 竞争条件:丢失的成交
上面输出的最后两行演示了 RACE-EXAMPLE 在交易系统里的原型:两个线程各自处理 2000 笔成交回报,每笔把持仓加 1。没有锁时,"读—改—写"被交错,记录的持仓只有 2000 股——一半的成交"丢了";加锁后正确得到 4000 股。代码里的 time.sleep(0) 人为放大了交错的概率,真实系统中这种 bug 可能几天才出现一次,却足以造成持仓、资金和风控计数的错误。
实盘系统常见的对策:
- 单写者原则:每份状态(某账户的持仓、某品种的订单簿)只由一个线程修改,其他线程通过消息队列提交请求。许多撮合引擎和交易网关采用单线程事件循环,正是为了从设计上消除竞争。
- 原子操作与锁:必须共享时,用原子指令或锁保护"读—改—写"的整个区间。
- 不可变数据与函数式更新:回测和研究代码中,让并行任务只读共享数据、各自返回结果、最后由主进程汇总——这正是 27.5.2 的写法,它在结构上满足原书"并行的链相互独立"的约定。
27.5.4 批处理流水线的跨度
日终批处理(行情清洗 → 复权 → 因子计算 → 风险模型 → 组合优化 → 生成订单)是一个任务 DAG。工作量是所有任务时间之和,跨度是关键路径(第 24 章 24.4 节)。加机器只能压缩工作量项 \(T_1/P\),压不动跨度。如果关键路径上有一个串行的 30 分钟任务(例如单线程的全市场协方差估计),无论加多少机器,整批任务都不可能在 30 分钟内完成——这时该做的是并行化关键路径上的那一步,而不是加机器。
本章小结
动态多线程只用 spawn、sync、parallel 三个关键字表达逻辑并行,删掉它们就是串行算法。计算被建模为由链组成的 DAG:工作量 \(T_1\) 是总时间,跨度 \(T_\infty\) 是关键路径长度,并行度 \(T_1/T_\infty\) 是可能的最大加速比。工作量定律 \(T_P\ge T_1/P\) 与跨度定律 \(T_P\ge T_\infty\) 给出下界,贪心调度达到 \(T_P\le T_1/P+T_\infty\),在最优的 2 倍以内;松弛度在 10 以上就接近线性加速。分析时串联跨度相加、并联取最大,并行循环的跨度为 \(\Theta(\lg n)\) 加单次迭代的最大跨度。必须避免确定性竞争。代表算法:矩阵乘法(三重循环并行度 \(\Theta(n^2)\),分治 \(\Theta(n^3/\lg^2n)\)),归并排序(串行合并只有 \(\Theta(\lg n)\),并行合并后 \(\Theta(n/\lg^2n)\))。量化中,参数扫描、蒙特卡洛是高并行度的 parallel for;任何满足结合律的累积运算(包括 EWMA 这样的线性递推)都能用并行前缀计算;交易系统中的共享状态更新是竞争 bug 的高发区。
| 概念 / 结论 | 公式或要点 |
|---|---|
| 工作量 / 跨度 | \(T_1\) = 所有链时间之和;\(T_\infty\) = 关键路径长度 |
| 工作量定律 | \(T_P\ge T_1/P\) |
| 跨度定律 | \(T_P\ge T_\infty\) |
| 并行度 | \(T_1/T_\infty\),加速比的上限 |
| 松弛度 | \(T_1/(PT_\infty)\),\(\ge10\) 时近似线性加速 |
| 贪心调度 | \(T_P\le T_1/P+T_\infty\le2T_P^*\) |
| 组合规则 | 串联:工作量、跨度都相加;并联:工作量相加、跨度取最大 |
| 并行循环跨度 | \(\Theta(\lg n)+\max_i iter_\infty(i)\) |
| 分治矩阵乘法 | \(T_1=\Theta(n^3)\),\(T_\infty=\Theta(\lg^2n)\) |
| P-MERGE | \(T_1=\Theta(n)\),\(T_\infty=\Theta(\lg^2n)\) |
| P-MERGE-SORT | \(T_1=\Theta(n\lg n)\),\(T_\infty=\Theta(\lg^3n)\) |
| 并行前缀(两遍法) | \(T_1=\Theta(n)\),\(T_\infty=\Theta(\lg n)\),要求运算满足结合律 |
| 确定性竞争 | 两条逻辑并行的指令访问同一位置,至少一条是写 |
练习
基础
- 把 P-FIB 中第二个递归调用也改成 spawn,对渐近工作量、跨度和并行度有什么影响?(原书 27.1-1。答:没有影响。)
- 画出 P-FIB(5) 的计算 DAG,求工作量、跨度和并行度,并给出 3 个处理器上的一个贪心调度。(原书 27.1-2。可用 27.5.1 的代码核对:工作量 29,跨度 10。)
- Karan 教授声称某程序 \(T_4=80\)、\(T_{10}=42\)、\(T_{64}=10\) 秒。用工作量定律、跨度定律和习题 27.1-3 的界 \(T_P\le(T_1-T_\infty)/P+T_\infty\) 证明这些数据不可能同时成立。(原书 27.1-5。)
- 写出原地转置 P-TRANSPOSE(两层 parallel for 交换 \(a_{ij}\) 与 \(a_{ji}\))的工作量、跨度与并行度;若内层改为普通 for 呢?(原书 27.1-7、27.1-8。)
- 设计工作量 \(\Theta(n^2)\)、并行度 \(\Theta(n^2/\lg n)\) 的矩阵-向量乘法。(原书 27.1-6。提示:内层用并行归约。)
进阶
- 多线程化 Floyd-Warshall:\(k\) 循环串行,\(i,j\) 两层并行。证明工作量 \(\Theta(n^3)\)、跨度 \(\Theta(n\lg n)\),并说明为什么 \(k\) 循环不能并行。(原书 27.2-6。)
- 证明思考题 27-1 中,按粒度 \(g\) 串行派生各段时跨度为 \(\Theta(n/g+g)\),并求最优粒度。在 27.5.2 的代码中把每个任务切成更小的块,观察加速比如何随块大小变化。
- 实现思考题 27-4 的 P-SCAN-2(递归求左右两半前缀,再用 parallel for 把左半最后值并入右半),统计工作量与跨度,与 P-SCAN-3 比较。
- 用并行前缀实现"带交易成本的累计净值":每日净值 \(V_t=V_{t-1}(1+r_t)-c_t\)。找出对应的结合运算,并用 27.5.1 的
p_scan验证。 - 思考题 27-6(b):1% 的情况 \(T_1=10^4\)、\(T_{10000}=1\),99% 的情况 \(T_1=T_{10000}=10^9\)。分别计算 \(E[T_1/T_P]\) 和 \(E[T_1]/E[T_P]\),说明为什么后者才是合理的加速比定义。
原书推荐习题:27.1-5、27.1-6、27.1-9、27.2-3、27.2-6、27.3-3,思考题 27-1、27-4(强烈推荐)、27-5。
原书对照
页码换算:原书页码 = PDF 页码 − 21。
| 本章小节 | 原书章节 | PDF 页码(原书页码) |
|---|---|---|
| 27.0 第七部分导言 | Part VII Introduction | p.789–792(768–771) |
| 27.1 动态多线程模型 | 第 27 章导言;27.1 The basics of dynamic multithreading | p.793–813(772–792) |
| 27.2 多线程矩阵乘法 | 27.2 Multithreaded matrix multiplication | p.813–818(792–797) |
| 27.3 多线程归并排序 | 27.3 Multithreaded merge sort | p.818–826(797–805) |
| 27.4 思考题选讲 | Problems 27-1 ~ 27-6;Chapter notes | p.826–833(805–812) |