計算流體力學程序單核指令級優化方法

劉闖,何峰,肖兮,董小社,張興軍

(西安交通大學計算機科學與技術系,710049,西安)

通過計算流體力學程序模擬真實情景是提高流體機械的制造精度和節約成本的重要手段,計算流體力學程序具有計算密集、處理數據量大的特性,因此常需要大量的計算資源。盡管計算機的計算能力有了飛速的發展,計算流體力學程序對計算性能的利用率仍然存在不充分的現象,導致計算資源的浪費。

目前對于流體計算程序性能提升的研究主要體現在并行化或異構加速上[1-3],例如采用多核系統并行處理或采用協處理器加速計算。超級計算機上部署的大型計算流體力學程序大多都采用了上面兩種計算方式,但是作為一個重要的影響因素,程序的單核執行效率并沒有得到足夠的重視,很多程序的單核性能只能達到系統理論峰值的5%～10%[4],導致大量的計算資源被浪費,因此對程序單核性能的優化應引起重視。

針對流體計算程序的優化,一些文獻提出了綜述性的優化方案和思路[5-7],或針對具體程序進行了優化[8],但是這種優化主要集中于對云計算、網格計算等并行化計算的優化,對單核級優化的工作不夠深入;針對單核級的程序優化的研究中,主要優化方案是流水優化或訪存優化,優化步驟一般是首先分析程序特性或測試程序熱點以找到瓶頸,然后通過調整代碼結構流水性能[9]或調整數據的組織方式優化訪存性能[4,10-12];從處理器的角度,賀愛香、Oyarzun等根據ARM處理器的計算特性對問題提出優化方案[13-14],這些研究具有一定的通用價值,但是對于大型計算程序常用的Intel處理器并沒有深入研究;申小偉等通過對編譯器等運行環境的優化來提升程序性能[15],但是這種方法對于跨專業的程序優化人員具有一定的門檻,且在公用的計算環境中難以部署。根據上述背景,本文針對軸流壓氣機模擬的流體計算程序,提出一種針對計算流體力學程序單核的優化方法。總體優化思路有3個方面:①優化指令結構,減少由相關引起的指令流水斷流,使程序充分發揮CPU的計算能力;②調整存儲結構和訪存結構并進行容器優化,從而增強存儲和訪存的性能;③減少冗余指令,優化底層代碼數量和結構從而減少程序執行所需的時間空間。本文在軸流壓氣機模擬程序上實現優化方法,在TIANHE-1A超級計算機和商用服務器上進行測試實驗,實驗結果表明,本文提出的優化方法可以有效地提高程序的緩存命中率、優化流水連續性、減少目標代碼數,從而大幅提升程序單核執行效率。

1 流體計算程序簡介

1.1 流體計算程序背景

本文采用軸流壓氣機模擬程序[16-17]對軸流壓氣機轉子進行模擬,程序采用有限體積法進行空間離散,3步Runge-Kutta法進行時間推進,控制方程組采用笛卡爾坐標系下的相對流動控制方程組,形式如下

(1)

湍流模型采用SA湍流模型,為加速收斂程序采用三層網格,通過在不同網格層上輪流迭代從而使不同頻率的誤差分量得到有效地衰減。

程序模擬36葉片壓縮機工作時的流體狀態,采用36個進程,每個進程對一個葉片進行模擬,單葉片網格數約為42萬,總計約1 535萬。在邊界處理時需進行通信,通信采用MPI完成。程序主要開銷是在單個進程的計算上,通信等待時間相對整個程序時間占比較少,且優化不涉及通信時間的優化,因此并行開銷可忽略不計。

1.2 流體計算程序實現細節及性能分析

考慮程序的緩存命中率對程序的存儲和計算特性進行分析。存儲方面,程序中網格的存儲采用STL庫中的Vector容器進行存儲,采用該容器的主要原因是Vector容器是順序存儲容器,有較好的隨機存取性能,以及可以利用STL中提供的豐富的操作,功能強大。根據物理網格模型,容器的維度為四維,存儲維度順序為x、y、z、n,其中x、y、z為空間維度,n為多層網格的層號。

計算方面,程序計算主體的運行時間約占整個運行時間的95%以上,剩余時間做程序流程控制、數據的通信與保存等操作,因此程序優化的主體在計算部分。在計算部分中,程序對整個網格層采用for循環進行更新,大多數循環深度為3層,由于邊界上存在虛擬網格,循環次數約為43.5萬次,略多于網格數,每次迭代中更新一個或多個變量。此外,為實現復雜計算,程序還調用了math庫中的部分計算函數,考慮到計算函數計算開銷較大,在某些調用中可采用等價短指令替換。

為進一步分析程序計算性能,本文采用性能分析軟件對程序進行分析。結果顯示,程序的流水線性能未充分發揮,雙向中斷率均超過50%,暫停周期數較高,導致執行某個程序的指令平均時鐘周期(CPI)與理論值有差距。因此,在優化時應考慮流水優化;緩存命中率約為50%,訪存上也有優化空間;此外程序執行指令數為38 534×109條,在優化中需尋找潛在的指令優化來減少冗余代碼。

2 優化策略研究

本文所述方法應用在TIANHE-1A超級計算機系統中Intel Xeon X5670 CPU上,采用Intel CPU作為示例是考慮到大多數的高性能計算環境都采用Intel處理器,且Intel CPU中采用的加速部件在其他處理器中有類似部件設計,因此以Intel CPU為基礎設計的優化方案具有較好的普適性。

圖1 Intel CPU架構示意圖

圖1給出了Intel CPU的架構示意圖,該示意圖基于Intel Xeon X5670 CPU所采用的Nehalem架構,主要體現CPU針對單核程序的加速方式,對其余部件做了省略和簡化。圖1中CPU包含兩個核,每個核內設置獨立的計算部件,從邏輯上分為取值部件、解碼部件、順序回退部件、計算單元,圖1中的箭頭代表了部件間的數據通路。基于上述架構,Nehalem CPU架構部署了取指(IF)、譯碼(ID)、執行(EXE)、寫回(WB)四階段流水機制。此外,Nehalem架構還采用了超標量技術,在CPU中建立了4條流水線,因此理論上可以同時執行4條流水指令。每個核設置獨立的L1和L2緩存,多個核間共享L3緩存,方便進行數據緩存,提升訪存性能。其中一級緩存分為數據緩存和指令緩存,每個緩存容量為32 KB,采用8×64組相連的映射方式,訪存速度為4個時鐘周期;L2緩存容量為256 KB,采用8×512組相連映射方式,訪存速度為10個時鐘周期;L3緩存容量為8 MB,采用16×8 192組相連映射方式,訪存速度為40～75個時鐘周期。

基于對Nehalem架構的分析,程序優化的思路主要集中在使程序充分發揮流水性能和訪存性能上。

除針對CPU架構的優化外,還需考慮程序本身的簡潔性。對于不同形式的程序,通過編譯器編譯得到的目標代碼有很大差距,若有冗余和過多同義代碼會導致性能下降。為保證程序性能,在編碼時需考慮高級語言被編譯后得到的目標代碼是否簡潔高效。

結合上述分析,本文提出的流體計算程序單核優化方案主要考慮3個方面:①優化存儲結構和訪存結構,發揮高速緩存效率;②調整程序的指令結構,發揮流水性能;③優化程序指令,使程序底層指令數得到削減和優化。前兩個方面面向CPU對單核指令執行的高速訪存和流水優化部件展開,目的是充分發揮CPU加速部件的加速效果。第3個方面針對程序代碼本身優化,考慮到不同指令在底層執行時采用不同的匯編指令翻譯,程序的目標代碼應足夠精簡高效。

2.1 訪存優化

訪存優化的目的是使程序在數據訪問時緩存的利用率更高,主要優化手段是調整數據的存儲結構和訪存結構,首先考慮存儲結構。

(a)訪存在內存上連續

(b)訪存在內存上離散

程序的訪存方式是通過循環遍歷的方式對網格進行訪問,優先訪問網格層號n,其次為x、y、z軸,這種存儲結構和訪存方式的不匹配造成了緩存效率較低。圖2給出了不同存儲方式下訪存情況的對比,對于M×N矩陣,規定訪存順序為從(0,1)到(0,N),若矩陣采取aM,N的形式存儲,則矩陣在內存中的物理分布如圖2a所示,在訪存時數據分布在連續的內存中,緩存會將連續內存存入,在下次訪存時可直接取用緩存,減少內存訪問。若矩陣采取aN,M的形式存儲,則矩陣在內存中分布如圖2b所示,在訪存時會產生跳躍,因此緩存無法命中,需要進行內存訪問,影響程序運行效率。

為解決訪存問題,將數據的存儲結構調整為an,z,y,x,這樣調整的好處有:①存儲結構更改后x網格內相鄰各點在內存中存儲的更近,更有利于程序的連續訪問;②對于多維容器來說,由于x>y>z>n,更改后的存儲形式減少了容器在較低維度所創建的數量,減少了目標代碼的指令數。例如,對于容器a20,10,需先創建一個容量為20的容器,每個容器中再創建一個容量為10的容器,總共創建21個,若改為a10,20,則僅需創建11個容器。

下面根據存儲結構調整訪存結構,訪存順序的控制由for循環語句進行。為使訪存順序更適合存儲結構,將循環結構和存儲維度調整為一一對應的關系。即對于an,z,y,x形式的數據,處理循環為最外層循環處理n、次外層處理z、內層處理y、最內層處理x。對不滿足這種處理方式的操作,考察操作的相關性,若無相關性可通過拆分循環來優化訪存順序。

訪存方式的優化和存儲結構的優化有類似的優點:調整訪存結構有利于高速緩存;調整后循環的控制變量初始化次數減少,消減底層指令數。

2.2 流水優化

考慮影響流水線性能的問題,在程序一級影響因素主要有:由分支轉移引起的控制相關會導致流水失效;由寫后讀引起的數據相關,這些相關會導致流水暫停,使流水性能下降。Nehalem CPU架構部署了四階段流水機制,雖然該架構加入了分支預測和數據旁路的方式減少流水線中斷,但是仍不足以完全消除相關性引起的性能損失,因此在優化時仍可消除相關性。

考慮計算中的循環部分,由于循環在每次迭代時需進行一次判斷,若循環內指令較少會導致計算和分支的比例過低,產生控制相關導致性能下降,同時循環中對同一變量的反復讀寫可能引起數據相關。針對控制相關和數據相關問題,一般采用循環展開并引入多中間變量的方法進行優化,循環展開可以在以下方面得到優化:①減少分支數;②減少循環控制的指令開銷;③減少循環內計算的數據相關。循環展開示例如下

for(i=0;i

s+=a[i];

改為:

for(i=0;i

s1+=a[i];

s2+=a[i+1];

}

s=s1+s2;

展開因子k=2,即每次循環計算次數為兩次。兩種代碼的流水情況如圖3所示,流水性能得到了提升。

(a)循環拆分前

(b)循環拆分后

理論上,在展開足夠多循環的情況下可以完全消除流水中斷,展開因子應大于硬件標量數與執行周期的乘積,如對于具有3個浮點乘器件,每個浮點乘需5個時鐘周期的指令集,展開因子應大于15;但是實際操作中需考慮寄存器溢出,即中間變量過多導致寄存器不夠,觸發訪存。為此還需保證最優循環因子k滿足

max{k};F(k)≤r,k∈Z

(3)

式中F(k)為k次展開需要的寄存器數。

此外,解除數據相關還可以采用長計算拆分為短計算、創建臨時變量的方式解除相關、多個連續的寫后讀相關語句交叉等方法。

對if語句引起的控制相關,分支預測部件往往不能正確預測走向,預測錯誤的懲罰往往達到十幾個到幾十個時鐘周期,帶來流水性能損失。在X86指令集中,條件數據傳輸指令可以被編譯為普通指令流水線的一部分,在判斷失敗時沒有懲罰,因此采用條件數據傳送代替條件控制轉移能提升程序性能。常用的條件數據傳輸指令是三目運算符,條件傳輸形式如s=a

x=a

x′=b

boolT=(a

if (!T)x=x′

可以看出,這種編譯方式將分支轉移操作放在偽代碼最后一行的if語句中,在目標代碼中這個if語句通過條件拷貝指令cmov實現,分支的結果只在這一條指令中體現,與上下文無關,從而消除了分支開銷。通過這種轉換可以提升分支在流水線中的穩定性,減少預測錯誤的懲罰。

2.3 指令優化

指令優化的目的是從指令一級減少冗余代碼。冗余代碼主要指程序編譯后的目標代碼中產生的冗余,一般主要產生于如下兩種情況。

(1)完成相同功能的不同語句在經過編譯器編譯后會產生不同的目標代碼,若選擇不合適的語句會引起冗余的目標代碼。這種冗余包括兩種,一種是目標代碼量的增加,指對于同樣功能的代碼由于實現方法不同導致編譯后代碼條數不同;另一種是代碼復雜度的增加,CISC指令集中的指令執行時間不同,因此同樣的指令條數下指令執行時間有所不同。

(2)為實現方法調用,操作系統需保存當前運行時信息,并在運行時棧開辟新的空間,在調用結束時恢復信息并釋放空間。這種調度工作會體現在目標代碼中,雖然一般認為正常的調用過程對程序性能影響不大,但當調用在程序中占比過大時會影響程序性能。調用引起代碼冗余的一種示例如下。

非調用版指令如下

s=a+b;

目標代碼為

1 mov %rdi,%rbp

2 mov %rax,rbp

3 add %rbx,%rax

調用版指令如下

Sum (inta,intb)

{returna+b;}

目標代碼為

1 push %rbp

2 push %rbx

3 sub %8,%rsp

4 mov %rdi,%rbp

5 mov %rax,rbp

6 add %rbx,%rax

7 add %8,%rsp

8 pop %rbx

9 pop %rbp

10 ret

調用版的第1、2行用來實現現場保存,第3行為新調用棧分配空間,然后執行加法功能,第7行釋放調用棧,第8、9行恢復現場,最后返回。可以看出,當實現同樣加法功能的s=a+b被替換為方法調用Sum(a,b)后,運行時棧調度的開銷超過了方法本身。雖然編譯器中有對高級語言的自動優化功能或調用產生冗余代碼的問題,如自動將短調用進行內聯處理,但是編譯器優化的前提是保證程序的安全,因此面對具有內存別名使用或變量類型未知等問題的代碼不能進行自動優化,需手動對高級語言進行優化。

分析本程序得出,導致性能下降的主要問題是:①部分庫函數調用時產生的隱藏代碼增加了代碼量,增加包括調用棧引起的代碼開銷和庫函數本身復雜性引起的代碼開銷;②代碼及代碼結構本身的計算復雜度較高,因為有較多耗時的除法計算需求。為解決上述問題,考慮使用簡單指令代替復雜指令,并結合訪存優化及流水優化考慮匯編指令在系統中的運行狀態進行優化。

下面討論指令優化的具體實現。針對庫函數調用的問題,結合程序分析發現,程序主要調用的庫函數為數學庫和容器類庫。數學庫的乘方、開方等操作一般不可替代,但是對于低階乘方,可以直接使用數乘代替。如pow(a,2)可以改寫為a*=a,這種替換可以減少調用開銷,同時簡化方法本身的復雜度。

對容器類庫函數的修改可考慮對C++語言自帶容器數組進行二次封裝代替庫容器,在本文的程序中為采用數組代替Vector容器。從可行性分析,Vector和數組有較為類似的存儲結構,因此在存取功能上具有可替代性。Vector相對數組的優勢在于靈活的操作方式和內存分配方式,程序中對容器的主要操作方式為隨機訪問,此外有少量的動態內存分配操作,因此在二次封裝數組時需要為數組增加隨機訪存函數,從而實現對Vector的替代,所以用數組代替Vector具有可行性。從必要性分析,Vector是STL封裝的容器庫,在使用時涉及庫調用會增加底層代碼量,容器的存儲結構上,Vector為實現動態增長,在分配空間時會分配額外的內存,在多維空間上,這種存儲空間的浪費變得更加明顯,在訪存時也會增加額外開銷,高速緩存的命中率也會下降。

針對代碼結構結合訪存特性考慮編譯器對語句的處理方式,由于對變量的訪問常被編譯為寄存器訪問操作,對引用的訪問被編譯為內存訪問操作,因此為優化訪存性能,需將循環中頻繁引用訪問修改為變量訪問。針對代碼本身,考慮科學計算中常出現的乘除法指令,在Intel提供的參考機中部分算數指令性能見表1,其中延遲表明實際運算周期數,發射表明兩次運算最小間隔周期數,容量表明該功能的運算單元數。

表1 Intel參考機算數指令性能

由表1可以看出,乘除法運算執行時間較長,除法尤為明顯,因此對乘除法的優化也應加以考慮。針對非浮點數的乘除法指令,采用移位運算代替乘除法的方法進行處理。該方法在乘除法數量較多且乘數或除數較為規律時可產生較好的優化效果。對于浮點數除法,可采用公共除法移出的方法將除法轉換為乘法,減少執行時間。

3 實驗結果

對第一節中的流體計算程序采用上述優化方案進行優化,并在商用服務器對各版本程序采用perf性能檢測工具對程序進行小規模性能分析。商用服務器采用Intel Xeon CPU E7-8850 @ 2.00 GHz(8×10核)CPU,內存125 GB,編譯器采用gcc 4.8.5。在三層網格上迭代輪次為5、5、10。程序優化前后的主要指標的分析結果和執行時間見表2。由于對容器以及對指令的優化,程序所使用的運行時空間的大小也有所下降,內存占用實驗結果見表3。

表2 程序主要指標優化前后的比較

表3 優化前后占用內存情況

由表2、3可以看出,程序性能的主要提升在于優化后導致指令數的減少,在系統級優化后緩存性能和流水性能也有所提升。在規模為5、5、10的程序中,經過優化的程序每周期指令數提升24.24%,指令暫停周期下降25.86%,流水中斷率下降6%,緩存命中率提升14%,執行時間下降70.56%,空間消耗下降49.97%。

將程序部署在TIANHE-1A超級計算機系統中測試性能,采用4節點運行36進程。采用的CPU為Intel Xeon X5670(2×6核),內存24 GB,編譯器為icc_xe_2013.0.07。三層網格迭代輪次采用50、50、10 000作為程序完整規模運行的實驗,測試結果見表4。

由表4結果可以看出,在TIANHE-1A系統上,三層網格迭代輪次分別為50、50、10 000計算規模時,程序耗時下降68.34%,空間消耗減少55.43%,與商用服務器上測試結果類似。根據測試結果可以看出,本文的優化方法能有效提升程序對CPU計算能力的發揮,大幅減少程序的底層代碼,在測試程序中取得較好的優化效果。

表4 程序完整規模運行的實驗測試結果

4 結 論

根據計算機系統單核級指令執行的加速特性和編譯特性,提出一種針對計算流體力學程序的單核優化方案,主要優化方面包括:①優化指令結構,減少指令流水的斷流,使程序充分發揮CPU的計算能力;②調整存儲結構和訪存結構,增強存儲和訪存的性能;③減少冗余指令,優化指令結構,從而減少程序執行所需的時間開銷和空間消耗。

將本文方法應用在壓氣機轉子模擬程序中,測試結果顯示優化方案能有效減少程序指令數、增加緩存命中率,減少指令流水斷流,在TIANHE-1A計算機系統中進行測試,結果顯示優化后時間性能提升約68.34%,空間性能提升約55.43%,加速比提升3倍以上,具有較好的優化效果。