从两个阵列的点积测量内存带宽

两个数组的点积

for(int i=0; i<n; i++) {
    sum += x[i]*y[i];
}

不重用数据,所以它应该是一个内存绑定操作。 因此,我应该能够从点积测量存储器带宽。

使用为什么向量化循环没有提高性能的代码我的系统带宽为9.3 GB / s 。 但是,当我尝试使用点积计算带宽时,我使用多线程(我的系统具有四个内核/八个超线程)获得了单线程速率的两倍以及三倍速率。 这对我来说没有意义,因为内存绑定操作不应该受益于多线程。 以下是以下代码的输出:

Xeon E5-1620, GCC 4.9.0, Linux kernel 3.13
dot 1 thread:      1.0 GB, sum 191054.81, time 4.98 s, 21.56 GB/s, 5.39 GFLOPS
dot_avx 1 thread   1.0 GB, sum 191043.33, time 5.16 s, 20.79 GB/s, 5.20 GFLOPS
dot_avx 2 threads: 1.0 GB, sum 191045.34, time 3.44 s, 31.24 GB/s, 7.81 GFLOPS
dot_avx 8 threads: 1.0 GB, sum 191043.34, time 3.26 s, 32.91 GB/s, 8.23 GFLOPS

有人可以向我解释为什么我使用多个线程获得两倍于一个线程的带宽和三倍以上的带宽?

这是我使用的代码:

//g++ -O3 -fopenmp -mavx -ffast-math dot.cpp
#include <stdio.h>
#include <string.h>
#include <stdlib.h>
#include <stdint.h>
#include <x86intrin.h>
#include <omp.h>

extern "C" inline float horizontal_add(__m256 a) {
    __m256 t1 = _mm256_hadd_ps(a,a);
    __m256 t2 = _mm256_hadd_ps(t1,t1);
    __m128 t3 = _mm256_extractf128_ps(t2,1);
    __m128 t4 = _mm_add_ss(_mm256_castps256_ps128(t2),t3);
    return _mm_cvtss_f32(t4);
}

extern "C" float dot_avx(float * __restrict x, float * __restrict y, const int n) {
    x = (float*)__builtin_assume_aligned (x, 32);
    y = (float*)__builtin_assume_aligned (y, 32);
    float sum = 0;
    #pragma omp parallel reduction(+:sum)
    {
        __m256 sum1 = _mm256_setzero_ps();
        __m256 sum2 = _mm256_setzero_ps();
        __m256 sum3 = _mm256_setzero_ps();
        __m256 sum4 = _mm256_setzero_ps();
        __m256 x8, y8;
        #pragma omp for
        for(int i=0; i<n; i+=32) {
            x8 = _mm256_loadu_ps(&x[i]);
            y8 = _mm256_loadu_ps(&y[i]);
            sum1 = _mm256_add_ps(_mm256_mul_ps(x8,y8),sum1);
            x8 = _mm256_loadu_ps(&x[i+8]);
            y8 = _mm256_loadu_ps(&y[i+8]);
            sum2 = _mm256_add_ps(_mm256_mul_ps(x8,y8),sum2);
            x8 = _mm256_loadu_ps(&x[i+16]);
            y8 = _mm256_loadu_ps(&y[i+16]);
            sum3 = _mm256_add_ps(_mm256_mul_ps(x8,y8),sum3);
            x8 = _mm256_loadu_ps(&x[i+24]);
            y8 = _mm256_loadu_ps(&y[i+24]);
            sum4 = _mm256_add_ps(_mm256_mul_ps(x8,y8),sum4);
        }
        sum += horizontal_add(_mm256_add_ps(_mm256_add_ps(sum1,sum2),_mm256_add_ps(sum3,sum4)));
    }
    return sum; 
}

extern "C" float dot(float * __restrict x, float * __restrict y, const int n) {
    x = (float*)__builtin_assume_aligned (x, 32);
    y = (float*)__builtin_assume_aligned (y, 32);
    float sum = 0;
    for(int i=0; i<n; i++) {
        sum += x[i]*y[i];
    }
    return sum;
}

int main(){
    uint64_t LEN = 1 << 27;
    float *x = (float*)_mm_malloc(sizeof(float)*LEN,64);
    float *y = (float*)_mm_malloc(sizeof(float)*LEN,64);
    for(uint64_t i=0; i<LEN; i++) { x[i] = 1.0*rand()/RAND_MAX - 0.5; y[i] = 1.0*rand()/RAND_MAX - 0.5;}

    uint64_t size = 2*sizeof(float)*LEN;

    volatile float sum = 0;
    double dtime, rate, flops;  
    int repeat = 100;

    dtime = omp_get_wtime();
    for(int i=0; i<repeat; i++) sum += dot(x,y,LEN);
    dtime = omp_get_wtime() - dtime;
    rate = 1.0*repeat*size/dtime*1E-9;
    flops = 2.0*repeat*LEN/dtime*1E-9;
    printf("%f GB, sum %f, time %f s, %.2f GB/s, %.2f GFLOPSn", 1.0*size/1024/1024/1024, sum, dtime, rate,flops);

    sum = 0;
    dtime = omp_get_wtime();
    for(int i=0; i<repeat; i++) sum += dot_avx(x,y,LEN);
    dtime = omp_get_wtime() - dtime;
    rate = 1.0*repeat*size/dtime*1E-9;
    flops = 2.0*repeat*LEN/dtime*1E-9;

    printf("%f GB, sum %f, time %f s, %.2f GB/s, %.2f GFLOPSn", 1.0*size/1024/1024/1024, sum, dtime, rate,flops);
}

我刚刚下载,编译并按照Jonathan Dursi的建议运行了STREAM,结果如下:

一个线程

Function      Rate (MB/s)   Avg time     Min time     Max time
Copy:       14292.1657       0.0023       0.0022       0.0023
Scale:      14286.0807       0.0023       0.0022       0.0023
Add:        14724.3906       0.0033       0.0033       0.0033
Triad:      15224.3339       0.0032       0.0032       0.0032

八个线程

Function      Rate (MB/s)   Avg time     Min time     Max time
Copy:       24501.2282       0.0014       0.0013       0.0021
Scale:      23121.0556       0.0014       0.0014       0.0015
Add:        25263.7209       0.0024       0.0019       0.0056
Triad:      25817.7215       0.0020       0.0019       0.0027

这里有几件事情可以归结为:

  • 你必须努力工作才能在内存子系统中获得最后一点性能; 和
  • 不同的基准测量不同的东西。
  • 第一个有助于解释为什么您需要多个线程来饱和可用内存带宽。 内存系统中存在很多并发性,并且它利用了你的CPU代码通常需要一些并发性。 多线程执行帮助的一个主要原因是延迟隐藏 - 当一个线程停顿等待数据到达时,另一个线程可能会利用一些刚刚可用的其他数据。

    在这种情况下,硬件可以在单线程上帮助您 - 因为内存访问非常可预测,硬件可以在需要时提前预取数据,即使使用一个线程也可以提供延迟隐藏的优势; 但预取可以做什么是有限制的。 例如,预取器不会跨越页面边界。 其中大部分规范参考是每个程序员应该知道的有关Ulrich Drepper关于内存的知识,现在已经足够老旧,一些差距已经开始显现(英特尔Sandy Bridge处理器的热芯片概述在这里 - 特别注意集成度更高的内存管理硬件与CPU)。

    至于与memset,mbw或STREAM进行比较的问题,跨基准进行比较总是令人头痛,甚至是声称测量同一事物的基准测试。 特别是,“内存带宽”不是一个单一的数字 - 性能会因操作而有所不同。 mbw和Stream都执行一些版本的复制操作,STREAMs操作在这里拼写出来(直接从网页中获取,所有操作数都是双精度浮点):

    ------------------------------------------------------------------
    name        kernel                  bytes/iter      FLOPS/iter
    ------------------------------------------------------------------
    COPY:       a(i) = b(i)                 16              0
    SCALE:      a(i) = q*b(i)               16              1
    SUM:        a(i) = b(i) + c(i)          24              1
    TRIAD:      a(i) = b(i) + q*c(i)        24              2
    ------------------------------------------------------------------
    

    所以在这些情况下,大约1 / 2-1 / 3的内存操作是写操作(并且对于memset,所有内容都是写操作)。 虽然单个写入可能比读取慢一些,但更大的问题是使用写入来饱和存储器子系统要困难得多,因为当然,您无法执行预取写入的等价操作。 交织读取和写入有助于实现,但基本上所有读取的dot-product示例都将针对在内存带宽上钉住针头的最佳可能情况。

    另外,STREAM基准测试(有意)是完全可移植的,只有一些编译器编译指示可以提供向量化,因此,跳过STREAM基准测试并不一定是警告信号,特别是当您在做两个流式读取时。


    我做了我自己的内存基准代码https://github.com/zboson/bandwidth

    以下是八个主题的最新结果:

    write:    0.5 GB, time 2.96e-01 s, 18.11 GB/s
    copy:       1 GB, time 4.50e-01 s, 23.85 GB/s
    scale:      1 GB, time 4.50e-01 s, 23.85 GB/s
    add:      1.5 GB, time 6.59e-01 s, 24.45 GB/s
    mul:      1.5 GB, time 6.56e-01 s, 24.57 GB/s
    triad:    1.5 GB, time 6.61e-01 s, 24.37 GB/s
    vsum:     0.5 GB, time 1.49e-01 s, 36.09 GB/s, sum -8.986818e+03
    vmul:     0.5 GB, time 9.00e-05 s, 59635.10 GB/s, sum 0.000000e+00
    vmul_sum:   1 GB, time 3.25e-01 s, 33.06 GB/s, sum 1.910421e+04
    

    以下是1个线程的电流结果:

    write:    0.5 GB, time 4.65e-01 s, 11.54 GB/s
    copy:       1 GB, time 7.51e-01 s, 14.30 GB/s
    scale:      1 GB, time 7.45e-01 s, 14.41 GB/s
    add:      1.5 GB, time 1.02e+00 s, 15.80 GB/s
    mul:      1.5 GB, time 1.07e+00 s, 15.08 GB/s
    triad:    1.5 GB, time 1.02e+00 s, 15.76 GB/s
    vsum:     0.5 GB, time 2.78e-01 s, 19.29 GB/s, sum -8.990941e+03
    vmul:     0.5 GB, time 1.15e-05 s, 468719.08 GB/s, sum 0.000000e+00
    vmul_sum:   1 GB, time 5.72e-01 s, 18.78 GB/s, sum 1.910549e+04
    
  • 写入:将常量(3.14159)写入数组。 这应该像memset
  • 复制,缩放,添加和三元组的定义与STREAM中的相同
  • mul: a(i) = b(i) * c(i)
  • vsum: sum += a(i)
  • vmul: sum *= a(i)
  • vmul_sum: sum += a(i)*b(i) //点积
  • 我的结果与STREAM一致。 我获得了vsum的最高带宽。 vmul方法当前不起作用(一旦值为零,它会提前结束)。 使用内在函数并展开稍后将添加的循环,我可以获得稍好的结果(大约10%)。

    链接地址: http://www.djcxy.com/p/36345.html

    上一篇: Measuring memory bandwidth from the dot product of two arrays

    下一篇: optimal in this memcpy implementation?