前言

GEMM（General Matrix Multiply，通用矩阵乘法）是深度学习、AI 推理、科学计算中最核心的计算原语——没有之一。几乎所有主流模型（Transformer、CNN、MoE）的计算瓶颈，最终都落在这一步。

昇腾 CANN 的 ops-blas 是第 2 层 AOL（Ascend Operator Library）中的高性能线性代数算子库，提供 BLAS 级别的基础算子和轻量化 GEMM 调用接口 aclBLASLt。与 ops-nn 中的融合 MatMul 不同，ops-blas 的设计目标是极致的硬件贴合：把每个 Tile 的大小、每条数据的预取时机、每次 Cache 行的复用率都算清楚，不留浪费。

这篇文章就来拆解 ops-blas 中 GEMM 的分块策略，以及它如何在昇腾达芬奇架构的 L1/L2 Cache 层级上把算力压榨到极限。

---

1. GEMM 的计算特征：M×K × K×N 的分层拆分

GEMM 计算 $C=A\times B+D$，其中 $A\in {\mathbb{R}}^{M\times K}$，$B\in {\mathbb{R}}^{K\times N}$，$C\in {\mathbb{R}}^{M\times N}$。表面上看，这只是一个三重循环：

for (int m = 0; m < M; ++m)
  for (int n = 0; n < N; ++n)
    for (int k = 0; k < K; ++k)
      C[m][n] += A[m][k] * B[k][n];

但这三重循环的顺序直接决定了数据搬运路径。实际上，GEMM 在硬件上被拆解为两级分块（Blocking）：

1.1 外层分块：MR×NR 大小的 Micro-Tile

昇腾 AICore 的 Cube 计算单元一次能吞吐一块 $16\times 16$ 或 $32\times 16$ 的矩阵块。整个 $M\times N$ 的结果矩阵被划分为若干 $MR\times NR$ 的子块，每个子块交给一个计算单元完成。这叫 Outer Blocking。

1.2 内层 K 循环：累积计算的分块

每次只取一小段 K（记为 $BK$）来计算 partial sum：

$$C_{M_i\times N_j}+=A_{M_i\times K_s}\times B_{K_s\times N_j}$$

随着 K 循环推进，$C_{M_i\times N_j}$ 不断累加，最终得到完整结果。内层 K 循环的分块大小 $BK$ 直接决定了 $A$ 和 $B$ 的数据复用窗口——$BK$ 越大，同一块 $B$ 被重读次数越多，但 Cache 容量的压力也越大。

ops-blas 的 GEMM 实现遵循这一分层思想，只不过每一层的分块大小都需要根据昇腾硬件的 Cache 层级精确推导。

---

2. 分块策略的硬件推导：L1/L2 Cache 层级适配

2.1 昇腾达芬奇架构的存储层次

昇腾 AICore 的存储层次如下（从快到慢）：

层级	类型	典型容量	带宽（相对）
寄存器文件	RF	~2KB/Thread	1x（基准）
L1 Cache	WMMA/SRAM	64KB~256KB（可配置）	~10x HBM
L2 Cache	L2 SRAM	8MB~16MB	~3x HBM
HBM	全局显存	32GB~256GB	1x（基准）

核心矛盾：Cube 单元的算力远高于 HBM 的带宽。Ascend 910 的 Cube 算力标称 256 TFLOPS（FP16），但 HBM 带宽只有 1.2 TB/s。这意味着如果每个浮点运算都要从 HBM 取数据，算力利用率连 10% 都不到。

分块策略的本质，就是用更大的局部计算换更少的数据搬运。

2.2 L1 Cache 约束下的 M/N 分块推导

L1 Cache 的容量限制了单个 Tile 的数据占用。对于 $MR\times BK+BK\times NR$ 的数据总量，必须满足：

$$\text{TileSize}\le \text{L1CacheSize}\times \text{UtilizationFactor}$$

典型配置下，昇腾 AICore L1 可用空间约 64KB（可调整），ops-blas 中一个 FP32 GEMM Tile 的典型配置为：

// blasLt/matmul_fp32/arch35/matmul_fp32_kernel.cpp（典型结构示意）
constexpr int M_TILE = 64;   // M方向分块
constexpr int N_TILE = 64;   // N方向分块
constexpr int K_BLOCK = 32;  // K方向内块大小

其中 $A$ 子块 $64\times 32$（FP32，4B/elem）= 8KB，$B$ 子块 $32\times 64$ = 8KB，加上累加器缓冲区 $C$ 子块 $64\times 64$ = 16KB，单个活跃 Tile 的在片数据约 32KB，远低于 L1 容量，留出了一部分空间给指令预取和中间变量。

2.3 L2 Cache 约束下的 K 分块推导

L2 Cache 的容量（通常 8MB~16MB）决定了在片（on-chip）可容纳的 $A$ 和 $B$ 的总量。当矩阵很大时，整个 $A$ 的一个 $M\times K$ 分块和 $B$ 的 $K\times N$ 分块可能都放不下——这时需要把 $K$ 方向再切一刀，形成 K-Split：

for (int k_offset = 0; k_offset < K; k_offset += K_OUTER) {
    // 将K分成若干 K_OUTER 大块，每块内再细分为K_INNER
    for (int ki = 0; ki < K_OUTER; ki += K_INNER) {
        ComputeTile(A[M_TILE][ki:ki+K_INNER], 
                    B[ki:ki+K_INNER][N_TILE],
                    C[M_TILE][N_TILE]);
    }
    // K_OUTER 块完成后才写回C
}

ops-blas 的 blasLt/epilogue 模块处理 K 分块间的累积语义，确保 partial sum 在不同 K 块之间正确累加。

---

3. ops-blas 中的分块实现

3.1 Tile 大小的硬件自适应选择

ops-blas 通过 aclblasLtMatmulPreference_t 提供 Heuristic 搜索，自动在不同硬件形态（Ascend 910 / 910B / 950）上选择最优 Tile 配置：

// include/cann_ops_blasLt.h
aclblasLtMatmulHeuristicResult_t heuristicResults[resultCount];
aclblasLtMatmulPreference_t preference;

// 偏好workspace大小（bytes），限制搜索空间
aclblasLtMatmulPreferenceSetAttribute(
    preference,
    ACLBLASLT_MATMUL_PREF_MAX_WORKSPACE_BYTES,
    &maxWorkspaceBytes, sizeof(maxWorkspaceBytes));

// 启发式搜索：在给定workspace约束下找出最优算法配置
aclblasLtMatmulAlgoGetHeuristic(
    handle, layoutA, layoutB, layoutC, layoutD,
    preference, resultCount, heuristicResults, &returnedCount);

aclblasLtMatmulHeuristicResult_t 中包含 wavesCount——这是设备利用率指标，描述该算法配置在目标硬件上需要多少waves来填充所有计算单元。wavesCount 越小，说明硬件利用率越高。

3.2 数据预取：双缓冲（Double Buffer）流水线

光分块不够，分块之间如果串行执行，Cube 单元在等待数据加载时会空闲。ops-blas 引入 双缓冲流水线：

// 预取逻辑示意（伪代码）
void GemmWithDoubleBuffer(..., int K_BLOCK) {
    // Buffer A[0], B[0] 已就绪，计算 C[0]
    ComputeTile(A_buf[0], B_buf[0], C);

    // 与计算并行：预取下一块数据到 A[1], B[1]
    __ascend_load_async(A_buf[1], A_next_ptr, ...);
    __ascend_load_async(B_buf[1], B_next_ptr, ...);

    // 第一块计算完成后，交换 buffer，重复
    Swap(A_buf[0], A_buf[1]);
    Swap(B_buf[0], B_buf[1]);
}

__ascend_load_async 是昇腾异步加载指令，触发 DMA 传输后立即返回，不阻塞 Cube 计算。两者流水重叠，消除了 $O(K)$ 次的串行等待开销。

3.3 Packing 布局：列优先的连续内存访问

GEMM 的 $A$ 矩阵在内存中通常按行主序存储，但 Cube 计算单元在读取 $A$ 的某一列片段（$A[m][k0:k0+BK]$）时，如果数据不连续，会产生大量非合并访问（non-coalesced access），导致 HBM 带宽利用率骤降。

ops-blas 引入 Packing：在计算前将 $A$ 和 $B$ 的待计算 Tile 重新排列为计算友好的内存布局：

// Packing A: 从 row-major 转换为 tile-column-major
void PackMatrixA(float* src, float* packed, 
                 int M_TILE, int K_BLOCK, int K_STRIDE) {
    for (int m = 0; m < M_TILE; ++m) {
        for (int k = 0; k < K_BLOCK; ++k) {
            // 按K方向连续存储，匹配Cube的列读取模式
            packed[m * K_BLOCK + k] = src[m * K_STRIDE + k];
        }
    }
}

// Packing B: 类似，转换为 tile-row-major
void PackMatrixB(float* src, float* packed,
                 int K_BLOCK, int N_TILE, int N_STRIDE) {
    for (int k = 0; k < K_BLOCK; ++k) {
        for (int n = 0; n < N_TILE; ++n) {
            packed[k * N_TILE + n] = src[k * N_STRIDE + n];
        }
    }
}

Packing 后，$A$ 的一个 Tile 变成 $M\_TILE\times K\_BLOCK$ 的连续内存块，$B$ 的 Tile 变成 $K\_BLOCK\times N\_TILE$ 的连续块。Cube 单元每次读取 $K\_BLOCK$ 个元素时，地址连续，HBM 的 Burst Read 效率达到最高。

Packing 的代价是额外的内存拷贝开销（约 5%~15% 的数据复制），但这笔开销远小于它换来的 Cache 命中率和带宽利用率提升，净收益为正。

3.4 Epilogue：融合后处理的分块策略

ops-blas 的 blasLt/epilogue 模块在每个 Tile 计算完成后立即执行后处理（ReLU、GELU、Bias、Quantize），而不需要单独的 Kernel 启动：

// 典型的 epilogue 融合配置
typedef enum aclblasLtEpilogue {
    ACLBLASLT_EPILOGUE_DEFAULT = 1,    // 仅缩放/量化
    ACLBLASLT_EPILOGUE_RELU = 2,       // +ReLU
    ACLBLASLT_EPILOGUE_BIAS = 4,       // +Bias
    ACLBLASLT_EPILOGUE_GELU = 32,      // +GELU 激活
    ACLBLASLT_EPILOGUE_GELU_BIAS = 36, // Bias + GELU
    // ...
} aclblasLtEpilogue_t;

这意味着对一个典型Transformer层的 GEMM计算，不需要额外启动独立的 Bias-add Kernel 或 ReLU Kernel，一次 Tile 计算完成，数据还在 L1 Cache 里，后处理就做完了——省掉一次完整的数据回写 HBM 再重新读取的过程。

---

4. Cache 利用率分析

4.1 Cache Miss Rate 的量化

分块策略的效果最终体现在 Cache Miss Rate 上。以 FP32 GEMM 为例，典型场景下的数据访问模式分析：

场景：$M=1024,N=1024,K=1024$，Tile 大小 $64\times 64\times 32$

数据块	无分块 HBM 访问次数	有分块 HBM 访问次数	Cache 命中率
A 矩阵	$M\times N\times K/M_{elem}$	$M\times (K/BK)\times MR\times BK/M_{elem}$	~95%
B 矩阵	$M\times N\times K/N_{elem}$	$N\times (K/BK)\times NR\times BK/N_{elem}$	~92%
C 矩阵	—	写回 $M\times N/(MR\times NR)$ 次	—

B 矩阵的 K 方向复用最关键：同一个 $B$ 的 $K\_BLOCK\times NR$ Tile，在整个 K 循环中被同一组 $N$ 列反复使用。一旦这个 Tile 落在 L1 Cache 里，$K$ 次乘法操作中只有第 1 次触发 HBM 读，后 $K-1$ 次全是 L1 Cache Hit。

4.2 算力利用率对比

以下为理论估算（假设 HBM 带宽利用率从 30% 提升到 80%，其他条件不变）：

指标	无分块基线	分块优化后	提升倍数
HBM 有效带宽	~360 GB/s	~960 GB/s	2.67x
Cube 算力利用率	~25%	~75%	3x
端到端 GEMM 吞吐	64 TFLOPS	192 TFLOPS	3x

ops-blas 的 aclblasLtMatmulHeuristicResult_t.wavesCount 字段可直接反映这一指标：当 wavesCount 从 4.0 降到 1.0，意味着所有计算单元都被填满。

4.3 不同矩阵形状下的 Cache 表现

Cache 利用率并非恒定——它高度依赖矩阵形状：

• $K$ 远大于 $M,N$（如 Transformer 的 Attention QK^T）：B 矩阵的 K 维度复用窗口极大，Cache 命中率高，接近理论最优
• $M,N$ 远大于 $K$（如小 batch 推理）：每个 Tile 的 K 循环次数少，B 的复用次数低，Cache 效果有限
• 方阵（$M=N=K$）：分块策略效果最稳定，$M/N/K$ 三个方向的分块都能充分复用

---

5. 两个关键陷阱

5.1 陷阱一：块大小过大导致 Cache Thrashing

症状：矩阵稍大时性能突然下降，profiling 显示 L1 Miss Rate 反而比小矩阵高。

根因：当 $MR\times K\_BLOCK+K\_BLOCK\times NR$ 的总数据量超过 L1 Cache 容量时，同一 Tile 的不同部分会相互驱逐。A 的前 $K\_BLOCK$ 行刚被加载进来，B 的 $K\_BLOCK$ 列就把它们覆盖了——数据还没用完就被挤出去了，下次用时只能重新从 HBM 读。

ops-blas 的应对：Heuristic 搜索会检测矩阵形状，自动降低 Tile 大小（如 $64\times 64$ 降到 $32\times 32$）来避免 thrashing。开发者也可以通过 aclblasLtMatmulPreferenceSetAttribute 设置 ACLBLASLT_MATMUL_PREF_MAX_WORKSPACE_BYTES，以 workspace 换 Cache 空间：

// 强制使用更小的 tile，以 workspace 内存换 Cache 命中率
size_t maxWorkspace = 32 * 1024 * 1024;  // 32MB workspace
aclblasLtMatmulPreferenceSetAttribute(
    preference,
    ACLBLASLT_MATMUL_PREF_MAX_WORKSPACE_BYTES,
    &maxWorkspace, sizeof(maxWorkspace));

判断方法：观察 wavesCount 是否异常增大，或者通过 CANN Profiling 工具查看 L1 Hit Rate 曲线——thrashing 时 L1 Hit Rate 会随时间剧烈振荡。

5.2 陷阱二：小矩阵的开销失衡

症状：矩阵规模较小时（$M,N,K<128$），分块策略反而比不分块慢。

根因：Packing 需要额外的内存拷贝，Tile 调度有固定的 Kernel Launch 开销。当矩阵足够小，数据本身就能放进 L2 Cache 甚至 L1 Cache，分块带来的"减少 HBM 访问"收益，抵不上 Packing + 调度带来的固定开销：

$$T_{\text{gemm}}=\underset{\text{固定开销}}{\underbrace{T_{\text{pack}}}}+\underset{\text{计算开销}}{\underbrace{N_{\text{tiles}}\times T_{\text{compute}}}}+\underset{\text{启动开销}}{\underbrace{N_{\text{tiles}}\times T_{\text{launch}}}}$$

矩阵很小时，$N_{\text{tiles}}$ 很小，但 $T_{\text{pack}}$ 和 $N_{\text{tiles}}\times T_{\text{launch}}$ 的占比反而很高。

ops-blas 的应对：aclblasLtMatmulAlgoGetHeuristic 内部有一个矩阵规模阈值判断：对于小矩阵，跳过 Heuristic，直接选择一个"零开销"的朴素路径（$MR=M,NR=N,BK=K$，整个矩阵一个 Tile 算完）。这一逻辑封装在 blasLt/utils/ 中的启发式规则里：

// utils/shape_adaptor.cpp（典型逻辑示意）
bool IsSmallMatrix(int M, int N, int K) {
    // 矩阵总元素数量小于阈值时，直接走朴素路径
    return (static_cast<int64_t>(M) * N * K) < SMALL_MATRIX_ELEM_THRESHOLD;
}

判断方法：对不同矩阵规模跑基准测试，观察 $M\times N\times K<128\times 128\times 128\approx 2M$ 的性能曲线是否有拐点。

---

6. 代码示例：端到端调用 aclBLASLt GEMM

以下为完整的 FP32 GEMM 调用流程，展示了 ops-blas 的接口设计：

#include "cann_ops_blasLt.h"

aclblasLtHandle_t handle;
aclblasLtMatmulDesc_t matmulDesc;
aclblasLtMatrixLayout_t layoutA, layoutB, layoutC;
aclblasLtMatmulPreference_t preference;

// 1. 初始化 handle
aclblasLtCreate(&handle);

// 2. 创建 Matmul 描述符
aclblasLtMatmulDescCreate(&matmulDesc, ACLBLASLT_COMPUTE_F32, ACL_ERROR_F32);

// 3. 配置 epilogue：GELU + Bias（典型 Transformer 场景）
aclblasLtMatmulDescSetAttribute(
    matmulDesc, ACLBLASLT_MATMUL_DESC_EPILOGUE,
    &ACLBLASLT_EPILOGUE_GELU_BIAS, sizeof(ACLBLASLT_EPILOGUE_GELU_BIAS));

// 4. 创建矩阵布局（COL32 格式，最适合昇腾Cube的内存布局）
aclblasLtMatrixLayoutCreate(&layoutA, ACL_DATA_TYPE_FLOAT32, M, K, K);
aclblasLtMatrixLayoutCreate(&layoutB, ACL_DATA_TYPE_FLOAT32, K, N, N);
aclblasLtMatrixLayoutCreate(&layoutC, ACL_DATA_TYPE_FLOAT32, M, N, N);

// 5. 创建 preference 并设置 workspace 上限
aclblasLtMatmulPreferenceCreate(&preference);
size_t maxWorkspace = 64 * 1024 * 1024;  // 64MB
aclblasLtMatmulPreferenceSetAttribute(
    preference,
    ACLBLASLT_MATMUL_PREF_MAX_WORKSPACE_BYTES,
    &maxWorkspace, sizeof(maxWorkspace));

// 6. Heuristic 搜索最优算法
aclblasLtMatmulHeuristicResult_t results[4];
int returned = 0;
aclblasLtMatmulAlgoGetHeuristic(
    handle, layoutA, layoutB, layoutC, matmulDesc,
    preference, 4, results, &returned);

// 7. 分配 workspace
void* workspace = nullptr;
if (returned > 0 && results[0].state == ACLBLAS_STATUS_SUCCESS) {
    size_t wsSize = results[0].workspaceSize;
    aclrtMalloc(&workspace, wsSize, ACL_MEM_MALLOC_HUGE_FIRST);
}

// 8. 执行 GEMM
aclblasLtMatmul(
    handle, matmulDesc,
    &alpha,               // 标量乘子
    A_dev, layoutA,
    B_dev, layoutB,
    &beta,                // D矩阵乘子（用于 C = A*B + D 场景）
    D_dev, layoutD,
    C_dev, layoutC,
    &results[0].algo,
    workspace, results[0].workspaceSize,
    stream);

---

7. ops-blas vs ops-nn MatMul：跨代设计对比

维度	ops-nn MatMul	ops-blas aclBLASLt
定位	神经网络融合算子	通用线性代数库
接口	AscendCL 单算子调用	cuBLAS-like API + 启发式搜索
Tile 配置	编译时固定	运行时 Heuristic 自适应
Epilogue	预定义融合路径	可配置（Bias/ReLU/GELU/Quantize）
精度	FP16/BF16/INT8 为主	FP32/FP16/MXFP8/MXFP4
最优场景	大 batch、融合后处理	小 batch、轻量化调用、多精度
硬件感知	通用路径	按 arch35/arch32 分架构特化

ops-nn MatMul 适合模型构建阶段——算子粒度大、融合度高、开销摊销好；ops-blas aclBLASLt 适合底层调用、第三方 BLAS 迁移、以及需要精确控制 Tile 和精度的场景。

---

结尾

GEMM 的分块策略，本质上是在存储层次和计算能力之间做空间换时间的博弈。ops-blas 通过 Tile 大小的硬件自适应选择、双缓冲流水线消除数据等待、Packing 布局优化 HBM 访问效率，以及 Heuristic 搜索规避 Thrashing/小矩阵陷阱，把昇腾达芬奇架构的 Cache 层级潜力压榨到了接近理论极限。

如果你对 ops-blas 的 GEMM 实现感兴趣，推荐进一步阅读：

• ops-blas 仓库首页：查看完整的 blasLt 模块源码和 matmul_fp32/arch35 的分块实现
• catlass 算子模板库：深入了解矩阵乘模板的硬件特化设计
• ops-nn MatMul 算子：对比融合算子与通用 BLAS 的设计取舍

昇腾的算力从来不缺，缺的是让算力真正跑起来的分块策略。