Volume Compute In SIMT Hardware Architecture

--轉載請註明出處 

 

 

 在進行基於體資料的PDE計算時,總是會涉及到鄰接單元(neighbir cell)的訪問,想要提高計算效率就必須盡量共用鄰域資料,減少訪問全域記憶體的次數.不同於二維的情況,尤其是需要多次迭代計算時,三維紋理的效率很多時候差強人意,且需要
在下一步迭代開始前進行大量資料的複製操作.而如果用二維紋理,cache命中率讓人不敢恭維,且同樣需要記憶體複製操作.雖然線型層次化紋理不需要複製操作,但最好的情況下也只能獲得二維紋理同樣的緩寸獎勵. 這裡介紹基於共用緩寸提高鄰域資料共用率以減少全域記憶體訪問次數的共用記憶體和寄存器迴圈切片法（SRC:Slice of Shared Memory-Refister Cycle Method）. 由於我沒有有效畫圖功能工具，所以雖然不情願而且可能還會影響你對我在這裡所說的話的理解，但是也只能如此，如有疑問可以聯絡我--
QQ:295553381
Email:cyrosly@163.

 

 不多廢話。首先，讓我們在腦中設想一個立方體作為對體資料的抽象，然後將該立方體在三個維方向上均勻的分割，然後我們 將所有分割的單元（cell)影射到平鋪的thread grid上。注意：這裡cell = sub volume.

grid policy: 每個線程塊（thread-block)處理一個cell,所以這裡整個grid的線程的數量並不等於整個體資料的元素數量。
從每個cell 的第一層開始，首先將第一遍的沿Z軸正反方向的鄰接點載入寄存器。然後在後面的迴圈中迴圈交換前中後
處理單元的次序，並更新共用記憶體。 共用記憶體的大小=(blockDimY+2)*(blockDimX+2)。
不同的機器最高效的塊尺寸的選擇不同。在我的機器上16x8可以獲得最好的效能。
雖然如何分配共用緩寸可能有N種策略。當然如果你願意。不過這樣的劃分可以最小化使用量（當然可以用register代
替，但是那樣可並發的thread-block數量就會因為平均可用寄存器的數量增加而減少反而影響效能（實驗過）。
另外也不要疑惑，在這裡（後面的程式）的共用記憶體訪問不會有bank conflicts(所以不用在X維度使用多餘的輔助空
間來避免衝突）。 為了便於計算線程對資料的索引（注意，由於每個thread-block-plane對應的是一個體資料，而不是簡單的平面映射
所以計算索引時要格外小心，你 也可能會因為在大腦中組織空間的邏輯結構時的”糾纏不清“而崩潰o(^!^)o。所以不要抓狂，
你可以簡單的將它設定成對應的體資料的X-Y層大小，而在程式中沿Z軸順序後推，事實上很多時候迴圈體的層數對效能影響
不大，例如256x256x256與512x512x64相同大小的體資料但決定不同的grid劃分策略時，效率幾乎沒有差別,當然也只是有些時候.後面的測試你將會看到差距.
然後的計算中每個迴圈遍中都有一slice共用記憶體和寄存器中的的資料被替換和交換：計算完後當前slice的後置slice
成為下一遍的當前slice，而前置slice將連同後置slice的鄰接單雲的的資料一併更新。並成為後置slice。而當前slice
將成為前置slice.

 

最後還需要注意邊界條件的處理，需要一個單獨的核心來計算邊界單元的值。

演算法描述：
<1>
global/tex
load------------>register:front neghbors
global/tex
load------------>register:middle
global/tex
load------------>register:back neghbors
global/tex
<2>
__LOOP__[ 0 : layers-1 ]
<2|0>
    smem
    store<--------registers
          /registers
    compute
          /shared mem
    sync
<2|1>
    update shared mem slice
    sync
    swap order of slice
<2|2>
    store the result to global memory which computed in step <2|0>

 

 

kernel code:

 

 __constant__ uint gridSizeU;<br />__constant__ uint gridSizeV;<br />__constant__ uint slice;<br />__constant__ uint subvolume;<br />__constant__ uint sublayers;</p><p>extern "C"<br />__global__<br />void kernel_volflo_pressure( float* target,<br /> const float* source,<br /> const float* div,<br /> float coeff )<br />{<br />__shared__ float cell[ DIMY+2 ][ DIMX+2 ];</p><p>const uint tidx=IM( DIMX, blockIdx.x )+threadIdx.x;<br />const uint iU=tidx&( gridSizeU-1 )+1;<br />const uint iV=IM( DIMY, blockIdx.y )+threadIdx.y+1;<br />const uint subdomainID=tidx/gridSizeU;</p><p> uint idx=slice+IM( subdomainID, subvolume )+IM( gridSizeU, iV )+iU;</p><p>float cell010=source[ idx ];<br />float cell100=source[ idx-slice ];<br />float cell001=source[ idx+slice ];</p><p>const uint slotU=threadIdx.x+1;<br />const uint slotV=threadIdx.y+1;<br />const uint cc0=threadIdx.x==0;<br />const uint cc1=threadIdx.y==0;</p><p>for( uint layer=0; layer<sublayers; ++layer )<br />{<br /> cell[ slotV ][ slotU ]=cell010;</p><p>if( cc0 ){<br />cell[ slotV ][ 0 ]=source[ idx-1 ];<br />cell[ slotV ][ DIMX+1 ]=source[ idx+DIMX ];<br />}</p><p>if( cc1 ){<br />cell[ 0 ][ slotU ]=source[ idx-gridSizeU ];<br />cell[ DIMY+1 ][ slotU ]=source[ idx+IM( DIMY, gridSizeU ) ];<br />} __syncthreads();</p><p>target[ idx ]=0.16666667f*( cell[ slotV ][ slotU-1 ]+<br /> cell[ slotV ][ slotU+1 ]+<br /> cell[ slotV-1 ][ slotU ]+<br /> cell[ slotV+1 ][ slotU ]+<br /> cell100+cell001-coeff*div[ idx ] );</p><p>cell100=cell[ slotV ][ slotU ];<br />cell010=cell001; idx+=slice;<br />cell001=source[ idx ];<br />}<br />}<br />

 

 

效能比較：在512x512x64的網格上（我的機器上一次可以分配的顯存的最大尺寸）
配置：
GPU:8800 GTS 320MB
CPU:DualCore Intel Core 2 Duo E6650,2333 MHz(7x333)
效率對比：
CPU:只計算一次大約14秒，
GPU:迭代計算500次大概3秒
Profile片段:
GTS 8800 320MB:

Reference:
<<Acceleration of a 3D Euler Solver Using Commodity Graphics Hardware>>
--Tobias Brandvik and Graham Pullan
--Whittle Laboratory , Department of Engineering,
--University of Cambridge , Cambridge , CB3 0DY,UK
<<Fluid Simulation>>:SIGGRAPH 2007 Course Notes.
--Robert Bridson1:University of British Columbia
--Matthias M¨ uller-Fischer: AGEIA Inc.