向量化友好的循環分塊因子選擇算法

柴曉菲，劉 松，屈 彬，王 倩，伍衛國

西安交通大學 電子與信息工程學部，西安 710049

1 引言

在許多計算密集型應用程序中，特別是科學和工程計算應用程序，嵌套循環會耗費大部分的運行時間，成為亟待解決的程序熱點[1]。循環分塊[2-7]是一種應用廣泛的循環優化技術，通過仿射變換對程序的嵌套循環部分進行代碼轉換，一方面增強程序的數據局部性，降低cache失效率；另一方面開發循環代碼的粗粒度并行性，充分利用多核處理器的計算性能。分塊后的循環迭代根據分塊因子大小重置訪存順序，從而減小數據重用距離。因此，分塊因子大小的選擇對循環分塊代碼的性能有著重要的影響。近年來，隨著SIMD擴展部件在微處理器和協處理器中的發展，向量寄存器的位數逐漸增加，使得自動向量化技術在開發嵌套循環的細粒度并行性方面得到有效提高[8-10]。但是，循環分塊的分塊因子大小不僅影響程序的局部性，也影響程序的自動向量化收益。因此，如何選擇合適的分塊因子，在保持程序訪存局部性的同時充分利用向量化收益，對程序性能的提高具有積極的意義。

早期比較成熟的循環分塊因子選擇算法（Tile Size Selection，TSS）主要以改善程序局部性和開發程序粗粒度并行性為優化目標，而近年來提出的TSS算法不僅考慮了以上目標，還將開發程序的細粒度并行性作為一個優化目標。Feld等人[11]通過限定可向量化循環層的分塊因子為向量因子（即向量寄存器可存放的最大操作數數目）的倍數關系來避免非對齊數據訪問，并通過分析程序的內存訪問模式，以最大化利用L1級cache為目標來確定最外層循環的分塊因子。Mehta等人[12]通過對建立的TSS模型進行分析，指出在某些嵌套循環中，循環分塊與向量化存在相互制約的關系，即循環分塊會對循環迭代進行分割而向量化則需要較長的循環迭代，因此對可向量化循環層不進行分塊處理。此后，Mehta等人[13]針對訪存數組規模較大的嵌套循環，提出了一種改進的TTS模型，以充分利用微處理器的數據預取功能。雖然以上這些算法考慮了分塊因子對向量化的影響，但是并沒有量化地分析不同分塊因子對向量化收益的影響。此外，以上算法研究的對象是訪存數組大小為向量因子整數倍的嵌套循環，對于不滿足該條件的嵌套循環（稱為病態規模循環），以上研究并未給出詳細分析。

為了解決這些問題，本文提出了一種針對可向量化嵌套循環的分塊因子選擇算法（VECtorization Tile Size Selection，VEC-TSS），能夠使分塊后的程序同時充分地利用編譯器的自動向量化技術和處理器的多級cache架構資源。該算法的創新點主要有：（1）量化分析了循環分塊后程序的向量化收益；（2）詳細分析了病態規模問題的分塊因子選擇過程。

2 VEC-TSS算法

通常情況下，當嵌套循環的最內層循環沒有跨迭代依賴時，主流編譯器會對該層循環進行自動向量化處理。本文提出的VEC-TSS算法僅適用于具有可向量化循環層的嵌套循環。對于正常問題規模的嵌套循環，只要分塊因子的大小是向量因子的整數倍，就可以避免地址非對齊數據訪問，不會影響向量化收益。但是對于病態規模問題，任何一種分塊因子都會導致向量化時的地址非對齊數據訪問。當嵌套循環為病態規模問題時，相較于其他TSS算法，VEC-TSS算法更能有效地降低向量化收益受分塊因子大小的影響。VEC-TSS算法分為兩個部分：可向量化循環層的分塊因子選擇和其他循環層的分塊因子選擇。

2.1 可向量化循環層的分塊因子選擇

為了有利于編譯器對嵌套循環進行自動向量化，可以借助代碼轉換工具（如PLuTo）將可向量化循環層置換到最內層，本節討論可向量化循環層的分塊因子選擇，即嵌套循環的最內層循環的分塊因子選擇。

程序能夠向量化的前提有兩點：（1）訪問地址連續的數據；（2）訪問地址對齊的數據[14]。經過循環分塊后的程序改變了循環邊界，進而改變了向量寄存器訪問數據的地址對齊性，但是并不會影響分塊內訪問數據的地址連續性，因此本文僅討論循環分塊對地址對齊性的影響。不同大小的分塊因子產生的循環邊界不同，使得分塊后的嵌套循環中非對齊數據的數量也不相同。目前的主流處理器對存在地址非對齊數據訪問的程序進行自動向量化時，會額外采用循環展開、循環剝離、數組填充等技術[15-16]，從而增加了開銷，降低了向量化收益。為了便于接下來的分析，定義可向量化數據塊為內存中連續的向量因子大小的數據，其第一個數據的地址對齊于向量寄存器；定義可向量化數據的數目（記為NUM_VEC）為所有可向量化數據塊中的數據個數?？上蛄炕瘮祿蕉?，向量化收益也越大。為了獲得最大化的向量收益，本文選擇的可向量化循環層分塊因子大小應使得循環體中所有可向量化數據的數目最大。

以圖1為例進行說明，（a）和（b）分別表示矩陣乘法分塊前后的核心代碼。其中，分塊因子為(I,J,K)，參與計算的數據類型為8 Byte的雙精度浮點型。目標機器的向量寄存器大小為256位，則向量寄存器一次能夠存取4個數據，要求存取的第一個數據的內存地址是32 Byte對齊的。假設數組的首內存地址為32 Byte對齊，當問題規模N=13時，數組C[i][j]在內存中的數據存儲情況如圖2所示。圖中每一個方格代表數組的一個元素。以分塊因子(I,J,K)對循環進行分塊，則影響數組C[i][j]的分塊因子是(I,J)，行分塊因子是I，列分塊因子是J。圖2中相同顏色（灰度）的方格表示在相同的分塊中。因為 j循環層可以向量化，此時可向量化循環層分塊因子為J。

圖1 分塊前后矩陣乘法matmul的核心代碼

圖2 內存中數組C的數據布局（分塊因子為6）

觀察圖2，數組的第一行C[0][:]分成三部分，位于三個塊中。第一塊中的C[0][0:5]首先被載入cache中參與計算，其中C[0][0:3]是可向量化數據，而C[0][4:5]大小不是向量寄存器大小，無法向量化，故記可向量化數據個數為4。第二塊中的C[0][6:11]被載入cache中參與計算時，C[0][6]的內存地址非對齊，C[0][6:7]無法向量化，而C[0][8:11]可以向量化，故記可向量化數據個數為4。同理分析數組的其余行，最終得到當 J=6時，NUM_VEC=104。

為了使可向量化循環層的收益最大，要求該層的分塊因子J對應的NUM_VEC最大。當數組大小為N時，遍歷所有的分塊因子 J(1≤J≤N)，并計算相應的NUM_VEC值，則可向量化循環層的最優分塊因子Best_J的計算如公式（1）所示：

2.2 其他循環層的分塊因子選擇

本節將詳細分析影響分塊因子選擇的另外兩個重要因素——程序局部性和并行粒度，并由此建立候選分塊因子的搜索空間，從中確定最優分塊因子。

2.2.1 候選分塊因子的搜索空間

（1）程序局部性

循環分塊技術通過減小塊內數據的重用距離，使得數據在下一次被訪問時沒有因為cache容量被替換出去而得到重用，從而帶來局部性收益。因此本文以數據重用距離（Reuse Distance，RD）來度量程序局部性。RD指同一個數據在兩次重用之間所訪問的其他不同數據的個數。以圖1（b）的代碼為例，分塊因子為(I,J,K)，對于數組B[k][j]來說，塊內循環的最外層i每迭代一次，數組B才能得到重用。數組B的兩次重用之間總共訪問了數組A[i][k]的K個不同元素、數組C[i][j]的J個不同元素以及數組B[k][j]本身的K×J-1個不同元素，故數組B的重用距離RDB=K+J+(K×J-1)。推廣到一般情況，下面給出數據重用距離RD的計算公式。

對于一個n層嵌套循環，分塊因子為T=(T1,T2,…,Tn)，其中Ti代表第i層循環的分塊因子。循環語句S訪問的數組有A1,A2,…,Am，則此嵌套循環中的數據訪問矩陣ACC的結構如公式（2）所示：

其中ai,j代表數組Ai在第 j層循環上的數據訪問情況，如公式（3）所示：

ai,j=1表示數組Ai在第 j層循環上具有數據重用。標記重用循環層的向量C=(C1,C2,…,Cm)表示訪問矩陣ACC每一行中第一次出現1的列號，代表的含義是數組Ai在第Ci層循環上具有數據重用。則數組Ai的重用距離和循環語句S的數據重用距離分別如公式（4）和（5）所示：

由公式（5）可知，循環語句S數據重用距離為語句中所有數組的最大重用距離。

通過以上分析可知，圖1（b）示例的循環語句S中RDB最大。當RDB小于某級cache的大小時，說明數組B[k][j]兩次重用之間所有訪問的數據都可以放入該級cache中，則數組 B[k][j]在該級cache中能夠得到重用。因為RDA和RDC都小于RDB，當數組B[k][j]得到重用時，A[i][k]和C[i][j]也一定得到重用。因此，定義循環語句S的數據重用距離為：RDS=max(RDA,RDB,RDC)。當RDS小于某級cache大小時，除去無法避免的第一次訪問數據冷失效外，語句S中的數據都能夠在該級cache中得到重用，從而獲得局部性收益。

當塊內循環語句S的數據重用距離RDS小于某級cache大小時，說明語句S中的所有數組均能在該級cache中得到重用，因此可以將整個分塊的工作集映射到該級cache，以發掘塊內的局部性收益。以圖1（b）的代碼為例來說明分塊內工作集大小的計算，循環層i、k、j構成一個分塊，塊內包含參與運算的數據總共有數組 A[i][k]的 I×K個元素，數組B[k][j]的K×J個元素和數組C[i][j]的I×J個元素，故一個分塊的工作集大小為：I×K+K×J+I×J。因此，可以建立一個基于程序局部性的分塊因子(I,K)的搜索空間，如公式（6），其中 SizeL2、SizeL2和 SizeL3分別為L1、L2和L3級cache的大小。

考慮到多級cache架構的特點，對于同時存儲數據和指令的cache，根據文獻[5]的處理方式，將可用于數據存儲的容量設為該級cache容量的75%。對于r個核共享的cache，將可用于數據存儲的容量設為該級cache容量的1/r。

（2）并行粒度

當采用多核對程序進行并行計算時，多核之間的負載均衡對并行收益有較大影響。文獻[5]提到當每個核處理的分塊數量大于2時，程序的并行收益能夠得到保證，且多個核之間也是相對負載均衡的。本文通過大量實驗也證實了這一點。

以圖1（b）的程序為例，最外層循環是iT，該層循環迭代的數目為。當使用r個核并行執行程序，granularity表示分配給每個核的負載數目（即分塊數目）時，為了保持良好的并行粒度，則應該滿足公式（7）：

2.2.2 優化目標函數

數據在cache中的重用次數（Reuse Count，RC）越多，局部性收益越大，因此本文采用數據的重用次數來度量嵌套循環的局部性收益。以圖1（b）的循環為例，通過公式（4）可以得到一個分塊內語句S中數組B的重用距離最大，在數組B得到重用的同時，數組A和數組C都得到了重用。定義主導數組為執行語句中重用距離最大的數組，圖1（b）中的主導數組為數組B，則主導數組帶來的局部性收益可以代表整個分塊的局部性收益。由于不同分塊內的數組重用情況相同，因此研究一個分塊內數組B的重用情況，可以代表整個程序的局部性收益。下面依然以圖1（b）的代碼來說明RC的計算過程。

假設公式（6）為L2 cache建立的搜索空間，即RDS

根據以上分析，當RCB取最大時，分塊能夠得到最大局部性收益。

通過以上分析，可以得到選擇最優分塊因子大小的過程。首先針對向量化收益，求出可向量化循環層的分塊因子J。然后通過公式（6）和（7）建立分塊因子(I,K)的搜索空間，并在搜索空間中找到使得公式（8）中RCB最大的分塊因子(I,K)。最后組合得到的(I,K,J)即為最優的分塊因子大小。其中在搜索空間中尋找分塊因子時，為了更好地利用空間局部性，K的取值全部為CLS的倍數，CLS表示一個cache行能包含數組元素的個數。

3 VEC-TSS算法的應用

VEC-TSS算法可以作為一個獨立的分塊因子選擇模塊，嵌入到其他針對嵌套循環優化的工具中，例如PLuTo。PLuTo[8]是一款優化嵌套循環的代碼轉換工具，它本身并不提供TSS策略，只采用默認32的分塊因子大小，但其為研究人員預留了重寫分塊因子的接口。本文將VEC-TSS算法在PLuTo中進行了集成和應用。圖3展示了PLuTo的框架和循環代碼轉換流程，以及VECTSS算法在PLuTo中的模塊位置，如圖中虛線框所示。

圖3 VEC-TSS算法在PLuTo中的應用

此外，相關的實現包含以下三部分。

（1）嵌套循環的可向量化識別。由于VEC-TSS算法的適用對象是可向量化的嵌套循環，因此第一步需要判斷嵌套循環是否具有可向量化循環層。在PLuTo的結構體hyperplane_properties中定義變量vec，用來標記每個循環層是否可以向量化。PLuTo中存在發掘分塊后程序可向量化循環層的函數pluto_pre_vectorize_band（）。借鑒以上函數可以實現函數pluto_pre_vectorize（），用來識別分塊前程序各循環層是否可以向量化。

（2）VEC-TSS算法模塊。如果嵌套循環存在可以向量化的循環層，則調用VEC-TSS算法。PLuTo工具利用Candl模塊和isl模塊，將程序的特征信息提取到多面體模型中，如問題規模、數組個數、數據類型、數組下標及循環層順序等。VEC-TSS模塊把以上信息作為輸入，計算得到最優分塊因子。為了使得計算結果能夠被調用，定義全局結構體vec_tile_size來存儲計算結果。

（3）調用VEC-TSS算法的結果。PluTo預留的重寫分塊因子接口為tile.sizes文件，PluTo通過函數read_tile_sizes（）讀取tile.sizes文件中的各循環層分塊因子。故可以通過修改函數read_tile_sizes（），使得PLuTo自動調用結構體vec_tile_size中的數據來進行代碼分塊。

4 實驗分析

為了驗證VEC-TSS算法的有效性，本文在一臺Intel Xeon E7-4820服務器上進行了實驗。該服務器具有4個8核的處理器，每個核擁有32 KB的L1級和256 KB的L2級cache，每個處理器具有16 MB的共享L3級cache；SIMD單元為256 bit；運行64位服務器版的CentOS 7.1.1503操作系統；采用icc 15.0.2編譯器，編譯選項為O3 parallel openmp。實驗使用的PLuTo版本為pluto-0.11.4。

本文從PLuTo的程序集PLuToBench和多面體編譯優化領域的PolyBench 3.2基準測試集[9]中選取了8個可以向量化的基準測試程序。其中測試程序matmul、dsyrk、dsyr2k、lu和trisolv的操作數為雙精度浮點型，測試程序tmm、corcol和covcol的操作數為單精度浮點型。

4.1 病態問題規模結果與分析

不同問題規模的嵌套循環會有不同的最優分塊因子。本文首先選取了一系列病態規模的程序進行測試，將VEC-TSS算法與目前先進的SICA算法[11]和TTS算法[13]進行比較。為了驗證算法的有效性，排除線程數目的影響，本組測試程序均采用8線程執行。圖4展示了分別采用三種算法的分塊程序相對于未分塊的原始程序的平均加速比。觀察圖4可知，隨著問題規模的擴大，VEC-TSS算法的優勢逐漸明顯。這是由于隨著問題規模的增大，可向量化的數據量也明顯增加，向量化技術帶來的優勢更加顯著，VEC-TSS算法能夠使嵌套循環的向量化收益最大，比其他兩個算法更具有優勢。

圖4 病態問題規模時VEC-TSS算法平均加速比

此外，表1展示了問題規模為3 199時各測試程序的分塊因子。觀察表1可知，各基準測試程序通過三種算法得到的分塊因子類似，是由于這些基準測試程序雖然具有不同的循環結構和循環體，但是數據重用距離大小和嵌套循環中工作集大小卻高度相似。

另外，通過VEC-TSS算法分塊后的測試程序tmm、corcol和covcol的加速比較為突出。這是因為這三個程序的操作數為單精度浮點型，相較于雙精度浮點型來說，向量寄存器中可以存放的數據更多，導致分塊不恰當時，更容易產生數據非對齊問題。VEC-TSS算法能夠更好地解決這個問題，相較于其他算法，性能更加優越。

表1 問題規模為3 199時的分塊大小

4.2 可擴展性結果與分析

循環分塊技術不僅能夠改善程序的局部性，還能夠開發程序的并行性。為了測試VEC-TSS算法對分塊后程序并行可擴展性的影響，本文在問題規模為3 199時，分別采用線程數為1、4、8、16和32測試了分塊后程序的加速比，結果如圖5所示。

圖5 VEC-TSS算法對程序并行可擴展性的影響

觀察圖5可以發現，隨著線程數的增多，分塊后程序的加速比近似表現為線性增長，這說明采用VEC-TSS算法的程序具有良好的并行可擴展性。究其原因有兩方面：第一，VEC-TSS算法考慮到了多核共享L3 cache的情況，將L3 cache的容量均分給多個核，故分塊后程序在多線程并行處理時，線程間不會競爭使用L3 cache；第二，VEC-TSS算法通過并行粒度分析限定最外層循環的分塊因子，使得每個核分配到的循環塊個數大于2，保證了多核之間的負載均衡。此外，測試程序lu和trisolv的并行可擴展性卻不理想，如表2所示。這是因為程序本身具有較好的局部性，在多線程并行執行時，反而容易造成偽共享，增加了數據一致性開銷。

表2 lu和trisolv的多線程加速比測試結果

本文僅采用Intel的CPU進行實驗驗證。但是本文提出的VEC-TSS算法同樣適用于與Intel具有相同SIMD擴展指令集的AMD CPU[17]，以及具有相似SIMD擴展的國產飛騰、申威等處理器[18-19]。

5 結束語

本文針對病態規模的嵌套循環經過循環分塊后會產生許多非對齊數據訪問，從而影響自動向量化收益的問題，提出了一種VEC-TSS算法。該算法對于大規模的嵌套循環具有較明顯的優勢，且采用該算法分塊后的循環具有良好的并行性。目前該算法僅適用于CPU處理器，對于具有向量化處理優勢的GPU并不適用。下一步將會考慮GPU架構的特點，并將該算法擴展到GPU架構中。