為什麼轉置512x512矩陣，會比513x513矩陣慢很多？

來源：evol128 的部落格

謹以此文，紀念剛退休的Professor Sibert以及Professor Goel。你們儘管已年過70，卻還仍然堅持在教導學生，實在令人欽佩。我今天所擁有的程式設計知識，經驗，技巧，很大一部分是從你們那兒學來的。謝謝你們。

問題的出處：Stackoverflow 問答貼

事情的起因是這樣的，先看下面這段程式碼：

很普通的一個求矩陣轉置的程式。但是，當MATSIZE取512和513的時候，出現了非常有意思的結果：

512 平均 2.19ms

513 平均 0.57ms

很讓人驚訝吧，513竟然比512快。更進一步的研究發現，size=512的時候，運算速度會比同數量級的其它數字慢很多很多。這是怎麼一回事呢？

stackoverflow上大牛給的解答非常正確，但是這次，我不想做翻譯了。我從Professor Sibert那裡，從Professor Goel那裡，學到的知識，足夠幫我解決這個問題了，我不是一個人。下面是我的解答：

很容易就聯想到，造成這個問題的原因是CPU cache，我們有很多種方式來儲存cache，具體可以參考這裡。

原作者沒有給出他的CPU型號，但是如今的pc幾乎都是採用的set associative的cache結構，下面我用2-way set associate來做例子，講解一下cache的工作原理。

（圖片取自Professor Sibert的講義，這可是純ascii畫的哦= =）
一個記憶體地址，可以劃分為block，tag，word，byte 4個部分。10bits的block，對應了1024個cache set，記憶體地址的block固定了，就必須儲存在相應的set裡面，這樣可以把查詢cache的事件從O(n)縮短為O(1)。

舉個例子，block是1023（1111111111），你的資料就放在第1023個set裡面。可能有人會覺得奇怪，為什麼block不是取的最前面的10bits，這當然是有道理的，通常在記憶體裡資料都是連續存放的，就是說，同一段程式用的資料，他們前10位幾乎都是一樣的，如果用前10位來定位block，那麼collision的發生率非常高，cache效率非常低下，所以才選了後面的10位來定位block。

當然，每個set裡面有多條記錄，2-way是2條，你得遍歷這兩條記錄，比較前面50位的tag，如果tag一樣，並且Valid bit(V)=1，那麼恭喜你，你的資料在cache裡面，接著就可以通過word和byte來取資料了。

如果遍歷完這兩條記錄，還是沒有找到tag的話，那麼很遺憾，你的資料不在cache裡，得從記憶體裡讀。從記憶體裡獲取相應的資料，然後把它存到對應的cache set裡，如果set裡有空位的話最好，如果沒有的話，用LRU來替換。因為一個set裡只有2條資料，所以實現LRU僅僅需要一個額外bit就可以了，非常高效。

好了，背景知識介紹的差不多了，讓我們回到這個問題上來。為什麼512大小的矩陣，會比其它數字慢那麼多？

讓我們來計算一下，512×512的int矩陣，在記憶體裡是連續存放的。每個cache line是16bytes，對應4個int，所以一個n階矩陣的row可以填充n/4個cache set。假設第一個資料a[0][0]正好對應cache set 0，那麼其中每一個資料a[i][j]對應的cache set是(512i+j)/4%1024=(128i+j/4)%1024。可以看到，前面的係數正好可以整除。很不巧的是，在進行矩陣轉置的運算時，在第2個for迴圈中，我們需要依次訪問每一個row中對應i的值。這樣會造成下面的結果：假設i=0，set(a[0][0])=0, set(a[1][0])=128, set(a[2][0])=256…set(a[7][0])=896,set(a[0][0])=0,後面開始重複了，到a[15][0]的時候剛好填完整個cache的所有128整數倍的set，當讀取a[16][0]的時候，將會發生replace，把a[0][0]從cache裡移除。這樣，當源程式的i=1時，將完全重複i=0的計算過程，每次取資料都需要先從memory讀到cache中來，cache的作用完全沒有體現。

而當size=513的時候，事情就不一樣了，mat[i][j]對應的cache set是(513i+j)/4%1024，前面的係數除不盡了，每遞增4次結果會比size=512時偏差1。例如：set(a[0][0])=0, set(a[1][0])=128, set(a[2][0])=256，set(a[3][0])=384, set(a[4][0])=513…這樣就很微妙的把cache set給錯開了。a[16][0]不在第0行而是第4行，不會覆蓋之前的資料。即使將全部的a[0-15][i]都讀入cache，也不會發生碰撞。之後，由於一個cache有4個word，a[0-15][i+1]，a[0-15][i+2]，a[0-15][i+3]也同時被讀進cache裡了，所以計算i+1,i+2,i+3時，僅僅需要讀對應行的資料就可以了，同一行的資料都是連續的，所以碰撞率很低。這個計算過程很好的利用了cache，如果不考慮其他因素（實際上，這個已經是影響執行時間的最大因素了），理論上我們可以節省75%的執行時間，可以看到，這個理論預測是和提問者給的資料相符合的。

總之，當你的data size是128的整數倍的時候，得特別小心，搞不好cache collision就把你的程式給拖慢了呢

Update 1: 原始碼有邏輯錯誤，這點大家都不要吐槽了，程式碼不是我寫的= =

Update 2：帥哥問我，為什麼可以加速這麼多。這個迴圈包括4次讀cache的操作，2次寫cache的操作，以及0-2次replace操作。每次replace操作會有一次memory read，有可能會有memory write（假設它是write back）。前面的讀寫cache時間和讀寫記憶體相比，幾乎可以忽略，對效率產生顯著影響的是後面的記憶體讀寫。如果cache的hit率高了，那麼記憶體讀寫的次數就少了，程式執行時間是會產生很大影響的

Update 3：當然，具體效果還視乎CPU架構而定，我自己試驗的只有節省25%左右時間

Update 4: 有人提出了用劃分矩陣（把大矩陣分成若干個小矩陣分別計算）的方法來求轉置。劃分矩陣可以解決類似的問題（譬如說求兩個矩陣乘積），但是對解決這個問題沒有任何幫助。因為求轉置的時候，每個資料只用到了一次，沒有重複訪問；即便劃分成更小的矩陣，在cache裡面的位置也沒有發生改變。

Update 5: 據說，Professor Goel只是因病休息幾個學期，沒有退休。。。（原來你還要回來教課！！！）