HPC

本文简要介绍通用矩阵乘（GEMM，General Matrix Multiplication）优化的基本概念和方法。GEMM 是 HPC 领域中最基础且计算密集型的工作负载之一。在人工智能、科学模拟和图像处理等领域，它的性能直接影响着整个应用程序的效率。虽然其数学概念简单，但高效的 GEMM 实现却需要对计算机体系结构有深刻的理解，包括缓存、SIMD 指令集和并行化技术。 Naive GEMM GEMM 通常定义为 $ C = A \times B $，对于矩阵 $ A \in \mathbb{R}^{M \times K} $，矩阵 $ B \in \mathbb{R}^{K \times N} $，其乘积矩阵 $ C\in \mathbb{R}^{M \times N} $ 可以表示为 $$ C_{i,j} = \sum_{k=0}^{K-1} A_{i,k}\times B_{k,j} $$ 对应的朴素代码通常如下（默认行主序存储）： 1 2 3 4 5 6 7 8 9 10 11 12 13 14 void gemm_naive(int M, int N, int K, const float* A, const float* B, float* C) { for (int i = 0; i < M; ++i) { for (int j = 0; j < N; ++j) { C[i][j] = 0.0f; // 初始化 C[i][j] } } for (int i = 0; i < M; ++i) { for (int j = 0; j < N; ++j) { for (int k = 0; k < K; ++k) { C[i][j] += A[i][k] * B[k][j]; } } } } 分析： ...

如何更好更快地访问内存是 HPC 中最大的瓶颈之一，仅仅了解 SIMD 或并行编程接口是不足够的，本文将梳理计算机的内存层次结构、缓存友好编程、内存墙现象、NUMA 架构以及预取技术。 Understanding Memory Hierarchy 为了充分利用现代 CPU 的性能，我们必须理解数据是如何在不同层级的内存组件之间流动的。 Registers, Caches, and Main Memory 寄存器 (Registers)： CPU 内置的、容量最小但速度最快的数据存储单元，用于存储正在被 CPU 活跃操作的数据。CPU 直接在寄存器上执行大部分计算。 缓存 (Cache)： 位于 CPU 和主内存之间的小容量、高速存储区域。它们的目的是通过存储最可能被 CPU 再次访问的数据来减少对主内存的访问延迟。 L1 缓存 (Level 1 Cache)：最小、最快，通常分为数据缓存 (L1d) 和指令缓存 (L1i)，每个 CPU 核心独有。其访问速度与 CPU 核心时钟周期相近。 L2 缓存 (Level 2 Cache)：比 L1 大且慢，每个 CPU 核心独有或由几个核心共享。 L3 缓存 (Level 3 Cache)：最大、最慢的缓存，通常由同一 CPU 插槽上的所有核心共享。 主内存 (Main Memory/RAM)： 容量远大于缓存，但访问速度慢得多。当数据不在任何缓存中时，CPU 必须从主内存中获取。 TLB (Translation Lookaside Buffer)： TLB 是一个专用的高性能缓存，用于存储虚拟地址到物理地址的转换映射。当 CPU 访问一个虚拟地址时，它首先检查 TLB。如果找到对应的物理地址（TLB 命中），则可以快速进行内存访问；如果未找到（TLB 未命中），则需要查询页表，这将导致显著的延迟。理解 TLB 对于优化内存页访问模式，尤其是在处理大型数据集时至关重要。 ...

去年超算队招新唯一没有解决的一道题，今在 Gemini 老师的帮助下成功解决，决定重写一份题解报告。 Description 题目传送门 题目要求参赛者优化一个 C 语言函数 rotate_the_bit_vector，函数功能是对一个 bit vector 中的一个指定子区间进行循环旋转操作。 具体而言，题目给的 bit_vector 是一种内存紧凑的数据结构，将 8 个 bit 打包存储到一个字节中。参赛者需要在只修改 submit_func.c 一个文件的前提下重写其中的 rotate_the_bit_vector 函数，使其在大规模数据时尽可能快。 最后评分程序会通过三个不同的 benchmark (-s, -m, -l) 来衡量，每个测试中的数据规模会随着层数的增加而几何级增长，最终的得分取决于规定时间内能到达的“层数”，层数越高。说明性能越好，最终分数也越高。 Analysis Three-Reversal Algorithm 假设要移动的区间长度为 $ n $ ，需要移动 $ k $ 位；由于向右旋转 $ k $ 位可以等效于向左旋转 $ n-k $ 位，因此只讨论向左的移动。 题目提供了一个初始的性能极差的实现，通过逐位移动 $ k $ 次来实现 $ k $ 位循环旋转，复杂度为 $ O(n^{2}) $。根据这个原始实现很容易想到一个初步的优化方案： 问题的核心是把数组 [A|B] 变成 [B|A] ，一个经典的算法是三步翻转法：如果我们把翻转操作记为 ' ，A' 表示数组 A 前后反转，那么可以发现原问题的操作实际上等价于三次翻转操作：[A'|B']' = [B|A] ，复杂度为 $ O(n) $ ，代码如下： ...

1. What is SIMD? SIMD，即 Single Instruction Multiple Data ，是一种并行计算的模式。传统的单指令单数据模型，也就是一条指令 CPU 只能处理一份数据，这在科学计算和图像渲染等大量数据密集的任务中是非常低效的。 SIMD 的核心思想是用一条指令同时对多个数据进行操作，现代的 CPU 为此设计了特殊的硬件单元，包括宽位（比如 128、256 或 512 位）的向量寄存器 (Vector Registers) 和能够操作这些寄存器的向量指令 (Vector Instructions)。一个向量操作可以同时完成多个标量操作，从而实现数据并行 (Data Parallelism)，提高效率。假设一个 256 位的向量寄存器可以容纳 8 个 32 位浮点数，一条向量加法指令就可以一次性完成 8 个浮点数的加法，理论上将这部分计算的吞吐量提升至原来的 8 倍；并且相比于执行 8 条独立的标量加法指令，CPU 只需要获取并解码一条向量加法指令，这降低了指令流水线的压力。 2. How SIMD Works 要理解 SIMD 的工作原理，需要了解两个核心概念：向量寄存器和向量指令。 2.1. Vector Registers 向量寄存器是 CPU 内部的特殊存储单元，其宽度远大于通用寄存器。不同的 Instruction Set Architecture (ISA, 指令集架构) 提供了不同宽度和名称的向量寄存器。 SSE (Streaming SIMD Extensions)：提供了 128 位的 XMM 寄存器。 AVX (Advanced Vector Extensions)：提供了 256 位的 YMM 寄存器。 ...

1. Why Memory Performance Matters in HPC? 在 HPC 领域，我们常常关注 CPU 的浮点运算能力 (FLOPS)，但真正的性能瓶颈往往不在于计算本身，而在于数据访问。现代 CPU 的计算速度远超于内存的访问速度，这种差距被称为内存墙 (Memory Wall)。程序的大部分时间可能都消耗在等待数据从内存加载到 CPU 寄存器的过程中。因此，优化内存访问模式，最大限度地利用 Cache，是提升 C/C++ 程序性能至关重要的一环。 2. Memory Alignment 内存对齐是指一个数据对象的内存地址是其自身大小或特定字节数（通常是 2 的幂）的整数倍。例如一个 4 字节的 int 类型变量，如果其内存地址是 4 的倍数（如 0x...00, 0x...04, 0x...08），那么它就是内存对齐的。 2.2 Why is Alignment Important? CPU 并不是逐字节地从内存中读取数据，而是以块（通常是缓存行 (Cache Line)，例如 64 字节）为单位进行读取。 性能提升：如果一个数据跨越了两个缓存行，CPU 就需要执行两次内存读取操作才能获取这一个数据，这会浪费一倍的时间。如果数据是对齐的，就可以保证它完整地落在一个缓存行内，CPU 只需一次读取操作。 硬件要求：许多现代 CPU 指令集，尤其是用于并行计算的 SIMD 指令强制要求操作的数据必须是内存对齐的，对未对齐的数据执行这些指令可能会导致程序崩溃或性能急剧下降。 2.3 How to Achieve Alignment in C/C++? C++11 alignas：这是 Modern C++ 的标准方式，可以指定变量或类型的对齐要求。 1 2 3 4 5 6 7 8 // 声明一个按 64 字节对齐的数组 alignas(64) float aligned_array[1024]; // 定义一个结构体，使其每个实例都按 32 字节对齐 struct alignas(32) MyData { float a; int b; }; GCC/Clang __attribute__((aligned(N)))：特定于编译器的扩展。 1 2 // 声明一个按 64 字节对齐的数组 float aligned_array[1024] __attribute__((aligned(64))); 动态内存对齐：标准的 malloc 不保证特定的对齐方式（通常只保证基本类型的对齐）。需要使用专用函数。 1 2 3 4 5 6 #include <stdlib.h> // C11 标准 // 分配 1024 个 float，并按 64 字节对齐 float* dynamic_array = (float*)aligned_alloc(64, 1024 * sizeof(float)); free(dynamic_array); // 必须用 free 释放 3. Data Locality 数据局部性是缓存工作的基本原理，也是性能优化的核心。描述了 CPU 访问内存地址的集中程度。 ...

HPC 领域中，除了基于共享内存的 OpenMP, 还有一种更广泛应用于分布式内存系统的并行编程范式——消息传递接口 (MPI)。MPI 不依赖于共享内存，而是通过进程间的显式消息传递来实现数据交换和同步，从而能支持更大规模的集群计算，是构建大规模 HPC 集群不可或缺的工具。 1. What is MPI? MPI (Message Passing Interface) 是一种用于分布式内存系统并行编程的标准化通信协议和库函数规范。它定义了一套可移植的函数接口，允许在并行计算环境中独立运行的进程之间进行消息传递，从而实现数据交换和协同工作。MPI 不指定如何启动进程，也不要求所有进程在同一台机器上，这使得它非常适合用于集群或多节点环境中的大规模并行计算。 2. The MPI Programming Model 分布式内存模型 在分布式内存模型中，各个处理节点可以独立运行自己的进程，使用自己的本地内存来存储和处理数据。每个进程的内存是私有的，其他进程无法直接访问它们。如果一个进程需要访问另一个进程的数据，就必须通过显式的消息传递机制将数据从一个进程发送到另一个进程。同一个节点（服务器）内部需要借助高速数据总线等硬件实现，而跨节点的通信通常由网络连接来实现，比如通过高速以太网、IB（InfiniBand）等。 核心概念 进程 (Process)：一个 MPI 程序由一个或多个独立的进程组成。这些进程通过调用 MPI 库函数来进行通信。 通信子 (Communicator)：一个通信子（MPI_Comm）定义了一个可以互相通信的进程组。最常用的通信子是 MPI_COMM_WORLD，它包含了程序启动时的所有进程。 秩 (Rank)：在同一个通信子内，每个进程都被赋予一个唯一的整数标识，称为秩。秩的范围是从 0 到 进程总数 - 1。 消息传递 (Message Passing)：进程间通信的核心机制，分为两大类： 点对点通信 (Point-to-Point)：在两个指定的进程之间进行。 集体通信 (Collective)：在一个通信子内的所有进程共同参与的通信。 通信协议：MPI 提供了多种通信协议，如阻塞通信（Blocking）、非阻塞通信（Non-blocking）、同步通信（Synchronous）等。 3. Basic Functions and Concepts 一个基础的 MPI 程序总是包含初始化、执行并行代码和结束这几个部分。 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 #include <mpi.h> #include <stdio.h> int main(int argc, char** argv) { // 1. 初始化 MPI 环境 MPI_Init(&argc, &argv); int world_size; int world_rank; char processor_name[MPI_MAX_PROCESSOR_NAME]; int name_len; // 2. 获取通信子信息 MPI_Comm_size(MPI_COMM_WORLD, &world_size); // 获取总进程数 MPI_Comm_rank(MPI_COMM_WORLD, &world_rank); // 获取当前进程的秩 // 获取处理器名称 (可选) MPI_Get_processor_name(processor_name, &name_len); // 3. 基于秩执行不同的代码 printf("Hello world from processor %s, rank %d out of %d processors\n", processor_name, world_rank, world_size); // 4. 结束 MPI 环境 MPI_Finalize(); return 0; } MPI_Init()：初始化 MPI 执行环境，必须是第一个被调用的 MPI 函数。 MPI_Comm_size()：获取指定通信子（这里是 MPI_COMM_WORLD）中的总进程数。 MPI_Comm_rank()：获取当前进程在指定通信子中的秩。 MPI_Finalize()：清理并结束 MPI 环境，必须是最后一个被调用的 MPI 函数。 4. Point-to-Point Communication 点对点通信是 MPI 中最基本的通信模式，用于在一个进程向另一个进程发送数据。核心操作是 Send 和 Recv。 ...

由于高性能计算场景下的并行编程任务的特性，OpenMP 可以通过简单受限的语法极大地化简了并行编程的复杂性，在普通的串行代码中添加一些指令就能够实现高效并行化。 1. What is OpenMP? OpenMP (Open Multi-Processing) 是一种用于共享内存多处理器系统并行编程的 API。它通过在 C, C++, 或 Fortran 代码中添加 #pragma 的方式，让开发者可以轻松地将串行代码并行化，而无需手动管理复杂的线程创建、同步和销毁过程。 2. The OpenMP Programming Model 共享内存模型：所有线程在同一个地址空间中共享数据。这意味着不同线程可以访问相同的内存位置，并且可以共享变量的值。 共享变量：在并行区域中，默认情况下，大多数变量是共享的，即所有线程都可以访问和修改这些变量的值。 私有变量：某些情况下，我们可能希望每个线程拥有变量的私有副本，这样不同线程之间不会相互干扰。OpenMP 通过 private 指令指定这些变量。 数据竞争（Race Condition）：由于多个线程同时访问和修改共享变量，可能会导致数据竞争问题。为了避免这种情况，OpenMP 提供了同步机制，如 critical 和 atomic 等。 并行区域（Parallel Region）：是 OpenMP 编程的核心概念。它是由编译器指令 #pragma omp parallel 指定的一段代码，告诉 OpenMP 在这段代码中创建多个线程并行执行。 Fork-Join 执行模型：从单线程开始执行，进入并行区域开始并行执行，在并行区域结尾进行同步和结束线程。 3. Core Directives and Constructs OpenMP 的功能主要是通过编译指令（Directives）和相关的子句（Clauses）来实现的。 parallel：用于创建一个并行区域。 1 2 3 4 5 #pragma omp parallel { // 这部分代码将由多个线程同时执行 printf("Hello from thread %d\n", omp_get_thread_num()); } for：用于并行化 for 循环，必须与 parallel 结合使用。它会自动将循环迭代分配给不同的线程，这是 OpenMP 最常用、最高效的指令之一。 1 2 3 4 #pragma omp parallel for for (int i = 0; i < n; i++) { // 循环的 n 次迭代会被分配给不同线程 } sections：用于将不同的、独立的任务代码块分配给不同线程。适用于任务并行而不是数据并行。 1 2 3 4 5 6 7 #pragma omp parallel sections { #pragma omp section { /* task A */ } #pragma omp section { /* task B */ } } 4. Data Scoping 数据作用域定义了并行区域中变量如何被线程共享或者私有，OpenMP 通过子句 clauses 来控制变量属性。 ...