矢量化和並行化至強Phi

Question

我正在尋找一個簡單的例子，其中使用至強Phi上的矢量化和並行化，這比只使用Xeon的性能更好。請問你能幫幫我嗎？

我正在嘗試下一個例子。我評論的線14，18和19只至強運行和uncoment這些對於至強披，但只有至強比至強-PHI有關自動向量化

1.void main(){ 
2.double *a, *b, *c; 
3.int i,j,k, ok, n=100; 
4.int nPadded = (n%8 == 0 ? n : n + (8-n%8)); 
5.ok = posix_memalign((void**)&a, 64, n*nPadded*sizeof(double)); 
6.ok = posix_memalign((void**)&b, 64, n*nPadded*sizeof(double)); 
7.ok = posix_memalign((void**)&c, 64, n*nPadded*sizeof(double)); 
8.for(i=0; i<n; i++) 
9.{ 
10. a[i] = (int) rand(); 
11. b[i] = (int) rand(); 
12. c[i] = 0.0; 
13.} 
14.#pragma offload target(mic) in(a,b:length(n*nPadded)) inout(c:length(n*nPadded)) 
15.#pragma omp parallel for 
16.for(i = 0; i < n; i++) 
17. for(k = 0; k < n; k++) 
18.  #pragma vector aligned 
19.  #pragma ivdep 
20.  for(j = 0; j < n; j++){ 
21.    c[i*nPadded+j] = c[i*nPadded+j] + a[i*nPadded+k]*b[k*nPadded+j]   
22.} 

Answer 1

一兩句話更好的性能。自動矢量化的優點是簡單。你需要設置一些關鍵字而不是發生魔法，編譯器會爲你製作快速代碼。如果你想這樣試試這個manual。

這種方法的缺點是沒有簡單的方法來理解編譯器如何工作。在矢量化報告中，您會看到「LOOP WAS VECTORIZED」或「LOOP WAS NOT VECTORIZED」。但是如果你想真正理解你的代碼是如何工作的，唯一的方法就是看你的程序集。這不是組裝的問題。你需要用-fcode-asm來編譯程序。但我認爲如果你需要閱讀程序集來檢查「簡單的自動插件」方法是如何工作的，那麼它並不那麼簡單。

替代自動矢量化的是內在函數（實際上，這不是單一的替代方法）。考慮像使用C函數包裝的組件一樣的內部函數。許多內部函數內部包裝單個裝配命令。

我推薦使用這個intrinsics guide。

所以我簡單的方法步驟：

1. 讓單線程參考實現。您將使用它來檢查內部版本的正確性。
2. 實施SSE內在版本。 SSE內部函數非常簡單，可以在Xeon上進行測試。
3. 實現Xeon Phi的AVX-512版本。
4. 衡量你的速度。

讓我們來做你的程序。 與您的程序有很多不同之處：

1. 我使用float而不是double。
2. 我使用_mm_malloc代替posix_memalign。
3. 我想n被16除以餘數（AVX-512向量寄存器中的16個浮點數）。在這個例子中，我不使用循環剝皮。
4. 我使用本地模式，而不是卸載模式。 KNL是可引導的，因此不再需要使用卸載模式。
5. 另外我認爲你的程序是不正確的，因爲它在一段時間內從多個線程修改c數組。但讓我們認爲這不重要，我們只需要一些計算工作。

我的代碼工作時間：

的Intel Xeon 5680

• 參考計算時間：97.677505秒
• 內在函數計算時間：6。189296秒

的Intel Xeon披（KNC）SE10X

• 參考計算時間：199.0秒
• 內在函數計算時間：2.78秒

代碼：

#include <stdio.h> 
#include <omp.h> 
#include <math.h> 
#include "immintrin.h" 
#include <assert.h> 

#define F_E_Q(X,Y,N) (round((X) * pow(10, N)-(Y) * pow(10, N)) == 0) 

void reference(float* a, float* b, float* c, int n, int nPadded); 
void intrinsics(float* a, float* b, float* c, int n, int nPadded); 

char *test(){ 
    int n=4800; 
    int nPadded = n; 

    assert(n%16 == 0); 

    float* a = (float*) _mm_malloc(sizeof(float)*n*nPadded, 64); 
    float* b = (float*) _mm_malloc(sizeof(float)*n*nPadded, 64); 
    float* cRef = (float*) _mm_malloc(sizeof(float)*n*nPadded, 64); 
    float* c = (float*) _mm_malloc(sizeof(float)*n*nPadded, 64); 
    assert(a != NULL); 
    assert(b != NULL); 
    assert(cRef != NULL); 
    assert(c != NULL); 

    for(int i=0, max = n*nPadded; i<max; i++){ 
    a[i] = (int) rand()/1804289408.0; 
    b[i] = (int) rand()/1804289408.0; 
    cRef[i] = 0.0; 
    c[i] = 0.0; 
    } 
    debug_arr("a", "%f", a, 0, 9, 1); 
    debug_arr("b", "%f", b, 0, 9, 1); 
    debug_arr("cRef", "%f", cRef, 0, 9, 1); 
    debug_arr("c", "%f", c, 0, 9, 1); 

    double t1 = omp_get_wtime(); 
    reference(a, b, cRef, n, nPadded); 
    double t2 = omp_get_wtime(); 
    debug("reference calc time: %f", t2-t1); 

    t1 = omp_get_wtime(); 
    intrinsics(a, b, c, n, nPadded); 
    t2 = omp_get_wtime(); 
    debug("Intrinsics calc time: %f", t2-t1); 

    debug_arr("cRef", "%f", cRef, 0, 9, 1); 
    debug_arr("c", "%f", c, 0, 9, 1); 

    for(int i=0, max = n*nPadded; i<max; i++){ 
    assert(F_E_Q(cRef[i], c[i], 2)); 
    } 

    _mm_free(a);     
    _mm_free(b); 
    _mm_free(cRef); 
    _mm_free(c); 
    return NULL; 
} 

void reference(float* a, float* b, float* c, int n, int nPadded){ 
    for(int i = 0; i < n; i++) 
    for(int k = 0; k < n; k++) 
     for(int j = 0; j < n; j++) 
     c[i*nPadded+j] = c[i*nPadded+j] + a[i*nPadded+k]*b[k*nPadded+j];   
} 

#if __MIC__ 

void intrinsics(float* a, float* b, float* c, int n, int nPadded){ 
    #pragma omp parallel for 
    for(int i = 0; i < n; i++) 
    for(int k = 0; k < n; k++) 
     for(int j = 0; j < n; j+=16){ 
     __m512 aPart = _mm512_extload_ps(a + i*nPadded+k, _MM_UPCONV_PS_NONE, _MM_BROADCAST_1X16, _MM_HINT_NONE); 
     __m512 bPart = _mm512_load_ps(b + k*nPadded+j); 
     __m512 cPart = _mm512_load_ps(c + i*nPadded+j); 
     cPart = _mm512_add_ps(cPart, _mm512_mul_ps(aPart, bPart)); 
     _mm512_store_ps(c + i*nPadded+j, cPart); 
     } 
} 

#else 

void intrinsics(float* a, float* b, float* c, int n, int nPadded){ 
    #pragma omp parallel for 
    for(int i = 0; i < n; i++) 
    for(int k = 0; k < n; k++) 
     for(int j = 0; j < n; j+=4){ 
     __m128 aPart = _mm_load_ps1(a + i*nPadded+k); 
     __m128 bPart = _mm_load_ps(b + k*nPadded+j); 
     __m128 cPart = _mm_load_ps(c + i*nPadded+j); 
     cPart = _mm_add_ps(cPart, _mm_mul_ps(aPart, bPart)); 
     _mm_store_ps(c + i*nPadded+j, cPart); 
     } 
} 

#endif