2010-09-22 98 views
5

我正在計算使用OpenCL的n維點之間的歐幾里德距離。我得到了兩個n維點的列表,我應該返回一個數組,其中只包含第一個表中每個點到第二個表中每個點的距離。使用OpenCL的累積數組求和

我的做法是做正規多布爾循環(在表1的每一個點{於表2 {每一點...}},然後做使用並聯每對點的計算。

的歐幾里得然後將距離分爲3部分: 1.取點 中的每個尺寸之差2.平方差(仍針對每個尺寸) 3.將所有在2中獲得的值相加。 4.取平方根的值(在本例中省略了此步驟)。

一切都像一個魅力,直到我嘗試累積所有差異的總和(即,執行上述過程的第3步,下面代碼的第49行)。

作爲測試數據,我使用描述符列表,每個列表都有2個點: 描述符列表1:001,002,003,...,127,128; (p1) 129,130​​,131,...,255,256; (p2)

描述符列表2:000,001,002,...,126,127; (p1) 128,129,130​​,...,254,255; (p2)

所以得到的向量應該具有值:128,2064512,2130048,128 現在我得到的隨機數隨每次運行而變化。

我感謝任何幫助或導致我做錯了什麼。希望一切都清楚我在工作場景

#define BLOCK_SIZE 128 

typedef struct 
{ 
    //How large each point is 
    int length; 
    //How many points in every list 
    int num_elements; 
    //Pointer to the elements of the descriptor (stored as a raw array) 
    __global float *elements; 
} DescriptorList; 

__kernel void CompareDescriptors_deb(__global float *C, DescriptorList A, DescriptorList B, int elements, __local float As[BLOCK_SIZE]) 
{ 

    int gpidA = get_global_id(0); 

    int featA = get_local_id(0); 

    //temporary array to store the difference between each dimension of 2 points 
    float dif_acum[BLOCK_SIZE]; 

    //counter to track the iterations of the inner loop 
    int loop = 0; 

    //loop over all descriptors in A 
    for (int i = 0; i < A.num_elements/BLOCK_SIZE; i++){ 

     //take the i-th descriptor. Returns a DescriptorList with just the i-th 
     //descriptor in DescriptorList A 
     DescriptorList tmpA = GetDescriptor(A, i); 

     //copy the current descriptor to local memory. 
     //returns one element of the only descriptor in DescriptorList tmpA 
     //and index featA 
     As[featA] = GetElement(tmpA, 0, featA); 
     //wait for all the threads to finish copying before continuing 
     barrier(CLK_LOCAL_MEM_FENCE); 

     //loop over all the descriptors in B 
     for (int k = 0; k < B.num_elements/BLOCK_SIZE; k++){ 
      //take the difference of both current points 
      dif_acum[featA] = As[featA]-B.elements[k*BLOCK_SIZE + featA]; 
      //wait again 
      barrier(CLK_LOCAL_MEM_FENCE); 
      //square value of the difference in dif_acum and store in C 
      //which is where the results should be stored at the end. 
      C[loop] = 0; 
      C[loop] += dif_acum[featA]*dif_acum[featA]; 
      loop += 1; 
      barrier(CLK_LOCAL_MEM_FENCE); 
     } 
    } 
} 

回答

7

你的問題就出在這行代碼:

C[loop] = 0; 
C[loop] += dif_acum[featA]*dif_acum[featA]; 

所有線程都在您的工作組(當然,實際上所有的線程,但讓我們稍後來談談)試圖同時修改這個數組的位置,而沒有任何同步。有幾個因素使這種真正問題的:

  1. 該工作組不能保證完全並行工作,這意味着對於某些線程C [循環] = 0可後其他線程已執行下一行稱爲
  2. 並行執行的所有程序都從C [loop]中讀取相同的值,用它們的增量修改並嘗試寫回到同一個地址。我不完全確定這個回寫的結果是什麼(我認爲其中一個線程成功寫回,而另一個線程失敗,但我不完全確定),但是它的錯誤方式都是錯誤的。

現在讓我們來解決這個問題: 雖然我們也許能夠得到這個使用原子能,它不會那麼快在全球內存工作,所以讓積聚在本地內存:

local float* accum; 
... 
accum[featA] = dif_acum[featA]*dif_acum[featA]; 
barrier(CLK_LOCAL_MEM_FENCE); 
for(unsigned int i = 1; i < BLOCKSIZE; i *= 2) 
{ 
    if ((featA % (2*i)) == 0) 
     accum[featA] += accum[featA + i]; 
    barrier(CLK_LOCAL_MEM_FENCE); 
} 
if(featA == 0) 
    C[loop] = accum[0]; 

中當然,您可以重用其他本地緩衝區,但我認爲這一點很明顯(順便說一句:您確定dif_acum將在本地內存中創建,因爲我認爲我在某處讀取了不會放入本地內存的地方會使預加載A成爲本地存儲器的一種毫無意義)。

這段代碼的其他一些觀點:

  1. 您的代碼似乎是面向只使用在工作組(不使用任何GROUPID也不全球ID看到工作哪些項目),爲您可能希望使用更多的最佳性能。
  2. 可能是個人偏好研究,但我對我來說似乎更好地使用get_local_size(0)爲workgroupsize比使用定義(因爲你可能在主機代碼改變它沒有意識到你應該改變你的OpenCL代碼)
  3. 代碼中的障礙都是不必要的,因爲沒有線程訪問由另一個線程寫入的本地內存中的元素。因此,您不需要爲此使用本地內存。

考慮到最後一顆子彈,你可以簡單地做:

float As = GetElement(tmpA, 0, featA); 
... 
float dif_acum = As-B.elements[k*BLOCK_SIZE + featA]; 

這將使代碼(不考慮前兩發子彈):

__kernel void CompareDescriptors_deb(__global float *C, DescriptorList A, DescriptorList B, int elements, __local float accum[BLOCK_SIZE]) 
{ 
    int gpidA = get_global_id(0); 
    int featA = get_local_id(0); 
    int loop = 0; 
    for (int i = 0; i < A.num_elements/BLOCK_SIZE; i++){ 
     DescriptorList tmpA = GetDescriptor(A, i); 
     float As = GetElement(tmpA, 0, featA); 
     for (int k = 0; k < B.num_elements/BLOCK_SIZE; k++){ 
      float dif_acum = As-B.elements[k*BLOCK_SIZE + featA]; 

      accum[featA] = dif_acum[featA]*dif_acum[featA]; 
      barrier(CLK_LOCAL_MEM_FENCE); 
      for(unsigned int i = 1; i < BLOCKSIZE; i *= 2) 
      { 
       if ((featA % (2*i)) == 0) 
       accum[featA] += accum[featA + i]; 
       barrier(CLK_LOCAL_MEM_FENCE); 
      } 
      if(featA == 0) 
       C[loop] = accum[0]; 
      barrier(CLK_LOCAL_MEM_FENCE); 

      loop += 1; 
     } 
    } 
} 
+0

我不得不說,首先,我是非常感激與灰熊的答案。我對OpenCL相當陌生,儘管我需要調整他給出的一些示例代碼,但它讓我直接朝着正確的方向發展。我注意到的重要事情(通過試驗和錯誤):無法處理陣列位置的線程需要丟棄; SCAN循環需要稍微調整,即使用輔助緩衝區累積部分結果並檢查邊界條件以找到要添加的項。再一次非常感謝你!我發佈了適用於我的代碼。 – SebastianP 2010-09-24 14:24:01

3

感謝灰熊,我現在有一個工作內核。我需要在灰熊的答案中修改一些東西:

我在例程的開頭添加了一個IF語句,以放棄所有不會引用我使用的數組中有效位置的線程。

if(featA > BLOCK_SIZE){return;} 

當複製第一個描述符本地(共享)內存(胃內到BS),該指數自該函數返回GetElement每次通話只有一個元素被指定(我跳過了,關於我的問題)。

Bs[featA] = GetElement(tmpA, 0, featA); 

然後,SCAN循環需要小的改變,因爲緩衝器每次迭代之後被覆蓋和數據第一個無法控制哪個線程訪問。這就是爲什麼我'回收'dif_acum緩衝區來存儲部分結果,並且以這種方式防止整個循環中的不一致。

dif_acum[featA] = accum[featA]; 

SCAN循環中還有一些邊界控件可以可靠地確定要加在一起的項。

if (featA >= j && next_addend >= 0 && next_addend < BLOCK_SIZE){ 

最後,我認爲在最後的IF語句中包含循環變量增量是有意義的,以便只有一個線程修改它。

if(featA == 0){ 
    C[loop] = accum[BLOCK_SIZE-1]; 
    loop += 1; 
} 

就是這樣。我仍然想知道如何使用group_size來消除BLOCK_SIZE定義,以及是否有更好的策略可以用於線程使用。

,使代碼看起來終於像這樣:

__kernel void CompareDescriptors(__global float *C, DescriptorList A, DescriptorList B, int elements, __local float accum[BLOCK_SIZE], __local float Bs[BLOCK_SIZE]) 
{ 

    int gpidA = get_global_id(0); 
    int featA = get_local_id(0); 

    //global counter to store final differences 
    int loop = 0; 

    //auxiliary buffer to store temporary data 
    local float dif_acum[BLOCK_SIZE]; 

    //discard the threads that are not going to be used. 
    if(featA > BLOCK_SIZE){ 
     return; 
    } 

    //loop over all descriptors in A 
    for (int i = 0; i < A.num_elements/BLOCK_SIZE; i++){ 

     //take the gpidA-th descriptor 
     DescriptorList tmpA = GetDescriptor(A, i); 

     //copy the current descriptor to local memory 
     Bs[featA] = GetElement(tmpA, 0, featA); 

     //loop over all the descriptors in B 
     for (int k = 0; k < B.num_elements/BLOCK_SIZE; k++){ 
      //take the difference of both current descriptors 
      dif_acum[featA] = Bs[featA]-B.elements[k*BLOCK_SIZE + featA]; 

      //square the values in dif_acum 
      accum[featA] = dif_acum[featA]*dif_acum[featA]; 
      barrier(CLK_LOCAL_MEM_FENCE); 

      //copy the values of accum to keep consistency once the scan procedure starts. Mostly important for the first element. Two buffers are necesarry because the scan procedure would override values that are then further read if one buffer is being used instead. 
      dif_acum[featA] = accum[featA]; 

      //Compute the accumulated sum (a.k.a. scan) 
      for(int j = 1; j < BLOCK_SIZE; j *= 2){ 
       int next_addend = featA-(j/2); 
       if (featA >= j && next_addend >= 0 && next_addend < BLOCK_SIZE){ 
        dif_acum[featA] = accum[featA] + accum[next_addend]; 
       } 
       barrier(CLK_LOCAL_MEM_FENCE); 

       //copy As to accum 
       accum[featA] = GetElementArray(dif_acum, BLOCK_SIZE, featA); 
       barrier(CLK_LOCAL_MEM_FENCE); 
      } 

      //tell one of the threads to write the result of the scan in the array containing the results. 
      if(featA == 0){ 
       C[loop] = accum[BLOCK_SIZE-1]; 
       loop += 1; 
      } 
      barrier(CLK_LOCAL_MEM_FENCE); 

     } 
    } 
} 
+0

我仍然認爲你不需要這些本地緩衝區(當然接受爲accum),因爲dif_acum和Bs只能在featA位置訪問,featA是線程的本地ID,因此數組的每個元素都可以被訪問只有一個線程。此外,爲掃描循環使用兩個緩衝區應該不是必須的,因爲一致性是由障礙保證的(迭代由障礙分隔,並且(至少對於我發佈的循環)在每次迭代中,每個元素只能被一個線程 – Grizzly 2010-09-24 15:25:54