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

Question

我正在尝试使用SIMD指令(_mm256_add_pd、store、load等）优化C中的2D矩阵加法。然而，我并没有看到一个大的加速。使用一些计时代码，我看到的加速范围是.8x-1.5倍的朴素解决方案）。我想知道这是否是典型的？我想这可能是一个内存瓶颈，因为在这种情况下，计算似乎很少。我相信这应该给我的速度提高大约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版本做得相当好，但当循环数较低时，朴素的解决方案往往做得更好。

Answer 1

CPU的计算速度比它们访问内存的速度快得多。

加双打是非常快的。你的核心很可能是内存瓶颈。

还有更多。

但最重要的是，你的基准测试是完全错误的。

>

要比较算法的不同实现，必须编译2个不同的程序。当您在一个程序中运行这两种实现时，会发生许多有趣的事情：缓存层次结构开始出现，CPU有另一个专用缓存用于解码指令，现代CPU在热节流方面速度极快，等等。

编译器并不愚蠢。当他们意识到你在不使用结果的情况下计算了一些东西时，他们可以丢弃代码。当他们意识到你在不改变输入数据的情况下多次计算某件事时，他们就可以消除冗余代码，只计算一次。