最大公因數範例
這是 Euclidean 演算法的 C 語言實作版本,可以找出兩個整數的最大公因數(greatest common divisor):
void gcd(int* p, int* q, int* r) { int a = *p; int b = *q; while (a != b) { if (a > b) a = a-b; else b = b-a; } *r = a; }
這裡我們以指標(pointer)的方式來傳遞參數,如果是一般的 C 語言程式是可以不需要這麼麻煩,直接以傳值呼叫即可,這裡這樣寫是為了後續我們要把這個函數改為平行化的向量版本。
若要讓這段城市可以在 GPU 上面執行,我們要使用 QPULib 的語法改寫這一個函數,以下是改寫後的版本:
#include <QPULib.h> void gcd(Ptr<Int> p, Ptr<Int> q, Ptr<Int> r) { // 取讀輸入向量,長度為 16 Int a = *p; Int b = *q; // 平行檢查 a 向量與 b 向量, // 若 a 與 b 中有任何元素不相等, // 則繼續執行 While 迴圈 While (any(a != b)) // 個別比較 a 向量與 b 向量中的元素, // 針對 a 大於 b 的那些元素進行運算 Where (a > b) a = a - b; End // Where 條件判斷式結尾 // 個別比較 a 向量與 b 向量中的元素, // 針對 a 小於 b 的那些元素進行運算 Where (a < b) b = b - a; End // Where 條件判斷式結尾 End // While 迴圈結尾 // 傳回計算結果 *r = a; }
QPULib 是以 C++ 物件導向的方式建立自己的程式語言語法,所以表面上看起來是 C++ 的程式,事實上在使用起來有很大的差異,用這種語法所撰寫的程式碼是放在 QPU 中執行的,所有的運算都是以長度為 16 的向量為基本單位,以平行化的方式在執行。
Int
是一個 QPULib 所定義的類別,代表長度為 16 的 32 位元整數向量。Ptr<Int>
則是代表Int
向量的記憶體位址向量。*p
代表記憶體中以p
開頭的Int
向量(跟 C 語言的陣列相同,第一個向量元素位址為p
)。a != b
代表兩個向量個別元素的比較,結果是長度為 16 的布林向量。any(a != b)
是判斷a != b
這個布林向量中,是否有任何值為true
。Where (a > b)
與End
兩個是一組的條件判斷式,向量中只有符合條件判斷的對應元素會執行情中的程式碼。
以上是 GPU 部分的核心函數,以下是 CPU 部分的 C/C++ 程式碼:
int main() { // 將 gcd 函數編譯成 QPU 核心函數 auto k = compile(gcd); // 配置 CPU 與 QPU 共享的陣列 SharedArray<int> a(16), b(16), r(16); // 設定亂數種子 srand(0); // 以亂數產生測試用資料 for (int i = 0; i < 16; i++) { a[i] = 100 + rand()%100; b[i] = 100 + rand()%100; } // 設定要使用的 QPU 數量 k.setNumQPUs(1); // 使用 QPU 進行平行運算 k(&a, &b, &r); // 輸出結果 for (int i = 0; i < 16; i++) printf("gcd(%i, %i) = %i\n", a[i], b[i], r[i]); return 0; }
這個範例就是前面我們在安裝 QPULib 函式庫測試時所用的 GCD
程式,編譯與執行方式請參考上面的說明。
回圈展開
回圈展開(loop unrolling)是程式碼最佳化常用的技巧,可以有效降低分支指令(branch instructions)的執行次數,加速整體程式的速度。
QPU 的分支指令會需要 3 個延遲間隙(delay slots),相當於 12 個時脈週期,目前 QPULib 在這個間隙中並沒有做其他用途(也就是浪費掉了),所以若要加速程式的執行,必須儘可能減少條件判斷出現的頻率。
QPULib 允許程式設計者以一般 C 語言的 for
迴圈,對特定的程式碼做回圈展開,例如:
#include <QPULib.h> void gcd(Ptr<Int> p, Ptr<Int> q, Ptr<Int> r) { Int a = *p; Int b = *q; While (any(a != b)) for (int i = 0; i < 32; i++) { // 迴圈展開 Where (a > b) a = a - b; End Where (a < b) b = b - a; End } // for 迴圈展開結尾 End *r = a; }
這裡加入的 for
迴圈並不是真正的迴圈,而是做迴圈展開,也就是讓編譯器產生多個重複的 QPU 指令,改善程式執行速度。
迴圈範例
這是一個旋轉座標的 C 語言函數:
void rot3D(int n, float cosTheta, float sinTheta, float* x, float* y) { for (int i = 0; i < n; i++) { float xOld = x[i]; float yOld = y[i]; x[i] = xOld * cosTheta - yOld * sinTheta; y[i] = yOld * cosTheta + xOld * sinTheta; } }
若以 QPULib 改寫後,會變成這樣:
void rot3D(Int n, Float cosTheta, Float sinTheta, Ptr<Float> x, Ptr<Float> y) { For (Int i = 0, i < n, i = i + 16) Float xOld = x[i]; Float yOld = y[i]; x[i] = xOld * cosTheta - yOld * sinTheta; y[i] = yOld * cosTheta + xOld * sinTheta; End }
Non-Blocking 記憶體存取
上面的程式碼只是最初淺的 GPU 版本,這段程式碼在執行時會因為記憶體搬移資料,造成程式產生多餘的等待時間,若要改善效能可以將計算與記憶體存取的時間區間重疊,讓程式一邊計算、一邊存取記憶體,以下是修改過後的程式碼:
void rot3D(Int n, Float cosTheta, Float sinTheta, Ptr<Float> x, Ptr<Float> y) { // index() 函數可產生 0, 1, 2, ..., 15 的向量 Ptr<Float> p = x + index(); Ptr<Float> q = y + index(); // 預先載入資料 gather(p); gather(q); Float xOld, yOld; For (Int i = 0, i < n, i = i+16) // 預先載入下一次迭代要用的資料 gather(p+16); gather(q+16); // 取得當次迭代要用的資料 receive(xOld); receive(yOld); // 將計算結果儲存回記憶體中 store(xOld * cosTheta - yOld * sinTheta, p); store(yOld * cosTheta + xOld * sinTheta, q); p = p+16; q = q+16; End }
QPULib 提供了 non-blocking 的記憶體載入方式:
gather(p)
可以預先載入p
中的每一個值。receive(x)
可以將gather
載入的值存入x
中,程式在這裡會等待所有的直都載入至x
中之後,才會繼續執行(也就是會 block 住)。
gather(p)
的語法與 x = *p
不同,它是針對每一個記憶體位址操作,所以要使用 index()
做一些調整。
記憶體的資料會在 gather
與 receive
之間進行載入,所以這段時間我們可以加入其他的計算,讓計算與資料載入同時進行。
QPU 中有 8 個 FIFO 緩衝區可以儲存 gather
載入的資料,而資料在載入 FIFO 之後,再以 receive
取出,也就是說同時間可以有 8 組 gather
與 receive
一起運作。
另外一個 non-blocking 操作是使用 store
將資料寫回記憶體,store(x, p)
可以將 x
向量寫入 p
這個記憶體位置,與 *p = x
不同的是,store(x, p)
是 non-blocking 的操作,不過若遇到下一個 store
時,就要等待前一個 store
執行完成才能接續執行。
多 QPU 平行運算
QPULib 可以透過 setNumQPUs(n)
函數來設定 QPU 數量,同時使用多顆 QPU 來做平行計算,而在使用多顆 QPU 時,程式設計者必須要能辨識其所使用的 GPU,這樣才能將計算分配到每一顆 QPU 中。
QPULib 的 me()
函數可以傳回 QPU 的 id(型態為 Int
),而 numQPUs()
則會傳回可以使用的 QPU 總數,QPU 的 id 值會介於 0
到 numQPUs - 1
之間。
以下是使用多顆 QPU 平行計算的版本:
void rot3D(Int n, Float cosTheta, Float sinTheta, Ptr<Float> x, Ptr<Float> y) { Int inc = numQPUs() << 4; // numQPUs() * 16 Ptr<Float> p = x + index() + (me() << 4); Ptr<Float> q = y + index() + (me() << 4); gather(p); gather(q); Float xOld, yOld; For (Int i = 0, i < n, i = i + inc) gather(p+inc); gather(q+inc); receive(xOld); receive(yOld); store(xOld * cosTheta - yOld * sinTheta, p); store(yOld * cosTheta + xOld * sinTheta, q); p = p + inc; q = q + inc; End }
以上就是簡單的 QPULib 介紹與教學,更詳細的語法與範例請參考 QPULib GibHub 網站上的說明文件以及其原始碼。
繼續閱讀: 12