使用 Intel SIMD 的段错误，即使空间非常大并且是 32 字节的倍数

Posted 2023-02-14

【中文标题】使用 Intel SIMD 的段错误，即使空间非常大并且是 32 字节的倍数

【英文标题】：Segment fault using Intel SIMD, even the space is very large and is of multiple of 32 bytes

【发布时间】：2022-01-10 13:43:20

【问题描述】：

在使用 SIMD 指令优化矩阵乘法时，我不断遇到段错误。

这里是核心计算部分。 矩阵存储如下： 分配了一个大小为 (3 * 1025 * 1025) 的大 vector<double> buf。矩阵A从buf[0]开始，矩阵B从buf[1025]开始，C从buf[1025*2]开始。我执行了大小为 4 到 1024 的各种矩阵乘法。所以它们都可以放入这个向量中。

#include <immintrin.h>
#define BLOCK_SIZE 4
/*
 * performs 4 * 4 matrix multiplication C=A*B
 * C is 4-by-4, A is 4-by-4, and B is 4-by-4, column major matrices
 * lda is the size of the large matrix.
 */
static void do_block(int lda4, double* A, double* B, double* C) 
    int n=4;
    for(int i=0; i<n; i++) // process i th column
      for(int j=0; j<n; j++)
        __m256d c = _mm256_load_pd(C+j*lda);
        c = _mm256_fmadd_pd(_mm256_load_pd(A+i*lda), _mm256_broadcast_sd(B+i+j*lda), c);
        _mm256_store_pd(C+j*lda, c);
      
    


/* This routine performs a dgemm operation
 *  C := C + A * B
 * where A, B, and C are lda-by-lda matrices stored in column-major format.
 * On exit, A and B maintain their input values. */
void square_dgemm(int lda, double* A, double* B, double* C) 
    for (int j = 0; j < lda; j += BLOCK_SIZE) 
        // Accumulate block dgemms into block of C
        for (int k = 0; k < lda; k += BLOCK_SIZE) 
            // For each block-row of A
            for (int i = 0; i < lda; i += BLOCK_SIZE) 
                do_block(lda, A + i + k * lda, B + k + j * lda, C + i + j * lda);
            
        
    

奇怪的是： 当我将向量从 (3 * 1025 * 1025) 的大小更改为 (3 * 1024 * 1024) 时，它给了我段错误。

我的问题是：

我了解到these instructions require aligned data。实际上，用 _mm256_loadu_pd 等未对齐的版本替换可以消除此错误。但是，由于 (3 * 1024 * 1024 * sizeof(double)) % 32 bytes == 0 的大小，它不是 32 字节对齐还是我误解了这个概念？ 我分配了非常大的连续空间，为什么在执行 small mat mul (4*4) 时它从一开始就崩溃？我以为只要我调用_mm256_load_pd(addr) 并从addr 开始分配至少32 个字节，它就不会崩溃，我错了吗？ 为什么它不会在 (3 * 1025 * 1025) 的 buf 上崩溃，但在 (3 * 1024 * 1024) 上崩溃？当大小为奇数（如 1025、1027、1029）时似乎不会崩溃，而当数字为偶数时（如 1024、1026）总是崩溃。

代码是使用 GCC 编译的，带有 -march=native 和 -O3。 CPU 支持 FMA、AVX 和 AVX2。这台机器是 Google Cloud VM，CPU 是 Intel Xeon，我无法获得确切的型号。谢谢你的建议！

【问题讨论】：

如果您查找sigaction（SO 上可能有示例），您可以为 SIGSEGV 安装一个处理程序，该处理程序提供了一个siginfo_t，它详细描述了故障的原因。您可能需要稍微熟悉一下/usr/include/<arch-toolset>/asm/siginfo.h 才能解析信息字段。

或者更简单地说，在调试器中运行程序，它会提供所有这些信息，将其与行号相关联，并提供一个方便的界面来检查程序在崩溃时的整个状态。

我认为vector<double> 不保证任何超出 8 字节对齐的内容，因此传入的指针 A,B,C 可能无法可靠地开始正确对齐。尝试在输入时打印出这些指针值。反正buf[0]和buf[1025]不能都对齐到32字节？

相关：How to solve the 32-byte-alignment issue for AVX load/store operations? - 如果您不将行步幅与实际矩阵几何分开填充为 4 个双精度数（32 个字节）的倍数，您的某些负载可能需要使用 _mm256_loadu_ps .

@PeterCordes - 因为这个问题不需要担心旧的编译器，所以不需要我的答案的更新版本。

【参考方案1】：

感谢 cmets，我想我已经弄清楚出了什么问题。

对 Q1 的回答：是的，我误解了“对齐”的概念。我必须确保地址是 32 的倍数。所以在手动（不是通过指针数学的优雅方式）从 36 字节对齐空间开始制作 3 个矩阵之后，它就可以工作了。

对Q2的回答：由于传递给_mm256_load_pd(addr)的地址不是32字节对齐的，即使它很大，也可能引发异常，根据Intel documentation

对 Q3 的回答：这很棘手：A 和 C 并非都是 32 字节对齐的。但是，更改总大小将导致 C 对齐 32 个字节（但 A 仍未对齐），并且不会崩溃。这只是偶然的，它不安全，因为 A 没有对齐。

所以矩阵乘法的向量化很麻烦，我需要

确保矩阵 A 和 C 以 32 个字节对齐，可以使用 cmets 中提到的 aligned_alloc 之类的内容。 将矩阵填充为 4 的倍数，以确保在调用 _mm256_load_pd(addr) 的每一步中，地址都正确对齐。