使用有限差分方法 - 拉普拉斯算子（第三部分）

Finite difference method - Laplacian part 3 — ROCm Blogs (amd.com)

2023年5月11日，由[Justin Chang](Justin Chang — ROCm Blogs)、[Rajat Arora](Rajat Arora — ROCm Blogs)、[Thomas Gibson](Thomas Gibson — ROCm Blogs)、[Sean Miller](Sean Miller — ROCm Blogs)和[Ossian O’Reilly](Ossian O’Reilly — ROCm Blogs)撰写。

在前两篇关于拉普拉斯算子的文章中，我们开发了一个基于HIP的有限差分代码，该代码围绕拉普拉斯算子设计，并应用了两种可能的代码优化，以优化L2缓存和全局内存之间的内存移动。第三部分将涵盖一些额外的优化和通用的微调性能的提示。快速回顾一下，记住拉普拉斯算子采用标量场u(x,y,z)的梯度的散度形式：

,

我们在[第一部分](Finite difference method - Laplacian part 1 — ROCm Blogs)开始的基线HIP实现的性能达到了理论峰值的50%左右。然而，基于初步的`rocprof`分析，我们预计有限差分内核至少应该达到71%的峰值。为了实现这一目标，我们应用了两种优化：

1. 引入循环块划分（loop tiling）来显式重用已加载的模板
2. 重新排序模板点的读取访问模式

请参阅[第二部分](Finite difference method - Laplacian part 2 — ROCm Blogs)中重新排序读取访问模式部分，以获取完整的代码实现。通过这些更改，我们已达到性能预测的95%，但仍有一些未解决的问题：

- 先前的优化需要手动调整某些参数，这可能会导致性能出现突然下降的大规模块划分。能否通过解决此性能下降的根本原因使我们更接近在[第一部分](Finite difference method - Laplacian part 1 — ROCm Blogs)中定义的评价指标（FOM）？
- 我们仅关注优化内核的缓存和读取操作的数据重用。通过改进写操作的同样方面，我们能否取得一些进展？
- 我们引入的优化需要非平凡的代码更改。是否有替代优化可以在不增加代码复杂性的情况下显著提高性能？

这篇博文将探讨其中一些剩余问题。接下来的几节将介绍和讨论以下概念：

1. 生成临时文件以了解寄存器使用和暂存内存
2. 应用启动界限（launch bounds）来控制寄存器使用
3. 应用非时态存储以释放更多缓存 

寄存器压力和溢出使用

在前一篇文章中描述的循环平铺优化中，当平铺因子 m=16 时，内核的性能指标 (FOM) 恶化。`rocprof` 指标 FETCH_SIZE 上升到超过理论限制的4倍，`L2CacheHit` 指标下降到50%以下。我们怀疑高平铺因子增加了寄存器的使用，导致溢出。为此，我们引入了一个新的编译标志 --save-temps，它告诉编译器生成关于寄存器使用、溢出、占用率等的重要信息，这些信息适用于每个 GPU 内核。它还包括主机和设备代码的指令集架构 (ISA) 转储。未来的实验笔记文章将详细介绍 AMD GPU ISA。

我们检查四个关键指标：

1. SGPR

2. VGPR

3. Scratch

4. Occupancy

SGPR 和 VGPR 指的是标量和矢量寄存器，Scratch 指的是可能表示寄存器溢出的临时存储器，Occupancy 代表在执行单元 (EU) 上可以运行的波前的最大数量。请注意，寄存器和溢出的统计数据可以直接从 rocprof 输出文件中找到，而其他细节如占用率和 ISA 汇编只能从临时文件中找到。用户可以简单地取消注释 makefile 中的 #TEMPS=true 行，以生成位于 temps/laplacian_dp_kernel1-hip.s 文件中的临时文件，包含这些信息。以下是基线 HIP 内核1 的一些示例输出：

  .section  .AMDGPU.csdata
; Kernel info:
; codeLenInByte = 520
; NumSgprs: 18
; NumVgprs: 24
; NumAgprs: 0
; TotalNumVgprs: 24
; ScratchSize: 0
; MemoryBound: 0
; FloatMode: 240
; IeeeMode: 1
; LDSByteSize: 0 bytes/workgroup (compile time only)
; SGPRBlocks: 2
; VGPRBlocks: 2
; NumSGPRsForWavesPerEU: 18
; NumVGPRsForWavesPerEU: 24
; AccumOffset: 24
; Occupancy: 8
; WaveLimiterHint : 1
; COMPUTE_PGM_RSRC2:SCRATCH_EN: 0
; COMPUTE_PGM_RSRC2:USER_SGPR: 8
; COMPUTE_PGM_RSRC2:TRAP_HANDLER: 0
; COMPUTE_PGM_RSRC2:TGID_X_EN: 1
; COMPUTE_PGM_RSRC2:TGID_Y_EN: 1
; COMPUTE_PGM_RSRC2:TGID_Z_EN: 1
; COMPUTE_PGM_RSRC2:TIDIG_COMP_CNT: 2
; COMPUTE_PGM_RSRC3_GFX90A:ACCUM_OFFSET: 5
; COMPUTE_PGM_RSRC3_GFX90A:TG_SPLIT: 0

除了 ISA 转储之外，还有很多信息需要解读。我们将所有感兴趣的读者推介给[这个演示文稿]以及这篇[关于寄存器压力的文章](Register pressure in AMD CDNA™2 GPUs — ROCm Blogs)，以了解有关寄存器、溢出、占用率等更多详细信息。

以下表格显示了基线和不同平铺因子 m 下优化内核的上述四个关键指标：

	SGPR	VGPR	Scratch	Occupancy
Kernel 1 - Baseline	18	24	0	8
Kernel 3 - Reordered loads m=1	24	18	0	8
Kernel 3 - Reordered loads m=2	26	28	0	8
Kernel 3 - Reordered loads m=4	34	54	0	8
Kernel 3 - Reordered loads m=8	52	90	0	5
Kernel 3 - Reordered loads m=16	90	128	180	4

寄存器/溢出使用与平铺因子有很强的相关性。当 m=16 时，寄存器使用量增大到“溢出”到临时存储空间，即寄存器在寄存器空间中不再适应，必须卸载到全局内存中。还有一个强逆相关性，定义为 EU 上的波前数（occupancy）与寄存器使用量 - 当寄存器使用量增加时，占用率下降。那么，如何防止溢出或增加占用率呢？ 

启动范围

一种快速控制寄存器使用的方法是应用启动范围（launch bounds）到内核。默认情况下，HIP编译器基于最大允许的线程块大小（1024个线程）限制每个线程的寄存器数量。如果在编译时知道线程块的大小，那么为内核设置启动范围是一个好习惯。设置启动范围采用以下参数：

__launch_bounds__(MAX_THREADS_PER_BLOCK,MIN_WAVES_PER_EU)

第一个参数 MAX_THREADS_PER_BLOCK 告知编译器线程块的维度，以便优化特定块大小的寄存器使用。第二个参数 MIN_WAVES_PER_EU 是一个可选参数，指定每个执行单元（EU）上需要激活的最小波阵面数量。默认情况下，这个第二个值设置为1，并且通常不需要修改，而默认的 MAX_THREADS_PER_BLOCK 值1024需要修改，因为我们不使用所有1024个线程。

到目前为止，我们一直在使用256×1×1的线程块大小，下面是如何为Kernel 3设置启动范围，`MAX_THREADS_PER_BLOCK = 256`：

template <typename T>
__launch_bounds__(256)
__global__ void laplacian_kernel(...) {

...

让我们用这个一行代码的修改定义“Kernel 4”，并检查其对寄存器和暂存空间使用的影响：

	SGPR	VGPR	Scratch	Occupancy
Kernel 1 - Baseline	18	24	0	8
Kernel 3/Kernel 4 - Reordered loads m=1	24/24	18/18	0/0	8/8
Kernel 3/Kernel 4 - Reordered loads m=2	26/26	28/28	0/0	8/8
Kernel 3/Kernel 4 - Reordered loads m=4	34/34	54/54	0/0	8/8
Kernel 3/Kernel 4 - Reordered loads m=8	52/52	90/94	0/0	5/5
Kernel 3/Kernel 4 - Reordered loads m=16	90/84	128/170	180/0	4/2

对于 m=4 及以下的瓦片因子，应用启动范围对寄存器使用没有影响。当瓦片因子为 m=8 时，仅VGPR略有增加，而对于 m=16，VGPR显著增加，暂存空间使用则完全消失。注意到占用率（occupancy）显著下降，这引起了关于性能是否会受到负面影响的问题。让我们来看一下FOM性能：

	Speedup	% of target
Kernel 1 - Baseline	1.00	69.4%
Kernel 3 - Reordered loads m=1	1.20	82.9%
Kernel 3 - Reordered loads m=2	1.28	88.9%
Kernel 3 - Reordered loads m=4	1.34	93.1%
Kernel 3 - Reordered loads m=8	1.37	94.8%
Kernel 3 - Reordered loads m=16	0.42	29.4%
Kernel 4 - Launch bounds m=1	1.20	82.9%
Kernel 4 - Launch bounds m=2	1.28	88.9%
Kernel 4 - Launch bounds m=4	1.34	93.1%
Kernel 4 - Launch bounds m=8	1.39	96.1%
Kernel 4 - Launch bounds m=16	1.34	93.2%

不出所料，对于启动范围没有影响寄存器、暂存位置或占用率的内核，其性能与之前相同。启动范围对SGPR、VGPR、Scratch以及Occupancy统计有影响的内核，其性能明显提升。瓦片因子为 m=8 和 m=16 的内核分别在其性能上有了改进。让我们来看相应的 rocprof 指标：

	FETCH_SIZE (GB)	Fetch efficiency (%)	L2CacheHit (%)
Theoretical	1.074	-	-
Kernel 1 - Baseline	2.014	53.3	65.0
Kernel 3 - Reordered loads m=1	1.347	79.7	72.0
Kernel 3 - Reordered loads m=2	1.166	92.1	70.6
Kernel 3 - Reordered loads m=4	1.107	97.0	68.8
Kernel 3 - Reordered loads m=8	1.080	99.4	67.7
Kernel 3 - Reordered loads m=16	3.915	27.4	44.5
Kernel 4 - Launch bounds m=1	1.346	79.8	72.0
Kernel 4 - Launch bounds m=2	1.167	92.1	70.6
Kernel 4 - Launch bounds m=4	1.107	97.0	68.8
Kernel 4 - Launch bounds m=8	1.080	99.4	67.3
Kernel 4 - Launch bounds m=16	1.094	98.2	66.1

NOTE: 虽然这些实验未展示 `WRITE_SIZE` 和写效率（%），但 `Kernel 3 - Reordered loads m=16` 的报告 `WRITE_SIZE` 和写效率（%）分别为2.547 GB和41.7%。没有暂存溢出的内核几乎在写效率方面达到100%。

现在，`m=16` 的启动范围不再将寄存器溢出到暂存空间中，因此我们看到了显著的性能提升、提取效率（fetch efficiency）和 L2 缓存命中率（L2CacheHit）的提高。使用启动范围可以快速恢复由于高寄存器使用导致的性能损失。虽然 m=8 具有比 m=16 更低的寄存器使用，但应用启动范围仍然对 VGPR 使用产生了影响，从而使总性能略有增加，使其成为目前表现最好的内核。新的 m=8 内核的 FOM 仍略低于预期目标，因此我们需要探索其他代码优化方案。

非临时存储器访问

我们的大多数优化工作都集中在提高空间局部性上，但我们尚未考虑的是时间局部性——即如何优先缓存变量在时间上的使用。加载各`u`元素以及存储各`f`元素都会占用缓存行。然而根据[第一部分](Finite difference method - Laplacian part 1 — ROCm Blogs)描述的数据布局，每个`u`元素理论上可能会最多被重用六次，而每个`f`元素仅被访问一次。因此，我们可以使用clang的内建非临时存储指令，让`f`跳过L2缓存，从而增加`u`条目的可用缓存。需要注意的是，这些内建函数是特定于AMD GPU的。

AMD clang编译器提供了两种重载的内建函数，允许生成非临时加载和存储：

T __builtin_nontemporal_load(T *addr);
void __builtin_nontemporal_store(T value, T *addr);

在Laplacian示例中，我们只需要非临时存储。首先将此内建函数应用于初始的基线内核：

Kernel 1 (Before)	Kernel 1 (After)
f[pos] = u[pos] * invhxyz2 + (u[pos - 1] + u[pos + 1]) * invhx2 + (u[pos - nx] + u[pos + nx]) * invhy2 + (u[pos - slice] + u[pos + slice]) * invhz2;	__builtin_nontemporal_store(u[pos] * invhxyz2 + (u[pos - 1] + u[pos + 1]) * invhx2 + (u[pos - nx] + u[pos + nx]) * invhy2 + (u[pos - slice] + u[pos + slice]) * invhz2, &f[pos]);

为了评估这一简单代码修改的影响，我们比较了基线实现（带有和不带有非临时存储）的性能，以及`m=1`时的Kernel 3性能：

	Speedup	% of target
Kernel 1 - Baseline	1.00	69.4%
Kernel 1 - Nontemporal store	1.19	82.5%
Kernel 3 - Reordered loads m=1	1.20	82.9%

对单行代码的改动与重新构建整个基线内核以利用循环分块和重新排序内存访问模式具有相似的改进。以下是`rocprof`统计数据：

	FETCH_SIZE (GB)	Fetch efficiency (%)	L2CacheHit (%)
Theoretical	1.074	-	-
Kernel 1 - Baseline	2.014	53.3	65.0
Kernel 1 - Nontemporal store	1.429	75.2	71.4
Kernel 3 - Reordered loads m=1	1.347	79.7	72.0

统计数据显示，当`m=1`时，非临时存储和重新排序加载的性能相当。这些结果表明，利用非临时内建函数进行内存访问实际上可能是用户首先应采取的优化措施，因为它只需要修改一行代码，这是一个“较低的果实”。显然，下一个问题是在结合使用非临时存储、循环分块因子`m=8`以及启动边界时会发生什么？让我们再次对Kernel 4进行一行修改：

Kernel 4 (Before)	Kernel 5 (After)
f[pos + n*nx] = Lu[n];	__builtin_nontemporal_store(Lu[n],&f[pos + n*nx]);

这个新的Kernel 5是结合了循环分块、重新排序加载、应用启动边界和利用非临时存储的所有优化。当循环分块因子`m=8`时的性能如下：

	Speedup	% of target
Kernel 1 - Baseline	1.00	69.4%
Kernel 3 - Reordered loads m=8	1.37	94.8%
Kernel 4 - Launch bounds m=8	1.39	96.1%
Kernel 5 - Nontemporal store m=8	1.44	100%

通过结合所有这些优化，我们实现了1.44倍的加速，并达到了100%的目标！让我们再来看看`rocprof`的度量数据：

	FETCH_SIZE (GB)	Fetch efficiency (%)	L2CacheHit (%)
Theoretical	1.074	-	-
Kernel 1 - Baseline	2.014	53.3	65.0
Kernel 3 - Reordered loads m=8	1.080	99.4	67.7
Kernel 4 - Launch bounds m=8	1.080	99.4	67.3
Kernel 5 - Nontemporal store m=8	1.074	100	67.4

在这些改进过程中，我们将从全局内存的加载次数减少了一半。测得的提取和写入大小已经达到理论极限，因此进一步的性能改进必须来自其他方向。由于报告的有效内存带宽已经达到预期目标，可能没有多少改进空间了。