Carved Marker

SSE 介紹 [Page 8]

前面介紹的計算都是蠻簡單的計算。現在我們來看一個比較複雜的例子。假設現在有一堆數字,都是單精確度的浮點數(float),而且範圍都在正負 0.5 之間。然後現在要計算它們的自然指數(exp)。當然我們知道 SSE 並沒有這個指令。事實上除了簡單的四則運算(加、減、乘、除)之外,SSE 和計算有關的指令,只有平方根而已。其它複雜的東西像是各種超越函數都是沒有的。

不過,這並不表示要算超越函數就不能用 SSE 了。這只是表示我們得自己算而已。以自然指數來說,最簡單的算法當然是用 Taylor series。事實上有更快的算法,不過我們為了舉例方便,所以就直接用 Taylor series 了。因為我們的數字範圍是在正負 0.5 之間,所以只需要取到第九項,也就是 x 的八次方項,就已經很接近單精確度浮點數所能表示的精確度了。下面是利用 Taylor series 算 exp(x) 的程式片斷:

  • int i;
  • float result = coeff[8] * x;
  • for(i = 7; i >= 2; i--) {
    • result += coeff[i];
    • result *= x;
  • }
  • return (result + 1) * x + 1;

上面的程式片斷中,coeff 陣列需要事先設定成 1/n! 的值。這可以事先算好,直接放在程式裡面,所以不需要花時間。這個程式片斷是利用 Horner 法來計算 Taylor series,以減少計算量。在上例中,每個 X 需要 8 個乘法和 8 個加法,共 16 個計算。

這個程式要改成用 SSE 來寫,並不會很困難。不過,coeff 陣列需要修改成 __m128 的陣列,而且每項存放四個相同的數字(也就是 1/n!)。這樣可以改成下面的程式:

  • int i;
  • __m128 X = _mm_load_ps(data);
  • __m128 result = _mm_mul_ps(coeff_sse[8], X);
  • for(i = 7; i >=2; i--) {
    • result = _mm_add_ps(result, coeff_sse[i]);
    • result = _mm_mul_ps(result, X);
  • }
  • result = _mm_add_ps(result, sse_one);
  • result = _mm_mul_ps(result, X);
  • result = _mm_add_ps(result, sse_one);
  • _mm_store_ps(out, result);

上面的程式片斷和 x87 的版本長得幾乎一樣,當然它是一次對四個浮點數做運算,所以應該會快四倍。下圖是測試結果:

SSE Taylor series test

這些測試結果是對 1,024 個隨機產生的數字進行的。第一個結果是直接用 exp 函式,它裡面是用 x87 的指令像 F2XM1FYL2X 等指令來完成計算的。當然這個函式並沒有範圍的限制,而且它得到的結果也是最精確的。不過,它的速度非常的慢,平均每個數字需要 128 個 cycles 才能算完。第二個結果是用 x87 版的 Taylor series 來計算,也就是第一個程式片斷。它很明顯要快得多,平均每個數字只需要 66 個 cycles,幾乎是快了一倍。不過,其實這裡就產生一個問題:為什麼 16 個計算會需要 66 個 cycles 才能完成?

第三個結果是用 SSE 計算 Taylor series,也就是第二個程式片斷。它比 x87 版本快了近四倍,平均每個數字只需要 17 cycles。可是,因為 SSE 指令一次對四個數字進行運算,所以其實每四個數字共花了 68 cycles。這裡也是同樣的問題:為什麼 16 個計算會需要 68 cycles?這表示平均一個 SSE 乘法指令或加法指令,花了四個 cycles 執行。但是,理論上 Pentium III 不是可以每 cycle 執行一個 SSE 指令嗎?

也許問題是出在上面的程式片斷裡面的迴圈(for(i = 7; i >=2; i--))。所以,也許可以把它展開試試。它會得到上面的第四個結果,實際上也是快了不少,平均每個四個數字花了約 43 cycles。也就是說,平均每個運算花了 2.7 cycles。這樣已經算是相當不錯了。可是,它和「每 cycle 執行一個 SSE 指令」還有一段距離。

事實上,這是因為這些指令都有相依性。在上面的程式片斷中,result 這個變數不但從頭用到尾,而且每個指令都用到這個變數,也會修改這個變數。因此,第二個指令要等到第一個指令執行完後,才能執行。在實際執行時,因為 Pentium III 的 instruction reorder buffer 夠大,所以可以達到一部份程度的平行執行,但是效率還不是非常理想。

要儘可能達到平行執行的效率,就要用另一個方法來展開迴圈。也就是說,不展開內部的迴圈,而是一次對多組資料進行計算。因為多組資料的計算之間並沒有相依性,要平行執行就容易多了。上圖中的第五個結果,就是一次對四組數字進行計算(也就是共 16 個數字),而內部的迴圈並沒有展開。這樣平均每四個數字只需要 34 cycles,也就是平均每個運算花了約 2.1 cycles。比起第四個結果,又要再快了一些。

那如果一次對八組數字進行計算,會再變快嗎?答案是不會,它反而會變慢。為什麼呢?因為 SSE 只有八個浮點數暫存器,所以如果一次對八組數字進行運算,那一定會有些東西無法放到暫存器中,而需要放到記憶體裡面。雖然這些資料幾乎是一定會放在快速的 L1 cache 中,但是它還是會需要一些額外的 load/store 指令。所以反而不會增加效率。不過,我們還是可以把內部的迴圈也展開,因為這個迴圈並不大,所以展開並不會太誇張。把內部的迴圈也展開後,得到的就是第六個結果,也就是平均每四個數字只需要 22.6 cycles,平均每個運算只需要 1.4 cycles!這已經離「每 cycle 執行一個 SSE 指令」不遠了!

還有辦法再加快嗎?恐怕不容易了。因為它已經達到約 2,837MFLOPS(在 Pentium III 1.0B Ghz 上執行)。事實上,這個結果比前面計算向量內積的速度還要快。所以,要再加快應該是很難了。

所以,這裡要討論的重點就是,當計算相當冗長,而且相依性很高的時候,光是展開內部的迴圈是不夠的,應該要考慮也展開外部的迴圈(也就是一次對多組資料進行同樣的計算)。這樣通常比展開內部的迴圈要更有效。當然,如果情況許可的話,也可以連內部的迴圈也展開,通常可以達到更高的效率。

另外,可能會有人問:這樣計算自然指數,精確度有多高呢?和用 exp 函數比較又是如何?在這個例子中,1,024 個隨機數字中,誤差都小於 1.2E-7,也就是單精確度浮點數的 epsilon。這表示有誤差的情形,都只有最後一個位元不同。因此,這應該已經是非常高的精確度了。

[Part 1] [Page 2] [Part 3] [Part 4] [Part 5] [Part 6] [Part 7] [Part 8]

11/29/2001, Ping-Che Chen


Sorry, Traditional Chinese only. This page is encoded in UTF-8.

Copyright© 2000, 2001 Ping-Che Chen