1. 程式人生 > >程式效能優化探討(6)——矩陣乘法優化之分塊矩陣

程式效能優化探討(6)——矩陣乘法優化之分塊矩陣

        有一種性格叫做偏執,有一種矩陣優化運算叫做分塊。實話說,也許我這輩子也用不上這種隨牛B但很複雜的演算法,有些版本的教材直接刪除這個內容。但越是這樣我越想不過,因此借寫這篇部落格,把分塊矩陣乘法徹底分析清楚。

         把矩陣乘法進行分塊優化,很奇妙的演算法,假設我們要做11X11矩陣乘法,AXB = C原本是這樣算的:

        比如我們要計算C1.2,就要把上圖的兩條射線,交叉乘加。把A1.0~A1.10和B0.2~B10.2遍歷完,也就是說要呼叫121個元素,還要進行11次乘法和10次加法運算,才能得出答案!如果按照上圖中分出四個粗線塊來進行分批交叉乘加,就能提高效率。

        關鍵是,怎麼理解這種演算法的轉變呢?如何證明兩種演算法是等價的呢?這裡我肯定沒辦法做完整的論證。不過我可以簡單做個定性的判斷。既然矩陣乘法的本質就是交叉乘加,既然最外層都是加法,那麼根據加法結合律,把射線分段進行乘加,如果射線拼起來剛好是原來整個射線,那麼兩種演算法等價了。

                  (1)                    (2)                (3)

        

        像上圖這樣,把計算C1.2這個元素的大的交叉乘加,分成了三個段組。而什麼時候計算哪組射線,這個順序是沒有影響的,只要小的交叉組配對正確(比如第一個橫條對應第一個豎條),再把交叉乘加的結果都加到C1.2的臨時值中,就能正確得出C1.2的最終結論。

        好了,下面把分塊矩陣乘法的完整程式碼一段段貼出來:

vord brck(array A, array B, array C, int n, int bsize) 
{
    int r, c, k, kk, cc;
    double sum;
    int en = bsize * (n/bsize); /* Amount that frts evenly into blocks */
  
    for (r = 0; r < n; r++)
        for (c = 0; c < n; c++)
            C[r][c] = 0.0;




    for (kk = 0; kk < en; kk += bsize) { 
        for (cc = 0; cc < en; cc += bsize) {
            for (r = 0; r < n; r++) {
                for (c = cc; c < cc + bsize; c++) {
                    sum = C[r][c];
                    for (k = kk; k < kk + bsize; k++) {
                         sum += A[r][k]*B[k][c];
                    }
                    C[r][c] = sum;
                }
            }
        }

        好,第一段程式碼截止在這裡。咋一看很複雜,我們來一條條分析。以我的習慣是從多層迴圈的最裡層進行分析。

    sum = C[r][c];
    for (k = kk; k < kk + bsrze; k++) {
        sum += A[r][k]*B[k][c];
    }
    C[r][c] = sum;

        我們先分析kk=0的情況,而bsize是塊的大小4,通過這個,把迴圈限制在某個分塊裡。程式碼中的C[r][c]儲存的是Cr.c臨時值,通過sum不斷增加C[r][c]的累積加法結果。

        假設kk = 0,r = 1、c = 2。那麼程式碼執行的就是下一段:

    for (k = 0; k < bsize; k++) {
        sum += A[1][k]*B[k][2];
    }

        這就相當於在做第一條紅橫線和第一條藍豎線的交叉乘加!能理解到這裡,那麼我們就能理解更多的程式碼了。

        假設kk = 0,cc = 0,r = 0,那麼以下程式碼:


for (c = 0; c <  bsize; c++) {
    sum = C[r][c];
    for (k = 0; k < bsize; k++) {
        sum += A[r][k]*B[k][c];
    }
    C[r][c] = sum;
}

        這相當於在做以下的交叉乘加遍歷:

        

        如果把外層迴圈for (r = 0; r < n; r++) {}加上,那麼他們的迴圈就是:

        

        我們看到r把所有行遍歷完,但由於k被限制在bsize範圍內,因此每行只處理四個列。

        好了,接下來剩下兩個外層迴圈:

    for (kk = 0; kk < en; kk += bsize) { 
        for (cc = 0; cc < en; cc += bsize) {

        首先是cc,cc被限制在en裡,而en的計算是限制在四個小分塊中,好,那麼我們先分析cc。我們看到,cc影響的是B[k][c],因此能把藍橫豎線遍歷的範圍擴大:

        

        我們看到,由於迴圈引數的限制,for (cc = 0; cc < en; cc += bsize) 的迴圈只增加了四條小射線,真正擴大的因素在 for (kk = 0; kk < en; kk += bsize) 。關於這個kk的變化,會同時影響到A[r][k]和B[k][c],而由於en本身是bsize的整數倍結果,那麼kk的遞增的擴充套件結果:

        

        好了,我們發現,第一段程式碼能遍歷的範圍似乎並不全,這時候我發現,在外層迴圈for (kk = 0; kk < en; kk += bsize)裡,還有一段程式碼忘分析了:

/* $end bmm-rck */
/* Now finish off rest of c values  (not shown in textbook) */
    for (r = 0; r < n; r++) {
       for (c = en; c < n; c++) {
           sum = C[r][c];
           for (k = kk; k < kk + bsize; k++) {
                 sum += A[r][k]*B[k][c];
           }
           C[r][c] = sum;
       }
    }
/* $begin bmm-rck */
}

        我們到,這段程式碼的for (r = 0; r < n; r++) 和A[r][k]為遍歷A的行提供了條件,而for (c = en; c < n; c++)和B[k][c]就是遍歷A的列。

        這段程式碼相當於遍歷矩陣B中分塊en以外的那部分列,要是單獨拿來分析,當kk = 0時的遍歷範圍就應該是:

        

        再加上外層迴圈for (kk = 0; kk < en; kk += bsize),遍歷範圍應該是:

        

        經過對迴圈for (kk = 0; kk < en; kk += bsize)內的兩段程式碼分析,我們找出了他們實現的交叉乘加遍歷範圍:

        

        剩餘未納入計算的區域,矩陣A和矩陣B中一目瞭然,接下來分析第三段程式碼:


    /* $end bmm-rck */
    /* Now finish remaining k values (not shown in textbook) */ 
    for (cc = 0; cc < en; cc += bsize) {
        for (r = 0; r < n; r++) {
            for (c = cc; c < cc + bsize; c++) {
                sum = C[r][c];
                for (k = en; k < n; k++) {
                    sum += A[r][k]*B[k][c];
                }
                C[r][c] = sum;
            }
        }
    }

        很明顯,這段程式碼和第一段程式碼的唯一區別就是for (k = en; k < n; k++) ,因此交叉乘加的遍歷範圍很容易得出:

        


        好了,就只剩下B的最後9個元素沒有被交叉乘加過,剩下的這點程式碼,和第三段程式碼唯一區別就是for (c = en; c < n; c++):


    /* Now finish off rest of c values (not shown in textbook) */
    for (r = 0; r < n; r++) {
        for (c = en; c < n; c++) {
           sum = C[r][c];
           for (k = en; k < n; k++) {
               sum += A[r][k]*B[k][c];
           }
           C[r][c] = sum;
        }
   }
    /* $begin bmm-rck */
}
/* $end bmm-rck */

        那麼這段程式碼實現的交叉乘加遍歷範圍是:

        

        至此,我們終於將所有交叉乘加遍歷乾淨!

        費了那麼多功夫,有個問題是,幹嘛我們要用這種分塊來實現演算法優化呢?主要考慮在於,當陣列大小增加時,時間區域性性會明顯降低,快取記憶體中不命中數目增加。當使用分塊技術是,時間區域性性由分塊大小來決定,而不是陣列總大小來決定。另外,分塊雖然能明顯提高效率,但使得程式碼灰常難懂,因此一般用於編譯器或者頻繁執行的庫例程……所以說嘛,反正我這輩子是用不上的……

         那麼我們就具體來實踐一下,把之前的六個版本,和分塊的兩個版本b_rck、b_rkc(我只分析了前一個),結果如下:

rck crk ckr kcr krc rkc b_rck b_rkc
25 0.12 0.02 0.01 0 0 0.31 12.38 14.48
50 0.09 0.01 0.01 0.1 0 0.02 0.01 0.07
75 0.02 0.01 0.02 0.05 0.01 0 0 0.15
100 0.07 0.07 3.71 1.78 4.48 0.11 1.37 1.6
125 0.67 1.88 6.13 4.23 0.67 0.66 0.66 2.25
150 1.06 1.44 9.8 8.16 1.77 1.79 1.78 2.47
175 1.66 1.66 13.05 12.32 1.95 2.05 1.94 2.37
200 1.59 1.97 17.39 15.43 2.21 2.16 1.87 3.09
225 1.97 2.17 20.16 17.38 2.57 2.38 2.22 2.82
250 2.86 2.71 22.07 21 2.6 2.41 2.17 3.16
275 4.45 3.65 22.49 20.3 2.49 2.39 2.5 3.48
300 4.74 5.44 22.78 21.42 2.34 2.42 3.04 3.51
325 6.95 7.13 23.76 21.08 2.74 2.3 2.95 3.22
350 7.75 7.57 23.97 21.55 2.3 2.95 2.69 3.65
375 7.38 8.53 23.75 21.13 2.91 2.64 2.96 3.99
400 8.99 8.95 23.98 23.09 3.07 3.09 2.93 3.92
425 9.83 9.87 24.79 21.97 3.26 2.72

相關推薦

程式效能優化探討6——矩陣乘法優化矩陣

        有一種性格叫做偏執,有一種矩陣優化運算叫做分塊。實話說,也許我這輩子也用不上這種隨牛B但很複雜的演算法,有些版本的教材直接刪除這個內容。但越是這樣我越想不過,因此借寫這篇部落格,把分塊矩陣乘法徹底分析清楚。          把矩陣乘法進行分塊優化,

程式效能優化探討5——快取記憶體、儲存器山與矩陣乘法優化

        這一節內容將綜合(3)和(4),討論快取記憶體相關的程式優化。 一、牛B完了的儲存器山         一個程式從儲存系統中讀資料的速率被稱為讀吞吐量或讀頻寬。如果一個程式在s秒的時間段內讀n個位元組,那麼讀吞吐量就是n/s,一般用MB/s作為單位。  

程式效能優化探討3——儲存器層次結構與快取記憶體

        連外行都大概清楚,目前硬體速度的瓶頸在硬碟而不是CPU。為了有效的克服不同器件之間的速度差,從CPU到硬碟引入了多級快取機制。由於快取影響程式讀取速度,因此是實現優化時必須考慮的內容。 一、儲存器層次結構         快取的思想可能存在於任何有速度差的

程式效能優化探討4——直接對映快取記憶體命中率問題的模擬

        前一節初步介紹了快取記憶體的結構和地址劃分策略,以及快取記憶體“讀”處理規則,這一節從討論“寫”開始。 一、快取記憶體寫的處理         快取處理讀的過程是,根據編號查詢相應的值,如果不命中,就從下一集快取調入新的資料,再根據替換策略(不細數),將新

改善深層神經網路——優化演算法6

目錄 1.Mini-batch gradient descent 前我們介紹的神經網路訓練過程是對所有m個樣本,稱為batch,通過向量化計算方式,同時進行的。如果m很大,例如達到百萬數量級,訓練速度往往會很慢,因為每次迭代都要對所

微信小程式開發常用技巧6——列表上拉載入更多

微信小程式API提供了監聽頁面使用者下拉重新整理事件,但是沒有提供上拉監聽事件,實際開發過程中經常會用到上拉列表,載入更多的需求。那麼如何實現呢? 實現原理:利用onReachBottom監聽頁面滑動到底部,然後執行具體的函式方法,例如請求資料,然後將列表資料

一步一步跟我學習lucene6---lucene索引優化多執行緒建立索引

這兩天工作有點忙,部落格更新不及時,請大家見諒; 前面瞭解到lucene在索引建立的時候一個IndexWriter獲取到一個讀寫鎖,這樣勢在lucene建立大資料量的索引的時候,執行效率低下的問題; 磁碟空間大小,這個直接影響索引的建立,甚至會造成索引寫入提示完成,但是沒

MySQL5.7效能優化系列——SQL語句優化2——子查詢-派生表-檢視--概述

章節內容: 使用Semi-join連線優化子查詢、派生表、檢視 使用Materialization優化子查詢 優化派生表、檢視 使用Exist 策略優化子查詢 概述 in或者any子查詢 MySQL查詢優化器具有不同的策略來評估子查詢。對於IN(

EasyPR源碼剖析6:車牌判斷LBP特征

extend 順序 位置 feature tput ray bpf range str 一、LBP特征 LBP指局部二值模式,英文全稱:Local Binary Pattern,是一種用來描述圖像局部特征的算子,LBP特征具有灰度不變性和旋轉不變性等顯著優點。 原始的LBP

6jmeter的監聽器斷言結果

這個是一個監聽器,目前我是緊跟響應斷言後面的,如圖 直接上圖吧,這裡來個錯誤的示範,我把響應斷言用在了header上。然後結果如圖 這說明斷言失敗了。我現在換過來 沒有fail,這說明斷言成功了 上圖可以把斷言結果放進檔案裡面,OK,這裡我們試試把所有斷言結果寫

企業架構研究總結6——聯邦企業架構FEAF的出現和構成

      美國聯邦政府可以說是企業架構應用的先行者和最大倡導者。通過企業架構的發展歷史我們可以看出,早在上世紀九十年代以來,美國軍方就對這種全域性性的資訊共享的理論開始了研究,並開發出符合其特色企業架構框架理論(DoDAF)。除此之外,在Zachman框架引入到美國聯邦政府

Nginx十五:Nginx併發優化專題3——php-fpm優化

pm = dynamicpm.max_children = 300pm.start_servers = 50pm.min_spare_servers = 50pm.max_spare_servers = 300pm.max_requests = 10000rlimit_fil

OpenGL 三維變換模型檢視矩陣

計算機三維圖形學中,一個基本的任務是如何描述三維空間中一個物體位置的變化,也就是如何 描述物體的運動。通常情況下,物體位置的變化包含三個基本的變化:平移、旋轉和縮放,物體的運動也可以用這三個基本的運動

自定義控制元件6---PorterDuffXfermode圖形過濾器橡皮擦應用

activity_main.xml <RelativeLayout xmlns:android="http://schemas.android.com/apk/res/android" xmlns:tools="http://schema

python學習6:python爬蟲requests和BeautifulSoup的使用

前言: Requests庫跟urllib庫的作用相似,都是根據http協議操作各種訊息和頁面。 都說Requests庫比urllib庫好用,我也沒有體會到好在哪兒。 但是,urllib庫有一點不爽的

操作系統筆記內存管理頁,分段和段頁式

分段式內存管理 筆記 關系 代碼 保護 系統 長度 段頁式內存管理 bit 基本內存管理: 進程占用空間必須連續,導致外部碎片以及附加的compaction 整個進程的swap in 和 swap out十分耗時。 解決:分頁 ->內存空間不必連續,無外部碎片,

為什麼大型矩陣乘法要用

      對於矩陣乘法,我們一般會用三重迴圈來實現,但當矩陣維數相當大時,將矩陣分分塊分割成為近似CPU快取大小,會大大提高計算效率.原因就是直接三重迴圈會導致單個矩陣元素來來回回的從快取進出,而分塊後,進出的資料是以分塊矩陣的大小為單位的.另外,平時我們估算演算法的效率主

POJ 3233 Matrix Power Series(求矩陣冪的和——矩陣快速冪 or 二分遞迴+矩陣快速冪)

Matrix Power Series Time Limit: 3000MS Memory Limit: 131072K Total Submissions: 21451 Accepted:

形態形成場矩陣乘法優化dp

現在 字母 個數 def long string ace 取模 lse 形態形成場(矩陣乘法優化dp) 短信中將會涉及前\(k\)種大寫字母,每個大寫字母都有一個對應的替換式\(Si\),替換式中只會出現大寫字母和數字,比如\(A→BB,B→CC0,C→123\),代表

java程式設計師菜鳥進階 HTTP權威指南 HTTP連線管理及對TCP效能的考慮

分享一下我老師大神的人工智慧教程!零基礎,通俗易懂!http://blog.csdn.net/jiangjunshow 也歡迎大家轉載本篇文章。分享知識,造福人民,實現我們中華民族偉大復興!