目标
分布式 MPI PCG 求解器
一句话总结:实现一个分布式 MPI 预处理共轭梯度(PCG)求解器,求解 1D Poisson 方程 Lx = b。你只能修改 distributed_pcg.cpp,在 1-64 秩上测试扩展性(N = 2²⁰ – 2²⁶)。
这是课程中最硬核的一次作业——综合了分布式内存编程(MPI)、稀疏线性代数、数值方法和性能优化。PCG 是科学计算中最重要的迭代求解器之一,理解它对你未来的 HPC 工作大有裨益。
| 任务 | 对应课程 | 核心概念 |
|---|---|---|
| 稀疏矩阵格式转换 | Lecture 16 | map → CSR/CSC 转换 |
| 分布式布局设计 | Lecture 7-8 | 行划分、halo 区域 |
| 并行化 PCG 操作 | Lecture 7-8, 16 | SpMV、内积、向量更新 |
| 优化 | Lecture 7-8 | 中间值重用、通信优化 |
| 扩展性研究 | Lecture 7-8 | 1-64 秩强扩展实验 |
| 报告 | — | 策略、性能、挑战分析 |
MPI_Send / MPI_Recv 点对点通信MPI_Allreduce 全局归约(内积计算必需)MPI_Allgather / MPI_Scatter 数据分发MPI_Isend / MPI_Irecv(用于通信隐藏)先回顾串行 PCG 的伪代码:
r = b - A*x # 残差
z = M^{-1} * r # 预处理
p = z # 搜索方向
while not converged:
alpha = (r·z) / (p·A*p)
x = x + alpha * p
r_new = r - alpha * A*p
z_new = M^{-1} * r_new
beta = (r_new·z_new) / (r·z)
p = z_new + beta * p
r = r_new; z = z_new
每步的主要计算:1 次 SpMV(A*p)、2 次内积(dot)、3 次向量更新(axpy)。
输入可能是 map 格式(键值对 (i,j) → value),需要转为 CSR(Compressed Sparse Row):
1D Poisson 矩阵是三对角矩阵:每行最多 3 个非零元素。
将 N 个未知数划分到 P 个秩(N 必须是 P 的倍数):
MPI_Allreduce 求和设计实验矩阵:
N 不是 P 的倍数: 题目明确要求 N 必须是秩数的倍数。如果不是,行划分会不均匀,导致负载不均衡甚至错误。
Halo 交换遗漏: SpMV 时,边界行需要相邻秩的第一个/最后一个元素。忘记 halo 交换会导致错误结果。
Allreduce 性能瓶颈: 每次迭代有 2 次全局归约(内积),这是主要的通信开销。用 MPI_Allreduce 而不是手动 send/recv。
收敛条件错误: PCG 的收敛判断通常是 ||r|| / ||b|| < tol。注意计算全局范数(需要 Allreduce),而不是用局部残差。
内存泄漏: MPI 分配的请求对象(MPI_Request)和 buffer 必须正确释放,尤其在迭代循环中。
在 [0,1] 区间上用中心差分离散 -u'' = f,步长 h = 1/(N+1):
# 三对角矩阵:主对角线 = 2/h²,次对角线 = -1/h²
A = 1/h² × [ 2 -1 0 0 ... 0 ]
[-1 2 -1 0 ... 0 ]
[ 0 -1 2 -1 ... 0 ]
[ ... ]
[ 0 ... 0 -1 2 -1 ]
[ 0 ... 0 0 -1 2 ]
| 预处理器 | 说明 | 代价 |
|---|---|---|
| Jacobi (Diagonal) | M = diag(A),M⁻¹ = 1/diag(A) | O(N),简单但效果一般 |
| Block Jacobi | 每个秩的本地块作为预处理器 | 本地求解,无通信 |
| Incomplete LU (ILU) | 不完全 LU 分解 | 效果更好但通信复杂 |
| Multigrid | 多尺度方法 | 最优 O(N) 但实现复杂 |
# 每个秩 i 负责行 [start_i, end_i]
# 计算 A*p 时:
# 1. 接收左邻居的 p[start_{i-1}](如果存在)
# 2. 接收右邻居的 p[end_{i+1}](如果存在)
# 3. 本地计算 SpMV
# 4. 发送本地 p[start_i] 和 p[end_i] 给邻居
PCG 的主要计算是 SpMV,对于三对角矩阵:
HW4 是课程并行数值线性代数的集大成者。相关主题在后续课程中深入:
卡住时的正确顺序: 复习 Lecture 7-8, 16 → 在纸上画分布式布局 → 用 P=2 调试 → 验证数值正确性 → 再扩展到更多秩 → 分析性能。
调试黄金法则: 先保证正确(串行和分布式结果一致),再追求性能。一个快但错的结果毫无价值。
问题 1:1D Poisson 方程离散后,系数矩阵 A 是什么结构?每行有多少个非零元素?
问题 2:分布式 SpMV 时,为什么需要 halo 交换?
问题 3:PCG 迭代中计算 r·z 时,为什么需要 MPI_Allreduce?
问题 4:最简单的 Jacobi 预处理器是如何工作的?
问题 5:为什么 HW4 要求 N 必须是秩数的倍数?
分布式 MPI PCG 求解器
稀疏矩阵 + SpMV + 集合通信
重用中间值 + 通信隐藏
扩展图 + 性能分析 + 挑战讨论