為什麼轉置一個512×512的矩陣，會比513×513的矩陣慢很多？

http://blog.renren.com/share/236761508/14152672566

在校內上無聊的時候看到這個問題，感覺不錯和大家分享一下。原來自己對這方面一點都不瞭解，這裡面還是有好多不懂，先留文在此後面細品。

問題出處：http://stackoverflow.com/questions/11413855/why-is-transposing-a-matrix-of-512x512-much-slower-than-transposing-a-matrix-of

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

#define SAMPLES 1000 
#define MATSIZE 512 

#include<time.h> 
#include<iostream> 
int mat[MATSIZE][MATSIZE]; 

void transpose() 
{ 
   for(int i =0; i < MATSIZE ; i++) 
   for(int j =0; j < MATSIZE ; j++) 
   { 
       int aux = mat[i][j]; 
       mat[i][j]= mat[j][i]; 
       mat[j][i]= aux; 
   } 
} 

int main() 
{ 
   //initialize matrix 
   for(int i =0; i < MATSIZE ; i++) 
   for(int j =0; j < MATSIZE ; j++) 
       mat[i][j]= i+j; 

   int t = clock(); 
   for(int i =0; i < SAMPLES ; i++) 
       transpose(); 
   int elapsed = clock()- t; 

   std::cout <<"Average for a matrix of "<< MATSIZE <<": "<< elapsed / SAMPLES;

 

 

很普通的一個求矩陣轉置的程式。但是，當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)。

（ 在這裡上面的是Memory address從名字應該知道其實存放的是地址，即指標。對於指標一般都是不放實際內容的，一般是起個指示作用。在這裡Memory address的地區功能和大小是固定的，分別分配了指向block address、Block indx、wd、by的地址空間。

而下面的是The cache，其實存放的是實際內容，應該是和上面對應的，地址對應內容很合理。後面從下面的說明可以分析10bit的block對應1024的cache set,所以Memory address裡面的索引block是聯絡的cache裡面的set，這裡是2-way，就是麼個set中一共2個blk，即blk應該是0~2047。

cache裡面的tag是後面提及到的標識符，後面的Block是存放的實際內容word0、word1等，其長度是任意的，比如16bytes。

即儲存單位還是Block，其可以分成n個word；標誌單位是tag，其起識別作用。

）

舉個例子，block是1023（1111111111），你的資料就放在第1023個set裡面（這裡應該是1-way）。可能有人會覺得奇怪，為什麼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大小的矩陣，會比其它數字慢那麼多？

讓我們來計算一下，512x512的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。可以看到，前面的係數正好可以整除。這樣會造成下面的結果：假設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（512/4*16=2048 個cache set），當讀取a[16][]的時候，將會發生replace，把a[0][0]從cache裡移除。這樣，當來源程式的i=1時，將完全重複i=0的計算過程，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。

總之，當你的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%左右