diff --git "a/guide/ch4.Kernel\346\200\247\350\203\275\345\273\272\346\250\241.sgemm.md" "b/guide/ch4.Kernel\346\200\247\350\203\275\345\273\272\346\250\241.sgemm.md"
index 78e5612e07317a27a96d15976e60314948c90a48..f48216233fa5b3499858b5166d883cc052af9504 100644
--- "a/guide/ch4.Kernel\346\200\247\350\203\275\345\273\272\346\250\241.sgemm.md"
+++ "b/guide/ch4.Kernel\346\200\247\350\203\275\345\273\272\346\250\241.sgemm.md"
@@ -1,25 +1,31 @@
 # ch4.精度矩阵乘法的性能建模与优化
 
-矩阵乘法堪称 GPU 最为核心的应用领域。从早期图形渲染中顶点渲染器处理的四元数线性变换，到近来深度学习中的全连接层以及注意力机制，无一不是矩阵乘法的用武之地。可以毫不夸张地说，应用对于矩阵乘法的性能需求推动了 GPU 芯片架构的不断演进，而一代代 GPU 芯片持续提升的矩阵乘法性能又催生出人工智能应用的大爆发。
+矩阵乘法堪称GPU最为核心的应用领域。从早期图形渲染中顶点渲染器处理的四元数线性变换，到近来深度学习中的全连接层以及注意力机制，无一不是矩阵乘法的用武之地。可以毫不夸张地说，应用对于矩阵乘法的性能需求推动了GPU芯片架构的不断演进，而一代代GPU芯片持续提升的矩阵乘法性能又催生出人工智能应用的大爆发。
 
-由此可见，理解 GPU 的关键之处就在于洞悉 GPU 是如何高效执行矩阵乘法的。在 GPU 和 AI 公司的技术岗位面试中，矩阵乘法常常作为应用场景，用以考察候选人对 GPU 的理解程度。
+由此可见，理解GPU的关键之处就在于洞悉GPU是如何高效执行矩阵乘法的。在GPU和 AI 公司的技术岗位面试中，矩阵乘法常常作为应用场景，用以考察候选人对GPU的理解程度。
 
-在 GPU 领域，性能的重要性不言而喻，可谓 “性能就是金钱”。客户在购买 GPU 时，不仅会关注硬件所提供的峰值性能，还会留意实测矩阵乘法的实测性能利用率，即实测矩阵乘法性能与峰值性能的百分比。这个值永远小于 100%。AI 客户还会额外关注给定 AI 模型的 MFU（Model Flop Utilization），而 MFU 也与矩阵乘法的实测利用率呈正相关关系。
+而从GPU销售的视角，性能的重要性不言而喻，可谓 “性能就是金钱”。客户在购买GPU时，不仅会关注硬件所提供的峰值性能，还会留意实测矩阵乘法的实测性能利用率，即实测矩阵乘法性能与峰值性能的百分比。这个值永远小于 100%。AI 客户还会额外关注给定 AI 模型的 MFU（Model Flop Utilization），而 MFU 也与矩阵乘法的实测利用率呈正相关关系。
+
+并行编程是一门复杂且容易出错的技术，常被比作 “火箭科学”。CUDA 之所以能够成为最受欢迎的高性能并行编程语言之一，其核心优势在于它倡导一种渐进式的编程策略。这种策略允许程序员首先以串行编程的思维模式，聚焦于单个线程和单一地址空间，构建功能原型。随后，他们可以充分利用GPU的独特特性进行性能上的优化和提升。
+
+在本节中，我们将采用这种渐进式的编程策略来完成单精度矩阵乘法在GPU上的性能优化。首先，基于单一地址空间上的单线程编程抽象，实现一个基准版本。之后，再利用GPU独有的编程特性，并结合量化分析方法，有的放矢地逐步进行优化。
 
 然而，需要注意的是，当程序被优化得越快，代码的可读性往往会变差，同时可维护性、可扩展性、安全性等方面也需要相应地付出代价。性能调优工程需要在努力接近硬件峰值性能的同时，在其他软件工程方面进行权衡。
 
-并行编程是一门复杂且容易出错的技术，常被比作 “火箭科学”。CUDA 之所以能够成为最受欢迎的高性能并行编程语言之一，其核心优势在于它倡导一种渐进式的编程策略。这种策略允许程序员首先以串行编程的思维模式，聚焦于单个线程和单一地址空间，构建功能原型。随后，他们可以充分利用 GPU 的独特特性进行性能上的优化和提升。
+![本教程渐进式地在Nvidia A100和MetaX C500上提高性能](figure/ch4/step-by-step.png)
+<!-- ![性能上界模型示例](figure/ch4/bounds.png) -->
 
-在本节中，我们将采用这种渐进式的编程策略来完成单精度矩阵乘法在 GPU 上的性能优化。首先，基于单一地址空间上的单线程编程抽象，实现一个基准版本。之后，再利用 GPU 独有的编程特性，并结合量化分析方法，有的放矢地逐步进行优化。
+GPU矩阵乘法的性能优化是一个渐进式的过程，如上图所示。每一步优化都能带来明显的性能提升，这种"小胜利"形成了强大的正反馈循环，激励开发者继续深入学习和探索。从最初的朴素实现到最终接近cuBLAS的高性能版本，每一步都代表着对GPU架构理解的加深。这种渐进式方法不仅降低了学习门槛，也让复杂的性能优化变得可行和有趣。
 
+值得注意的是，矩阵乘法的性能表现高度依赖于矩阵维度（M、N、K）。在小矩阵（如1024×1024）情况下，启动开销占比较大，而在大矩阵（如4096×4096）情况下，内存带宽和计算效率成为主要瓶颈。不同维度组合下的最佳优化策略也各不相同：方形矩阵（M=N=K）适合均衡的分块策略，而细长矩阵（如M>>N）则需要特殊处理以维持良好的内存访问模式和负载均衡。
 
-希望有意掌握 GPU 的读者能够在实际的 GPU 上，跟随本教程，一步步亲自尝试编写矩阵乘法的性能优化代码。每一步优化所带来的性能提升，都代表着读者对 GPU 的理解更进一步。日拱一卒，每一小步都是一次小小的正反馈，让读者更有动力深入学习愈加艰深的教程。
+GPU厂商随芯片推出的cuBLAS, mcBLAS等矩阵运算库，之所以能在各种矩阵维度下都表现出色，正是因为它采用了复杂的算法选择机制。它不是单一的实现，而是包含数十种针对不同场景优化的kernel实现，通过启发式规则或查表方式在运行时选择最适合当前问题规模的算法。例如，对于小矩阵可能选择减少启动开销的实现，而对于大矩阵则选择最大化计算访存比的实现。真正高效的矩阵乘法库需要针对不同的M、N、K组合采用不同的算法。理解这种"一个尺寸不适合所有情况"的原则，是迈向高性能GPU编程的关键一步。
 
+本节介绍的方法，不局限在特定一款GPU上。希望有意掌握GPU的读者能够在手边的GPU上，跟随本教程，一步步亲自尝试编写矩阵乘法的性能优化代码。每一步优化所带来的性能提升，都代表着读者对GPU的理解更进一步。日拱一卒，每一小步都是一次小小的正反馈，让读者更有动力深入学习愈加艰深的教程。
 本节中所有代码及可运行的示例见[case/4.sgemm_tt](../case/4.sgemm_tt)
 
-![性能上界模型示例](figure/ch4/bounds.png)
+# 1.2 矩阵乘法的问题描述
 
-# 1.2 1.2.	矩阵乘法的问题描述
 在现代深度学习框架中，高达95%以上的计算时间都花在了矩阵乘法上。这个看似简单的数学运算，背后蕴含着丰富的计算机科学智慧。
 
 ## 1.2.1 问题的数学描述
@@ -28,15 +34,37 @@
 
 $$C_{ij} = \sum_{k=0}^{K-1} A_{ik} \times B_{kj}$$
 
-看起来很简单。但是当矩阵维度达到数千甚至数万时，计算量就会急剧增加。以一个常见的场景为例：
+对于矩阵乘法C = A × B，涉及三个关键维度：
+
+- A矩阵：M × K
+- B矩阵：K × N
+- C矩阵：M × N
+
+我们通常将这样的运算简记为"M × N × K矩阵乘法"。这里有一个重要的观察：K维度在输出矩阵C中消失了，这是因为K方向上进行了归约（Reduction）操作。归约维度和非规约维度的区分，在深度学习自动调优框架、分布式框架中，有重要的应用。
+
+![M x N x K的矩阵乘法](figure/ch4/naive.png)
+
+矩阵乘法C = A × B的计算量与访存量与矩阵维度M、N、K密切相关。对于M×K的矩阵A与K×N的矩阵B相乘得到M×N的矩阵C，计算量约2×M×N×K次浮点运算。
+
+读访存量方面，需要从内存读取A矩阵的M×K个元素和B矩阵的K×N个元素，总计(M×K+K×N)个元素，若每个元素为4字节的float类型，则读访存量为4×(M×K+K×N)字节。
+
+写访存量则是将计算结果C矩阵的M×N个元素写回内存，为4×M×N字节。
+
+随着矩阵维度增大，计算量以三次方增长(O(n³))，而访存量仅以二次方增长(O(n²))，这种计算访存比的提升使得大规模矩阵乘法更容易达到计算密集型状态，从而更充分地利用GPU的计算能力。
+
+当矩阵维度达到数千甚至数万时，计算量就会急剧增加。以一个常见的场景为例：
 
-假设我们要进行一次批量图像处理，输入是100张224×224的图片，要通过一个具有1024个神经元的全连接层。这个操作可以表示为一次矩阵乘法：(100×50176) × (50176×1024)。简单计算一下，这涉及到超过50亿次的乘加运算！
+假设我们要进行一次批量图像处理，输入是100张224×224的图片，要通过一个具有1024个神经元的全连接层。这个操作可以表示为一次矩阵乘法：(100 x 1024 x (224 x 224))。简单计算一下，这涉及到超过50亿次的乘加运算！
+
+TFLOPS（每秒万亿次浮点运算）定义了芯片每秒可以执行的浮点运算数量，数值越高越好。需要注意的是，不同数据类型（如 int8、bf16、tf32）的 TFLOPS 定义不同。通常，int8 的速度是 bf16 的 2 倍，而 bf16 又是 tf32 的 2 倍。用总计算量，除以芯片对应精度的实际TFLOPS，就能得到矩阵乘法计算时间。
+
+实测的 TFLOPS 取决于矩阵大小。即使加速器 A 的理论 TFLOPS 高于加速器 B，也不意味着 A 的实际性能一定更快。理论值在实践中永远无法完全达到，实际的 TFLOPS 效率（HFU，Hardware Floating-point Utilization）在不同厂商之间，甚至同一厂商的不同加速器架构之间可能存在显著差异。
 
 ## 1.2.2 BLAS：基础线性代数子程序库与SGEMM详解
 
 面对如此庞大的计算量，我们显然需要高度优化的实现。这就要提到BLAS（Basic Linear Algebra Subprograms，基础线性代数子程序库）。BLAS诞生于1979年，是计算机科学史上最成功的软件接口标准之一。
 
-让我们首先深入理解BLAS中最重要的操作之一：SGEMM。这个名字是"Single-precision GEneral Matrix Multiplication"的缩写，让我们逐字解析它的含义：
+让我们首先深入理解BLAS中最基本的算子之一：SGEMM。这个名字是"Single-precision GEneral Matrix Multiplication"的缩写，让我们逐字解析它的含义：
 
 1. **Single-precision (S)**：表示使用32位浮点数（float类型）。BLAS支持四种基本数据类型：
    - S：单精度浮点数（float）
@@ -44,10 +72,11 @@ $$C_{ij} = \sum_{k=0}^{K-1} A_{ik} \times B_{kj}$$
    - C：单精度复数（std::complex<float>）
    - Z：双精度复数（std::complex<double>）
 
-值得注意的是，随着深度学习的发展，特别是混合精度训练的兴起，新的数据类型不断涌现。现代GPU还支持：
-- HGEMM：半精度（16位）矩阵乘法
-- FP8GEMM：8位浮点数矩阵乘法
-等更多变体。
+    值得注意的是，随着深度学习的发展，特别是混合精度训练的兴起，新的数据类型不断涌现。现代GPU还支持：
+    - HGEMM：半精度（16位）矩阵乘法
+    - FP8GEMM：8位浮点数矩阵乘法
+
+    等更多变体。
 
 2. **GEneral (GE)**：表示处理的是一般性矩阵，即密集二维数组，不假定具有任何特殊结构。BLAS支持多种特殊矩阵类型，常用的一些矩阵类型还有：
    - SY：对称矩阵
@@ -62,18 +91,9 @@ $$C_{ij} = \sum_{k=0}^{K-1} A_{ik} \times B_{kj}$$
 
 在深入接口设计之前，我们需要理解两个关键概念：矩阵维度和内存布局。
 
-对于矩阵乘法C = A × B，涉及三个关键维度：
-- A矩阵：M × K
-- B矩阵：K × N
-- C矩阵：M × N
-
-我们通常将这样的运算简记为"M × N × K矩阵乘法"。这里有一个重要的观察：K维度在输出矩阵C中消失了，这是因为K方向上进行了归约（Reduction）操作。归约维度和非规约维度的区分，在深度学习自动调优框架
-
-![M x N x K的矩阵乘法](figure/ch4/naive.png)
-
-
 **内存布局：二维数据的一维存储**
 虽然矩阵是二维的，但计算机内存是一维的，这就引出了内存布局的问题。假设矩阵的大小为M x N，主要有两种布局方式：
+
 1. **行主序（Row-Major）**：
    - C语言/NumPy/PyTorch的默认顺序
    - 元素访问：C[i][j] = i*N + j
@@ -84,6 +104,8 @@ $$C_{ij} = \sum_{k=0}^{K-1} A_{ik} \times B_{kj}$$
    - 元素访问：C[i][j] = j*M + i
    - 在BLAS中记为N（Not Transposed）
 
+![一维内存中用行主序和列主序的方式存同一个矩阵](figure/ch4/major.png)
+
 值得思考的是：如果用列主序的方式读取一个以行主序存储的矩阵，实际上相当于读取了原矩阵的转置！这个特性在BLAS接口设计中得到了巧妙的利用。
 
 ## 1.2.4 BLAS接口的设计
@@ -125,6 +147,7 @@ void SGEMM_TT_GPU(cublasHandle_t handle, int m, int n, int k,
 ```
 
 这段代码有一些值得注意的细节：
+
 1. 所有操作都使用CUBLAS_OP_N（列主序），但实际数据是行主序
 2. A和B的调用顺序被交换了
 3. 维度参数的顺序也做了调整
@@ -147,13 +170,12 @@ void SGEMM_TT_GPU(cublasHandle_t handle, int m, int n, int k,
 
 1. **性能敏感**：在很多应用中，GEMM是性能的关键瓶颈。深入理解其实现原理，可以帮助我们做出更好的架构决策。
 
-2. **定制化需求**：标准BLAS库针对的是通用场景，在特定应用中（如稀疏矩阵、低精度计算、特殊硬件等），往往需要定制化的实现。
+2. **定制化需求**：标准BLAS库针对的是通用场景，在特定应用中（如稀疏矩阵、融合算子等），往往需要定制化的实现。
 
-3. **算法创新**：理解GEMM的优化思路对开发新的算法很有帮助。例如，Winograd算法就是通过减少乘法次数来加速卷积运算，其核心思想与矩阵乘法优化类似。
+3. **算法创新**：理解GEMM的优化思路对开发新的算法很有帮助。例如，FlashAttention算法就是通过减少数据搬移次数加速Attention算子在GPU上的性能，其核心思想就来源于对GPU上矩阵乘法优化的深入洞察。
 
 在接下来的章节中，我们将从最简单的三重循环实现开始，逐步引入各种优化技术，最终实现一个高性能的GPU GEMM。这个过程不仅能帮助你掌握性能优化的技巧，更能培养系统性的算法设计思维和硬件性能建模思维。
 
-
 # 1.3 简单矩阵乘法：三重循环
 
 在开始深入GPU优化之前,让我们先从最基础的实现开始。就像学习任何复杂的知识一样,理解基础版本的不足之处,才能更好地理解为什么需要后续的优化。我们就从最简单的三重循环矩阵乘法开始,看看它在GPU上是如何实现的,以及存在哪些问题。
@@ -218,12 +240,14 @@ __global__ void sgemm_naive(int M, int N, int K,
 ### 线程索引计算
 
 首先看这两行代码:
+
 ```cpp
 int row = blockIdx.x * blockDim.x + threadIdx.x;
 int col = blockIdx.y * blockDim.y + threadIdx.y;
 ```
 
 这是CUDA编程中最基本的线程定位公式。想象一个大礼堂中的座位安排：
+
 - `blockIdx`就像区号(比如第3区)
 - `blockDim`是每个区的行/列数(比如每区16排16列)
 - `threadIdx`则是区内的具体位置(比如第5排第6列)
@@ -296,19 +320,18 @@ B[k * N + col]
 
 在深入内存聚合之前,我们需要深入理解线程束（Warp）。如果说GPU是一支军队,那么Warp就是其中的一个班组。它具有以下特点:
 
-- 每个Warp固定包含32或64个线程,这些线程步调一致地执行指令。Warp的大小随GPU的不同而不同，可以在CPU上用以下代码查询。
-```
-    cudaDeviceProp prop;
-    cudaGetDeviceProperties(&prop, 0);  // 查询设备 0。多卡时需要更改
-    printf("Warp size: %d\n", prop.warpSize);
-```
+- 每个Warp固定包含32或64个线程,这些线程步调一致地执行指令。Warp的大小随GPU的不同而不同，可以在CPU上用以下代码查询。通常Nvidia的GPU上每个warpSize = 32，而沐曦的C500上warpSize = 64。
+    ```
+        cudaDeviceProp prop;
+        cudaGetDeviceProperties(&prop, 0);  // 查询设备 0。多卡时需要更改
+        printf("Warp size: %d\n", prop.warpSize);
+    ```
 - 每个SM(Stream Multiprocessor)配备了4个warp调度器,就像4个班长在管理自己的班组
 - 线程是如何被分到各个"班组"(Warp)的呢？这要从threadId说起
 - 对于多维的blockDim,threadId的计算遵循以下公式:
-
-```cpp
-threadId = threadIdx.x + blockDim.x * (threadIdx.y + blockDim.y * threadIdx.z)
-```
+    ```cpp
+    threadId = threadIdx.x + blockDim.x * (threadIdx.y + blockDim.y * threadIdx.z)
+    ```
 
 相邻threadId的线程会被分配到同一个Warp中,这个特性对我们后面的优化至关重要。
 
@@ -334,6 +357,8 @@ threadId = threadIdx.x + blockDim.x * (threadIdx.y + blockDim.y * threadIdx.z)
 
 这种访问模式就像是32个人要买东西,但他们要去的商店隔得很远。即使每个人只买一个东西,也无法合并成一次"团购"。
 
+为什么会有128字节的合并访问跨度的约束呢？熟悉硬件的同学可能已经看出来，这是GPU上Cacheline（缓存行）的大小。GPU缓存行大小设定为128B是片上存储与片外存储之间平衡的结果。DRAM作为片外存储，就像一个容量巨大但访问较慢的仓库；而SRAM作为片上缓存，则像是一个小型高效的便利店，其中128B的缓存行大小代表了一次合理的"购物量"。SRAM虽然访问延迟仅有1-2ns，速度远超DRAM，但其面积大、成本高，导致容量有限，必须精打细算地使用。128B的缓存行大小充分考虑了程序执行的时间局部性和空间局部性特征——近期访问过的数据很可能再次被访问，而访问某个数据时其附近数据也很可能被用到。缓存行过小会导致地址和控制电路的额外开销过大，浪费宝贵的芯片面积；而过大则会因预取过多无用数据而降低缓存效率。128B恰好能容纳32个单精度浮点数，既能充分利用空间局部性预取临近数据，又能与GPU的warp大小（32线程）完美匹配，使得一个warp的线程在访问连续内存时能够通过一次内存事务完成，实现最佳的内存访问效率。
+
 ## 1.4.3 优化实战 - 重新设计线程映射
 
 明白了原理,让我们动手优化代码。关键在于重新设计线程与矩阵元素的对应关系:
@@ -374,7 +399,6 @@ int x = blockIdx.y * blockDim.y + threadIdx.y;
 
 实际测试中，这种方法的性能提升与之前的方法一致。
 
-
 # 1.5 分块矩阵乘法与屋檐模型
 
 上一节我们通过内存聚合让GPU的访存更高效了,性能也确实提升了不少。但是,当你满怀期待地运行benchmark,看到的性能可能还是不尽如人意。距离cuBLAS的性能差了何止一个数量级！
@@ -394,6 +418,7 @@ for(int k = 0; k < K; k++) {
 这个比值究竟意味着什么？这里我们要引入性能分析领域的最重要的工具之一：**屋檐模型(Roofline Model)**。它就像是一把"尺子",能帮我们量出程序的性能上限在哪里。
 
 ### 1.5.1 屋檐模型:性能的"天花板"
+
 屋檐模型（Roofline Model）是一个用于分析和预测程序性能上限的可视化工具。它最早由 UC Berkeley 的 Samuel Williams 等人在2009年提出。这个模型之所以叫"屋檐"，是因为它的图形看起来就像房子的侧面屋檐 —— 一个斜坡加一个水平线。
 
 ![alt text](figure/ch4/roofline.png)
@@ -416,13 +441,13 @@ for(int k = 0; k < K; k++) {
 
 根据程序在屋檐图上的位置，我们可以将性能状态分为三种典型情况：带宽受限、算力受限和性能不足。让我们详细探讨这三种状态及其优化思路。
 
-![alt text](figure/ch4/region.png)
+![算力受限与带宽受限](figure/ch4/region.png)
 
 当程序的性能点落在斜线区域附近时，我们称之为带宽受限状态（Memory-Bound）。这种情况下，程序的计算访存比小于硬件的算力带宽比，表明内存带宽成为了性能提升的瓶颈。在带宽受限的情况下，性能会随着内存带宽的提升而线性提升。针对这种情况，优化的重点应该放在改善内存访问模式上。我们可以采用缓存友好的数据结构，调整数据布局以提高空间局部性，或者合并访存操作以减少内存事务。同时，减少不必要的数据移动也很重要，这包括消除冗余的内存读写、充分利用寄存器复用数据，以及在可能的情况下压缩数据格式。
 
 当程序的性能点接近水平线时，我们称之为算力受限（Compute-Bound）状态。此时程序的计算访存比大于硬件的算力带宽比，性能主要受限于处理器的计算能力。对于GPU应用而言，达到算力受限通常是一个好消息，意味着优化工作已经相当成功。在这种情况下，如果还需要进一步提升性能，就需要考虑从算法层面进行改进，或者使用更高算力的硬件。
 
-然而，在实际工作中，我们最常遇到的是第三种情况：性能不足（Poor Performance）。这种情况下，程序的性能点远离屋檐模型的边界（无论是斜线还是水平线），与理论预期之间存在明显的差距。造成这种情况的原因可能有很多：例如，指令级并行度不足、线程负载不均衡、同步开销过大、缓存命中率低，或者其他系统级开销。换句话说，系统中存在当前模型中尚未弄清楚的其他上界。
+然而，在实际工作中，我们最常遇到的是第三种情况：性能差。这种情况下，程序的性能点远离屋檐模型的边界（无论是斜线还是水平线），与理论预期之间存在明显的差距。造成这种情况的原因可能有很多：例如，指令级并行度不足、线程负载不均衡、同步开销过大、缓存命中率低，或者其他系统级开销。换句话说，系统中存在当前模型中尚未弄清楚的其他上界。
 
 例如，在最初的三重循环矩阵乘法，虽然理论分析认为它是带宽受限的状态，但实测性能距离表示带宽受限的斜线，依然有很大的距离。在上一节我们经过分析，了解到三重循环矩阵乘法，每个Warp内相邻线程的访问请求无法聚合成一笔128B的内存请求。每个线程的访问都会产生独立的内存事务。因此，实际的有效访存带宽，只有峰值带宽的1/32(因为每个线程访问4-Byte float，理想情况下每32个线程能聚合为一笔128B的内存请求；未聚合是内存请求的数量，是理想聚合的32倍；于是未聚合的有效带宽是理想峰值带宽的1/32)。这一点可以通过在之前的屋檐模型下，新增一个未聚合请求有效带宽的斜线，从而观察到性能受限与该带宽。
 
@@ -565,7 +590,7 @@ __global__ void sgemm_tiling(int M, int N, int K, float *A, float *B, float *C)
 
 ![alt text](figure/ch4/c500_memory_subsystem.png)
 
-现代GPU的存储系统就像一座金字塔：从延迟最低、带宽最大但容量最小的寄存器，到相对较慢但容量更大的L1缓存、共享内存，再到容量最大但延迟最大、带宽最小的主存（DRAM）。每一层都有其独特的带宽特征，而这些特征都可能成为程序性能的潜在瓶颈。分层屋檐模型正是为了刻画这种复杂的层次化结构而诞生的。
+现代GPU的存储系统就像一座金字塔：从延迟最低、带宽最大但容量最小的寄存器，到相对较慢但容量更大的L1缓存、共享内存，再到容量最大但延迟最大、带宽最小的主存（DRAM）。每一层都有其独特的带宽特征，而这些特征都可能成为程序性能的潜在瓶颈。分层屋檐模型正是为了刻画这种复杂的层次化结构而诞生的。芯片外的存储，使用DRAM实现，容量密度高但速率慢；芯片上的存储，使用SRAM实现，容量密度低但速度快。
 
 ![alt text](figure/ch4/hierarchical_roofline.png)
 
@@ -647,8 +672,6 @@ for (int k = 0; k < BK; k++) {
 } 
 ```
 
-![这张图需要重画](figure/ch4/outer_product.png)
-
 换一个视角可以更直观地感受到外积的优势。如果采用内积，A矩阵的每一行都需要跟TN个B矩阵的列做内积，而B矩阵的每一列都需要跟TM个A矩阵的行做内积。而采用外积，A矩阵的第k列仅仅需要跟B矩阵的第k行做一次外积，之后这一行和这一列都不需要被使用。
 
 理论分析表明，外积相比内积，可以减少线程对共享内存的访问次数，从而提高计算访存比。
@@ -683,7 +706,6 @@ r_c[1][0] += a10 * reg_b;
 r_c[1][1] += a10 * b01;
 ```
 
-
 ## 启示与思考
 
 这个有趣的例子告诉我们：在性能优化的道路上，理论分析固然重要，但也要充分考虑实际执行环境的特点。现代硬件和编译器的能力往往超出我们的想象，有时候"过度优化"反而可能适得其反。
@@ -692,7 +714,7 @@ r_c[1][1] += a10 * b01;
 
 # 1.8 深入理解GPU共享内存冲突与优化方案
 
-在我们的矩阵乘法优化之旅中，共享内存的访问效率是一个关键问题。本节我们将通过一个具体的例子，深入分析共享内存访问中的bank conflict问题，并学习如何解决它。
+在我们的矩阵乘法优化之旅中，共享内存的访问效率是一个关键问题。本节我们将通过一个具体的例子，深入分析共享内存访问中的bank 冲突问题，并学习如何解决它。
 
 ## 理解GPU共享内存的组织方式
 
@@ -723,7 +745,7 @@ for (int m = 0; m < TM; m++) {
 } 
 ```
 
-## 我们来分析哪里会产生Bank Conflict
+## 我们来分析哪里会产生Bank冲突
 
 ### 1. 访问s_a数组时的情况
 
@@ -736,7 +758,7 @@ for (int m = 0; m < TM; m++) {
   - 线程2计算：2 * 8 + 0 = 16
   - 线程3计算：3 * 8 + 0 = 24
   
-当k固定时，这些线程会访问s_a[0][k]、s_a[8][k]、s_a[16][k]、s_a[24][k]。由于共享内存是按照行优先存储的，这些访问会落在同一个bank上！这就造成了严重的bank conflict。
+当k固定时，这些线程会访问s_a[0][k]、s_a[8][k]、s_a[16][k]、s_a[24][k]。由于共享内存是按照行优先存储的，这些访问会落在同一个bank上！这就造成了严重的bank冲突。
 
 ### 2. 访问s_b数组时的情况
 
@@ -745,7 +767,7 @@ for (int m = 0; m < TM; m++) {
 - tx = threadIdx.y，同一warp中tx相同
 - comp_b_smem_n = tx * 8 + n
 - 当k固定时，同一warp的线程访问的是连续的内存位置
-- 这种连续访问模式不会造成bank conflict，因为连续的地址会被映射到不同的bank
+- 这种连续访问模式不会造成bank冲突，因为连续的地址会被映射到不同的bank
 
 ## 如何解决这个问题？
 
@@ -758,8 +780,9 @@ __shared__ float s_a[BM][BK + 1];
 ```
 
 增加一列不使用的padding数据后，原本会发生冲突的访问就被错开了。比如：
-- 原本：线程0、1、2、3分别访问第0、8、16、24个元素（都在bank 0）
+- 原本：线程0、1、2、3分别访问第0、8、16、24个元素，线程间访存间隔为8（都在bank 0）
 - 使用padding后：这些元素会落在不同的bank上，因为每行多了一个元素，改变了映射关系
+- 神奇地是，使用padding后，线程之间访存间隔BK为任意2的整数次幂时，都可以避免bank冲突。这是因为2的整数次幂，与bank数量32，它们的最小公约数为1。
 
 ### 方案二：使用Swizzle重排访问
 
@@ -770,31 +793,65 @@ int swizzled_k = k ^ threadIdx.x;
 r_c[m][n] += s_a[comp_a_smem_m][swizzled_k] * s_b[swizzled_k][comp_b_smem_n];
 ```
 
-这种方法通过异或运算，巧妙地将原本会冲突的访问分散到不同的bank上。这就像是在一个图书馆里，原本所有学生都按照书架的物理顺序（比如从左到右）找书，导致大家都挤在同一个书架前。我们通过重新排列书的位置（但不改变书架的数量），让原本会在同一书架前拥挤的学生自然地分散到不同的书架前。虽然每本书的新位置需要通过一个简单的公式来计算，但一旦熟悉这个规则，就能高效地找到需要的书，而且不会与其他同学产生拥堵。
+这种方法通过异或运算，巧妙地将原本会冲突的访问分散到不同的bank上。
+
+Swizzle初次接触时，很难理解。我们用一个理工科高校图书馆的例子，详细解释一下。
+
+图书馆中，书籍通常按照学科分类系统排列：
+TP类（计算机科学）的书籍全部放在TP书架
+I类（文学）的书籍全部放在I书架
+以此类推...
+
+假设计算机专业的学生很多，期末复习时大量学生需要查阅TP类书籍：
+
+- 学生A需要《操作系统》(TP316)
+- 学生B需要《计算机网络》(TP393)
+- 学生C需要《人工智能》(TP18)
+
+TP类（计算机科学）的书籍全部放在40号书架。
+结果所有学生都挤在40号书架前，造成严重拥堵，而其他书架却几乎无人问津。这就像GPU中的bank冲突，多个线程访问同一内存bank时必须串行化处理。
+
+图书馆管理员设计了一个新系统：书籍的实际存放位置由"科目编码"与"书籍编号的数字部分"进行异或运算决定：
+
+现在，同样的书籍分配情况：
+
+- 《操作系统》(TP316) → 存放位置：40 ^ (316 % 64) = 40 ^ 60 = 20号书架
+- 《计算机网络》(TP393) → 存放位置：40 ^ (393 % 64) = 40 ^ 9 = 33号书架
+- 《人工智能》(TP18) → 存放位置：40 ^ (18 % 64) = 40 ^ 18 = 58号书架
+
+结果原本会全部集中在40号书架的TP类书籍现在分散到了20、33、58等不同书架上，即使多名学生同时查找计算机类书籍，也不会在同一书架前拥挤。
+
+科目编号（如TP对应的40）可以看作是矩阵的行索引，而书籍的数字编号部分（如316、393、18）则对应于矩阵的列索引。传统的图书馆分类方式，就像是按行优先顺序访问矩阵元素，当多个读者同时访问同一"行"（同一科目）的不同"列"（不同书籍）时，就会造成拥堵。
+
+在GPU编程中，这种情况对应着多个线程访问共享内存中同一bank的不同地址，导致bank冲突。通过Swizzle技术（科目编号与书籍数字编号异或），我们实际上是在重新映射矩阵元素的物理存储位置，使得原本在同一"行"上的元素被分散到不同的"行"，从而使得并行访问时能够分散到不同的内存bank上。
+
+这种重新映射不改变数据的逻辑组织，只改变其物理存储位置，从而优化并行访问模式，减少冲突，提高吞吐量。这正是Swizzle技术在高性能计算中的精妙之处——通过简单的位操作，巧妙地解决了复杂的资源竞争问题。
 
 ## 如何选择优化方案？
 
 两种方案各有优势：
+
 - Padding方案实现简单，容易理解，但会多占用一些内存
 - Swizzle方案不需要额外内存，但实现较复杂，需要仔细处理访问模式
 
 在实际应用中，建议：
-1. 先用性能分析工具确认是否存在严重的bank conflict
+
+1. 先用性能分析工具确认是否存在严重的bank冲突
 2. 如果共享内存不紧张，可以先尝试简单的padding方案
 3. 如果内存资源紧张，或者追求极致性能，可以考虑使用swizzle方案
 4. 最重要的是要实际测试性能，选择最适合你的应用场景的方案
 
 记住，优化不是目的，提高程序性能才是。有时候，简单的方案可能比复杂的方案更好维护，更容易调试，这些都是我们需要权衡的因素。
 
-# 1.9  矩阵规模与性能的关系
-
 # 1.10 总结
 
 本章从GPU的并行性和局域性特性出发，系统地介绍了SGEMM优化的关键要素：从基础的三重循环实现到内存聚合优化，从分块矩阵乘法到层次化张量表示，从共享内存冲突处理到矩阵规模影响分析。这些优化技术构成了一个层次分明的知识体系，帮助我们深入理解现代GPU架构下的性能优化原理。
 
 近年来GPU呈现出一个显著趋势：计算能力（屋顶）快速提升，而内存带宽（屋檐斜率）增长相对缓慢。即硬件的算力带宽比在不断提高。这意味着充分利用计算能力需要我们更加注重数据的高效复用，这也是为什么本章详细讨论了分块策略、层次化优化等技术的原因。
 
-从认知科学的角度看，GPU性能优化涉及两种思维模式：一是利用屋檐模型进行快速的性能估算和决策（类比于人类的"快思考"系统），二是通过深入分析发现新的性能瓶颈和优化机会（类比于"慢思考"系统）。这种双系统思维框架不仅帮助我们更好地理解和应用本章介绍的各种优化技术，也为后续的性能优化工作提供了清晰的方法论指导。优秀的性能工程师正是需要在这两种思维模式之间灵活切换，既能快速决策，又能深入思考。
+从认知科学的角度看，GPU性能优化涉及两种思维模式：一是利用屋檐模型等已知性能上界，进行快速的性能估算和决策（类比于人类的"快思考"系统），二是通过深入分析发现新的性能上界和优化机会（类比于"慢思考"系统）。这种双系统思维框架不仅帮助我们更好地理解和应用本章介绍的各种优化技术，也为后续的性能优化工作提供了清晰的方法论指导。优秀的性能工程师正是需要在这两种思维模式之间灵活切换，既能快速决策，又能深入思考。
+
+![alt text](figure/ch4/bounds.png)
 
 # 参考资料
 
diff --git a/guide/figure/ch4/app.png b/guide/figure/ch4/app.png
new file mode 100644
index 0000000000000000000000000000000000000000..a6a4d30d1dea9fcbadf9c0601ff6937fa38889ac
Binary files /dev/null and b/guide/figure/ch4/app.png differ
diff --git a/guide/figure/ch4/major.png b/guide/figure/ch4/major.png
new file mode 100644
index 0000000000000000000000000000000000000000..08970bbffca00de979ce02f75c56fb3bbbf5ae20
Binary files /dev/null and b/guide/figure/ch4/major.png differ
diff --git a/guide/figure/ch4/region.png b/guide/figure/ch4/region.png
index d5ce4372f4c42b4be176a1c9cfb40c2d9601877e..a6a4d30d1dea9fcbadf9c0601ff6937fa38889ac 100644
Binary files a/guide/figure/ch4/region.png and b/guide/figure/ch4/region.png differ
diff --git a/guide/figure/ch4/step-by-step.png b/guide/figure/ch4/step-by-step.png
new file mode 100644
index 0000000000000000000000000000000000000000..a16735ab5e8937ad8a6bab66ca68fea42f56ff28
Binary files /dev/null and b/guide/figure/ch4/step-by-step.png differ
diff --git a/plot/roofline/roofline_MI300X.py b/plot/roofline/roofline_MI300X.py
new file mode 100644
index 0000000000000000000000000000000000000000..7d36ef2c2b55264d892105c74b46da161003bfb0
--- /dev/null
+++ b/plot/roofline/roofline_MI300X.py
@@ -0,0 +1,63 @@
+import numpy as np
+import matplotlib.pyplot as plt
+
+def plot_roofline(memory_bandwidth, arithmetic_throughput, tile_map, title):
+    # Find maximum bandwidth and maximum throughput
+    B_max = max(memory_bandwidth.values())
+    T_max = max(arithmetic_throughput.values())
+    
+    # Define operation intensity range with higher resolution
+    op_intensity = np.logspace(-2, 5, 1000)
+    
+    plt.figure(figsize=(8, 6))
+    
+    # Plot memory bandwidth lines, cut off to T_max
+    for mem_type, bandwidth in memory_bandwidth.items():
+        intercept_intensity = T_max / bandwidth
+        roof = bandwidth * op_intensity / 1e12  # Convert to TB/s
+        end_idx = np.where(op_intensity <= intercept_intensity)[0][-1]  # Find ending index based on intercept_intensity
+        plt.loglog(op_intensity[:end_idx+1], roof[:end_idx+1], label=f'{mem_type} ({bandwidth / 1e9:.0f} GB/s)', linestyle='-', linewidth=2)
+    
+    # Plot arithmetic throughput lines, starting at the intersection with the leftmost bandwidth line
+    for arith_type, throughput in arithmetic_throughput.items():
+        intercept_intensity = throughput / B_max
+        roof = np.full_like(op_intensity, throughput / 1e12)  # Convert to TFLOPS
+        start_idx = np.where(op_intensity >= intercept_intensity)[0][0]  # Find starting index based on intercept_intensity
+        plt.loglog(op_intensity[start_idx:], roof[start_idx:], label=f'{arith_type} ({throughput / 1e12:.0f} TFLOPS)', linestyle='-', linewidth=2)
+    
+    # Add vertical dotted lines for specific HGEMM tile sizes
+    for tile_name, tile_intensity in tile_map.items():
+        plt.axvline(x=tile_intensity, color='black', linestyle=':', linewidth=2)
+        plt.text(tile_intensity, 1, tile_name, rotation=90, verticalalignment='bottom', horizontalalignment='right')
+    
+    # Set plot labels and title
+    plt.xlabel('Operation Intensity (OPs:Bytes)')
+    plt.ylabel('Attainable Performance (TFLOPS)')
+    plt.title(title)
+    plt.legend(loc='lower right')
+    plt.grid(True, which="both", ls="--")
+    plt.xlim(1e-1, 1e4)
+    # plt.ylim(1, T_max / 1e12 * 1.1)  # Set y limit a bit above T_max for clarity
+    plt.ylim(1, 400)
+    plt.savefig('roofline_MI300X.png')
+    plt.show()
+
+# Example parameters for C500
+memory_bandwidth = {
+    'DRAM': 1.8e12,  # 1.8 TB/s
+    'L2': 4.608e12,     
+    'shared': 14.9e12  
+}
+
+arithmetic_throughput = {
+    'FP32 Vector': 15e12,  
+    'FP16 Tensor': 240e12  # 312 TFLOPS
+}
+
+tile_map = {
+    # '(128x128)': 128,
+    # '(64x128)': 85,
+    # '(64x64)': 64
+}
+
+plot_roofline(memory_bandwidth, arithmetic_throughput, tile_map, 'C500 Hierarichal Roofline')