1. 程式人生 > >CSI-S2:編寫快取記憶體友好的程式碼

CSI-S2:編寫快取記憶體友好的程式碼

       在CSI-VII一篇中,我們瞭解了儲存器系統的層次結構,並知道了層次結構自上而下使用了快取(cashing)技術,因此我們著重介紹了儲存系統中快取記憶體的工作原理。本篇內容,我們通過分析幾個程式碼例項來分析快取記憶體如何影響程式,並提出如何編寫快取記憶體友好程式碼的方法。

程式碼例項

1. 首先讓我們看第一個方法transpose,這是一個矩陣轉置函式,定義如下:

 Typedef int array[2][2];
 Void transpose(array dst,array src)
{
         Int I,j;
         For(i=0;i<2;i++)
             For(j=0;j<2;j++){
                  Dst[j][i]=src[i][j];
             }
}

通過這個方法我們先看如何分析快取的命中情況,在這之前我們假設程式碼執行的機器具有如下屬性:

         Sizeof(int)=4。

         只有一個L1資料快取記憶體,它是直接對映、直寫、寫分配的,塊大小為8位元組。

         快取大小為16個數據位元組,初始為空。

         Src陣列地址從開始,dSt陣列地址從16開始。

據此我們可以給出記憶體和快取的對映關係:


根據這樣的對映關係,我們可以看出每次對於讀src或者寫Dst陣列,都會載入資料到快取,並且剛好每行只能夠放2個元素。例如當我們讀元素src[0][0]時,會載入src[0][0],src[0][1]到同一快取記憶體行中,這是有其在主存中的地址所決定的。同樣,載入dst[0][0]時,也會將dst[0][1]載入進來。對於兩個陣列中的每個元素,下面的表格給出了其命中情況,m代表不命中,h代表命中。

Dst 陣列

列1

列2

行1

m

m

行2

m

                                  

src 陣列

列1

列2

行1

m

m

行2

m

h

        可以看到,對於dst中的每個元素寫都不命中,這個很容易分析出來,從內部迴圈第一次開始,src[0][0],src[0][1]首先被載入到快取記憶體行0中,而在寫Dst的時候也會從行0載入,而這時候是src的資料,所以不命中,根據寫分配的原則,這時會把dst[0][0],dst[0][1]載入到快取記憶體行0中,待到內部迴圈第二次載入src[0][1]時,由於這時快取行0中的資料變為了dst[0][0],dst[0][1],所以依舊不命中,就這樣依次分析,我們可以得到表格中的結果。

在這裡,如果我們將快取大小調大為32位元組,情況會不一樣?

       當然,如果快取的大小為32時,那麼足夠容納下兩個陣列的資料,記憶體行和快取記憶體的對映會是一一對應的情況,因此,所有的不命中都將會是開始的冷不命中。

Dst 陣列

列1

列2

行1

m

h     

行2

m   

h

src 陣列

列1

列2

行1

m

h

行2

m

h

2. 我們將要看到的第二個例子,問題是我們在白紙上打印出黃方格,在列印過程中,需要設定方格中每個點的CMYK(藍色,紅色,黃色,黑色)值,這裡有幾個演算法,我們試著分析他們在一個具有2048位元組、直接對映、塊大小為32位元組的資料快取記憶體上的效率。

Struct point_color

{

           Int c;

           Int m;

           Int y;

           Int k;

}

Strcut point_color sequare[16][16];

Int I,j;

同樣有下面的假設:

       Sizeof(int)=4;

       Sequare其實儲存地址0

       快取記憶體初始為空。

        唯一的儲存器訪問是對於sequare資料中的元素。變數i和j被放在暫存器中。

方法一:                  

  For(i=0;i<16;i++){
         For(j=0;j<16;j++){
               Sequare[i][j].c=0;
               Sequare[i][j].m=0;
               Sequare[i][j].y=1;
               Sequare[i][j].k=0;
         }
}

         首先,我們要說明,這裡衡量演算法的效率的標準就是在程式中快取記憶體的不命中率,所以我們需要明確在方法中的引用數量,和不命中數量。這裡的引用數量即寫數量=16*16*4=1024次,資料大小是16*16*sizeof(point_color)=16^3=4096顯然大於快取大小,且對於每個快取行都可以容納兩個point_color,即每兩迴圈兩次會有一次不命中,不命中數量為16*16/2=128,不命中率為128/1024=12.5%.

 方法二:                  

For(i=0;i<16;i++){
      For(j=0;j<16;j++){
            Sequare[j][i].c=0;
            Sequare[j][i].m=0;
            Sequare[j][i].y=1;
            Sequare[j][i].k=0;
      }
}

         方法二中改變了引用陣列的索引方式,採用了“列主元”的訪問方式,我們知道sequare[0][0]和sequare[1][0]相差16個point_color結構,快取記憶體總共有2048/32=64組,即64個快取行,每行兩個point_color結構,那麼資料在快取的儲存會是這個樣子:

Sequare[0][0]         Sequare[0][1]

Sequare[0][2]         Sequare[0][3]

 ……..                                     ….

Sequare[1][0]         Sequare[1][1]

Sequare[1][2]         Sequare[1][3]

…..                               ……

可以預見對於每次內部的每次迴圈的,都會有一次不命中發生,所以這時的不命中率為16*16/1024=25%

方法三:

For(i=0;i<16;i++){
        For(j=0;j<16;j++){
             Sequare[i][j].y=1;
        }
}
For(i=0;i<16;i++){
        For(j=0;j<16;j++){
              Sequare[i][j].c=0;
              Sequare[i][j].m=0;
              Sequare[i][j].k=0;
        }
}

         方法三寫操作分兩次執行,但這沒有降低快取的不命中率,在第一次迴圈中依然是每兩次迴圈發生一次不命中的情況,所以有128次不命中,第二個迴圈中,雖然寫的元素增多了,但是對於某個結構來說,一旦不命中,意味著寫的第一個元素寫不命中,而其他元素將會命中,一旦命中,則寫都命中,所以跟第一個迴圈的情況其實是一樣的,也有128次不命中,那麼總得不命中率為(128+128)/1024=25%

通過例二我們學到了分析和預測程式快取行為的能力,這裡我們通過計算程式執行過程中的不命中率來衡量利用快取記憶體的效能。

3.  考慮一個n*n矩陣相稱的問題:C=AB。例如,如果n=2,那麼

                 

        矩陣乘法函式通常需要三個巢狀的迴圈來實現,分別用索引I,j,k來標識。如果我們改變迴圈的次序,對程式碼進行一些其他的小改動,我們就能夠得到矩陣乘法的六個在功能上等價的版本,如下圖所示:

 

         對於這六個不同版本,我們內迴圈訪問的矩陣對來表示每個類,例如ijk和jik是類AB的成員,因為最內層的迴圈使用的是A和B的矩陣。接下來我們對每個類統計和分析每個內迴圈迭代中載入讀和儲存寫的數量,每次迴圈迭代中對A、B和C的引用在快取記憶體中的不命中數量,以及每次迭代快取不命中的總數。為了分析的目的,我們假設:

         <1>每個陣列都是double型別的陣列,sizeof(double)=8

         <2>只有一個告訴快取,塊大小為32位元組

         <3>陣列大小的n很大,以至於矩陣的一行都不能完全裝入快取記憶體行。

在類AB的方法中,內層迴圈每次有兩次載入讀操作,A以步長1讀取資料,因為每個快取行可以儲存4個數組元素,所以每次迭代A不命中0.25次。另外,B以步長n讀取資料,因為n很大,對於B的訪問都會不命中,所以每次迭代的總不命中次數為1.25.

           在類AC的方法中,內層迴圈每次有兩次載入讀操作和一次儲存寫操作,A和C的步長都為n,所以每次迭代的不命中次數都為1,總的不命中次數為2.

類似的我們分析類BC的方法,可以得到下面的資料:

矩陣乘法版本

   類

每次迭代使用的載入次數

每次迭代使用的儲存次數

每次迭代A不命中次數

每次迭代B不命中次數

每次迭代C不命中次數

每次迭代總不命中次數2

Ijk&jik(AB)

2

0

0.25

1.00

0.00

1.25

Jki&kji(AC)

2

1

1.00

0.00

1.00

2.00

Kij&ikj(BC)

2

1

0.00

0.25

0.25

0.50

        通過上面的分析我們明白通過重新排列迴圈可以提高空間區域性性從而使得程式變得快取記憶體友好的,尤其對於步長為1的引用,有更將良好的程式效能。但是,隨著陣列大小的不斷增大,時間區域性性會降低,而快取記憶體中不命中的資料也隨著增加。為了改進這個問題,我們使用一種分塊的技術,來提供其時間區域性性。

         分塊的思想是將一個程式中的資料結構組織成為塊的塊組,使得每個能夠載入到L1快取記憶體中,並在這個塊中進行所需的所有的讀和寫,然後丟掉這個塊,載入下一個塊,直到完成計算。

        因此對於矩陣的乘法的分塊我們可以這麼處理,將每個矩陣劃分成子矩陣,例如n=8,那麼我們可以將每個矩陣劃分成四個4*4的子矩陣。

分塊演算法的思想是將A和C劃分成1*bsize的行條,將B劃分成bsize*bsize的塊。最內層的迴圈對用B的一個塊去乘以A的一個行條,將結果放到C的一個行條中,用B中同一個塊,i迴圈迭代通過A和C的行條。

 

注:這個簡單的版本假設的陣列大小是塊大小的整數倍

分塊的關鍵思想是它載入B的一個塊到告訴快取中,使用它,然後丟棄它,對的引用由很好的空間區域性性,因為是以步長1來訪問每個行條的。它也有很好的時間區域性性,因為是連續bsize次引用整個行條的。對B的引用有好的時間區域性性,因為是連續n次訪問整個bsize*bsize塊的。最後,對C的引用由好的空間區域性性,因為行條的每個元素時連續寫的,但對C的引用沒有好的時間區域性性,因為每個行條都只被訪問一次。

        

編寫快取記憶體友好的程式碼

         具有良好區域性性的程式應該更容易有較低的不命中率,而不命中率較低的程式往往比不命中率較高的程式執行的更快。所以,我們應該嘗試著去編寫快取記憶體友好的程式碼。下面我看下幾種常見的基本方法:

1) 讓最常見的情況執行得快。程式通常把大部分時間花在少量的核心函式上,而這些函式通常把大部分時間花在少量迴圈上。根據Amdahl定律,我們要把注意力集中佔系統大部分時間的核心函式中的迴圈上從而提升系統的整體效能。

2) 在每個迴圈內部使快取不命中數量最小。例如,對於下面的方法

Int sumvec(int v[N])

{

           Int I,sum=0;

           For(i=0;i<N;i++)

                    Sum+=v[i];

           Return sum;

}

我們注意到,對於區域性變數i和sum,迴圈體有良好的時間區域性性。考慮一下對向量v的步長為1的引用。一般而言,如果一個快取記憶體的塊大小為B位元組,那麼一個步長為k的引用模式平均每次迴圈迭代會有min(1,(wordsize*k)/B)次快取不命中。當k=1時,它取最小值,所以步長為1的引用是快取記憶體友好的。

(3)一旦從儲存器器中一個數據物件,就儘可能多的使用它,從而使得你程式中的時間區域性性最大。