1. 程式人生 > 其它 >矩陣相乘的strassen演算法_C++加速矩陣乘法的最簡單方法

矩陣相乘的strassen演算法_C++加速矩陣乘法的最簡單方法

技術標籤:矩陣相乘的strassen演算法

假設我們封裝了矩陣類,現在要實現矩陣乘法

為了討論方便,設所有矩陣都是 階矩陣。本文假定矩陣是按行儲存的。 最普通的實現方式如下:
for(int i=0;i<n;++i)
    for(int j=0;j<n;++j)
        for(int k=0;k<n;++k)
            C[i][j]+=A[i][k]*B[k][j];

這種方式,當然是速度不快的方式。對於矩陣乘法的加速,有如下兩種策略:

  • 基於演算法的優化
  • 基於硬體的優化

對於演算法優化,最廣為人知的是Strassen演算法,能達到

的時間複雜度,這甚至還不是漸進時間複雜度意義上最快的演算法。但在實際的庫中,沒有用Strassen演算法實現的,因為常數太大,要想超越樸素的演算法,矩陣的規模必須非常大。本文暫不探討演算法優化。

對於硬體優化,有通用矩陣乘法(GEMM),裡面有各種迴圈展開、基於記憶體佈局等的優化技巧。此外還可以多執行緒(C++11有<thread>庫)併發執行。這裡不說的那麼複雜,就說說最簡單的一個技巧:調換迴圈順序。例如將原來的

順序改成 順序:
for(int i=0;i<n;++i)
    for(int k=0;k<n;++k)
        s=A[i][k];
        for(int j=0;j<n;++j)
            C[i][j]+=s*B[k][j];

對於一千階的矩陣,速度提升5倍左右。下面分析一下為什麼效果這麼顯著。

首先從矩陣的儲存方式說起。一般而言,矩陣有兩種儲存方式:

  • 一維陣列
  • 二維陣列

對於矩陣乘法而言,一維陣列顯然比二維陣列好得多(下文的分析也能看出這一點)。但是如果用矩陣類實現其他功能的話,可能其他功能用二維陣列更方便一些。所以下面對於這兩種實現方式都加以分析。

造成矩陣乘法慢的原因,除了演算法上的

以外,還有 記憶體訪問不連續。這會導致cache命中率不高。所以為了加速,就要儘可能使記憶體訪問連續,即不要跳來跳去。我們定義一個概念: 跳躍數,來衡量訪問的不連續程度。

對於最普通的實現方式(順序:

),它是依次計算
中的每個元素。當計算 中任一個元素時,需要將 對應的行與 對應的列依次相乘相加。之前已假設過,矩陣是按行儲存的,所以在 相應行中不斷向右移動時,記憶體訪問是連續的。但 相應列不斷向下移動時,記憶體訪問是不連續的。計算完 的一個元素時, 相應列中已經間斷地訪問了 次,而 只間斷 次(這一次就是算完後跳轉回本行的開頭),故總共是 次。這樣計算完 中所有 個元素,跳轉了 次。但剛才沒有計數 的跳轉次數,加上以後是 。(注意,在計算完 中每行的最後一個元素時, 是從相應行末尾轉到下一行開頭。如果使用一維陣列實現的話,這是連續地訪問,要減掉這 次。同時, 沒有跳轉次數了,還要減掉 次。因此對於一維陣列,跳轉數是 次)

而如果以順序

實現,它將 中元素一行一行計算。當計算 中任一行的第一個元素時,先訪問 中相應行的第一列元素,和 中第一列的第一行元素,然後 不斷往右挪(不間斷),算完後跳轉到下一行(如果二維陣列則間斷一次,一維陣列不間斷),此時 往右挪一個元素(不間斷)。依次這樣挪動,這樣算完 的這一行元素後,恰好按順序將 遍歷一遍,間斷了 次(一維陣列是 次),且恰好從左往右遍歷了 的相應行,間斷了 次(一維陣列沒有間斷),加起來是 次(一維陣列是 次)。故算完 的所有 行後,跳轉了 次(一維陣列是 次)。剛才沒有算 的跳轉,算上後跳躍數是 次(一維陣列是 次)。

(我上面這塊寫的很亂,找時間修改一下)

由此可見:

  • 順序 的跳轉數漸進地少於順序 的跳轉數
  • 一維陣列比二維陣列好

下面是各個迴圈順序的跳躍數列表(下面是我寫文章時現計算的,可能會因為粗心犯錯)

  • 順序 —— (二維陣列)—— (一維陣列)
  • 順序 —— (二維陣列)—— (一維陣列)
  • 順序 —— (二維陣列)—— (一維陣列)
  • 順序 —— (二維陣列)—— (一維陣列)
  • 順序 —— (二維陣列)—— (一維陣列)
  • 順序 —— (二維陣列)—— (一維陣列)

因此從速度來說:

實測速度:(

階, 優化)

發現基本符合(

略微違反,其實它們很接近,回過頭看理論分析也表明很接近)