内存是矩阵加法（SIMD 指令）的瓶颈吗？

Posted 2023-02-19

技术标签:

【中文标题】内存是矩阵加法（SIMD 指令）的瓶颈吗？【英文标题】：Is memory a bottleneck in matrix addition (SIMD Instructions)? 【发布时间】：2021-03-05 00:14:22 【问题描述】：

我正在尝试使用 SIMD 指令（_mm256_add_pd、存储、加载等）优化 C 中的二维矩阵加法。但是，我根本没有看到很大的加速。使用一些时序代码，我看到在 0.8x-1.5x 范围内的加速是天真的解决方案）。我想知道这是否是典型的？我在想这可能是一个内存瓶颈，因为在这种情况下计算似乎很少。我相信这应该让我的速度提高 4 倍左右，因为我一次要进行 4 次加法，所以我不完全确定瓶颈是什么。

我编写了一些代码来演示我在做什么（测试并行 + SIMD 与仅 SIMD）：

#include <stddef.h>
#include <stdio.h>
#include <stdlib.h>
#include <omp.h>
#include <time.h>
#include <omp.h>
#include <string.h>

#if defined(_MSC_VER)
#include <intrin.h>
#elif defined(__GNUC__) && (defined(__x86_64__) || defined(__i386__))
#include <immintrin.h>
#include <x86intrin.h>
#endif

void add_matrix_naive (double **result, double **mat1, double **mat2, int rows, int cols) 
    int simdCols = cols / 4 * 4;
        if(simdCols > 0)
            for(unsigned int i = 0; i < rows; i++)
                for(unsigned int j = 0; j < simdCols; j += 4)
                    _mm256_storeu_pd(result[i] + j, _mm256_add_pd(
                        _mm256_loadu_pd(mat1[i] + j)
                        , _mm256_loadu_pd(mat2[i] + j)));
                
            
        

        //Handle extra columns
        if(simdCols < cols)
            for(unsigned int i = 0; i < rows; i++) 
                for(unsigned int j = simdCols; j < cols; j++)
                    result[i][j] = mat1[i][j] + mat2[i][j];
                
            
        


void add_matrix(double **result, double **mat1, double **mat2, int rows, int cols) 
    int simdCols = cols / 4 * 4;
    #pragma omp parallel if (rows*cols >= 2000)
    
        if(simdCols > 0)
            #pragma omp for collapse(2)
            for(unsigned int i = 0; i < rows; i++)
                for(unsigned int j = 0; j < simdCols; j += 4)
                    _mm256_storeu_pd(result[i] + j, _mm256_add_pd(
                        _mm256_loadu_pd(mat1[i] + j)
                        , _mm256_loadu_pd(mat2[i] + j)));
                
            
        

        //Handle extra columns
        if(simdCols < cols)
            #pragma omp for collapse(2)
            for(unsigned int i = 0; i < rows; i++) 
                for(unsigned int j = simdCols; j < cols; j++)
                    result[i][j] = mat1[i][j] + mat2[i][j];
                
            
        

    


int main() 
 
    omp_set_num_threads(8);
    //Allocate Matrices
    int rows = 200;
    int cols = 200;

    double **matrix_a = malloc(rows * sizeof(double *) + rows*cols*sizeof(double));

    double * dataStart = (double *) matrix_a + rows; //Offset row pointers
    for(unsigned int i = 0; i < rows; i++)
        matrix_a[i] = dataStart + i * cols;
        memset(matrix_a[i], 0, sizeof(double) * cols);
    

    double **matrix_b = malloc(rows * sizeof(double *) + rows*cols*sizeof(double));

    dataStart = (double *) matrix_b + rows; //Offset row pointers
    for(unsigned int i = 0; i < rows; i++)
        matrix_b[i] = dataStart + i * cols;
        memset(matrix_b[i], 0, sizeof(double) * cols);
    

    double **result = malloc(rows * sizeof(double *) + rows*cols*sizeof(double));

    dataStart = (double *) result + rows; //Offset row pointers
    for(unsigned int i = 0; i < rows; i++)
        result[i] = dataStart + i * cols;
        memset(result[i], 0, sizeof(double) * cols);
    

    //Assign random values to matrices.
    for(int i = 0; i < rows; i++)
        for(int j = 0; j < cols; j++)
            matrix_a[i][j] = rand();
            matrix_b[i][j] = rand();
        
    

    
    int LOOP_COUNT = 4;

    double prevTime = omp_get_wtime();
    for(int i = 0; i < LOOP_COUNT; i++)
        add_matrix(result, matrix_a, matrix_b, rows, cols);
        
    
    double endTime = omp_get_wtime();
    double firstTime = (endTime - prevTime)/LOOP_COUNT;
    printf("Took %f Seconds\n", firstTime);

    //Assign random values to matrices.
    for(int i = 0; i < rows; i++)
        for(int j = 0; j < cols; j++)
            matrix_a[i][j] = rand();
            matrix_b[i][j] = rand();
        
    

    prevTime = omp_get_wtime();
    for(int i = 0; i < LOOP_COUNT; i++)
        add_matrix_naive(result, matrix_a, matrix_b, rows, cols);
    
    endTime = omp_get_wtime();
    double secondTime = (endTime - prevTime)/LOOP_COUNT;
    printf("Took %f Seconds\n", secondTime);
    printf("Naive Time: %f Faster\n", firstTime/secondTime);

我注意到结果似乎完全取决于 LOOP_COUNT。高循环数的并行/SIMD 版本表现得相当好，但低循环数的朴素解决方案往往会做得更好。

【问题讨论】：

如果您在 gcc/clang 上使用 -O3 编译您的幼稚 C 代码，他们可能也能够对其进行矢量化（查看生成的汇编代码）。 “我不允许在网上发布我的代码”翻译为“我对这件事有这个问题”，这意味着我们可能无能为力。我们需要更多细节。我们需要我们可以用来重现问题的代码。但如果没有代码或任何细节描述可谈，这对于未来的读者来说不是一个有用的问题。 @tadman 有道理，我在帖子中添加了代码。呃，为什么要使用指向数组的指针数组，而不是单个高效的二维数组？ A different way to malloc a 2D array?。这将使编译器更难证明或检查没有别名（即没有输出行指向与某些输入行相同的存储）。 【参考方案1】：

CPU 计算事物的速度比它们访问内存的速度要快。

添加双打非常快。您的内核很可能出现内存瓶颈。

还有更多。