技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析

動圖

01 引言:CFD大規模仿真效率瓶頸

1.1 背景與目標

  隨著工業數值模擬在航空航天、能源、船舶及裝備制造等領域的廣泛應用,對計算流體力學(Computational Fluid Dynamics, CFD)軟件的性能要求日益提升。某基于有限體積方法的通用 CFD 仿真軟件,能夠模擬穩態與瞬態下的復雜流動,并支持旋轉機械流動分析與先進的共軛熱傳遞(CHT)模型。結合特定的前處理器,該軟件可以快速生成計算網格,從而為用戶提供從前處理到數值求解的完整解決方案。

  然而,隨著工程問題規模的不斷擴大,傳統代數求解器在處理大規模稀疏線性方程組時暴露出性能瓶頸。具體而言,現有求解器在并行計算中的通信開銷較大,尤其是代數多重網格(AMG)方法,難以充分發揮高性能計算平臺的并行潛力,從而限制了整體計算效率和擴展性。

  此次高性能改造的目標是對該CFD商軟中的代數求解器進行替換與優化,開發出性能更高、并行效率更好的數值求解模塊,以解決大規模并行通信開銷過大的問題。通過此次改造,期望顯著提升CFD軟件在大規模計算場景下的計算速度與可擴展性,為復雜工程問題的高效仿真提供支持。

  結合軟件的主要使用場景,本次改造確定弱擴展性并行效率測試的環境為: 400萬以上網格算例,以單核為基準計算10倍擴展時的并行效率。

1.2 高性能改造的必要性

  在 CFD 數值求解過程中,代數求解器往往占據總計算時間的很大一部分,是影響整體性能的關鍵環節。當前求解器在單機計算中能夠滿足一定規模的計算需求,但在大規模算例的并行環境下,通信與同步開銷逐漸成為性能瓶頸,尤其是在不可壓流動壓力方程的求解中,AMG作為經典的求解方法,存在著收斂性與并行效率之間的矛盾。這種性能限制不僅導致計算效率下降,還影響到軟件在高性能平臺上的可擴展性。

  因此,有必要對代數求解器進行高性能改造,采用先進的并行求解方法與優化策略,以實現以下目標:

  • 提高求解速度:通過優化算法和數據結構,減少單次迭代的計算與通信開銷;
  • 提升并行效率:增強在并行計算環境下的擴展性,充分利用多核和多節點計算資源;
  • 增強工程適用性:支持更大規模、更復雜的工業仿真計算需求。

02 解析:AMGCL效率衰減實證

2.1 原始算法架構

  最初階段神工坊?技術團隊采用了 AMGCL 作為代數求解庫。AMGCL 是一個基于 C++ 的開源代數多重網格(Algebraic MultiGrid, AMG)求解器框架,具有模塊化、可擴展和高性能的特點,廣泛應用于大規模稀疏線性方程組的求解。

  其基本架構可以分為以下幾個層次:

  • 前端接口層 

  提供用戶友好的 API,用于構建稀疏矩陣、設置預條件子與迭代方法,并驅動求解流程。

  • 預條件子模塊

  AMGCL 的核心在于代數多重網格預條件子構建。它通過稀疏矩陣自動生成層次化的粗網格與限制/延拓算子,從而加速 Krylov 子空間迭代方法的收斂。

  • 迭代求解器模塊

  支持多種 Krylov 子空間方法,例如 CG、BiCGStab 等,通常與多重網格預條件子結合使用,以實現高效收斂。

  • 后端執行層

  提供多種計算后端以適配不同硬件環境,包括:

  • CPU 并行后端(OpenMP、Threading Building Blocks 等)
  • GPU 加速后端(CUDA、OpenCL 等)
  • 分布式后端(MPI)
  • 并行化與擴展性

  AMGCL 通過模板化設計實現后端解耦,能夠在不同硬件架構和規模下快速切換。并行化策略主要體現在粗網格構建與迭代求解階段。

2.2 主要瓶頸分析

  整體上,AMGCL 提供了一個相對成熟且易于使用的求解框架,為后續的高性能改造打下了基礎。但在面對超大規模計算和高并行環境時,其在內存訪問效率、通信開銷等方面仍存在優化空間。

1. 計算性能瓶頸

  盡管AMGCL中使用了許多OpenMP指令進行線程級別優化,但受稀疏矩陣迭代法核心計算特點限制,反而可能產生一些同步/規約開銷,尤其是AMG方法中,圖遍歷 + 不規則訪問 + 分支判斷,OpenMP 細粒度任務化開銷大。使得通過OpenMP使用多個CPU核心進行計算的收益沒有那么明顯。

2. 內存帶寬瓶頸

  稀疏線性方程組的求解并不是算法強度(浮點運算數/內存訪問字節數)很高的運算,使得其主要瓶頸受限于訪存效率而不是CPU的計算能力。以迭代法的核心操作:稀疏矩陣-向量乘(SpMV)為例,假設用CSR存儲, 每個非零元需要:

  • 讀數值(8字節,double)
  • 讀列索引(4或8字節)
  • 讀對應的向量元素(8字節)
  • 計算只涉及一個乘法 + 一個加法(2 FLOPs)
  • 結果:每20-24 字節數據,只有2FLOPs → 算法強度 ≈ 0.1 FLOPs/Byte

  根據Roofline Model,該強度較低,位于帶寬瓶頸區域 Memory-Bound,內存帶寬(斜率)決定了FLOPS上限。

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖2

 此外,并行效率很容易受內存帶寬瓶頸影響,在某個并行規模開始產生效率急劇下降的現象。下面通過測試進行驗證:

 測試環境:12核CPU Intel? Xeon? Processor E5-2650 v4其  Roofline Model  圖如下所示,可見其峰值性能約206.3GFLOPS,內存帶寬約為29.4GB/s

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖3

  測試了400w網格壓力方程使用AMGCL的AMG預條件+BiCGStab并行求解,發現10進程的并行效率僅35.83%

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖4

  可以發現4進程開始并行效率下降十分迅速,分析其計算性能參數結果如下:

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖5

  分析數據可以得出以下結論:

  1. 單核性能受限于內存帶寬,單進程時算法強度僅 0.097 FLOP/Byte,遠小于現代 CPU 的計算能力所需的平衡點,和上文對SpMV的算法強度分析結果一致;
  2. 隨進程數增加,單核計算能力在下降:
  • 從 1 到 2 核,MFLOP/s 從 467 提升到 919,幾乎翻倍
  • 到 4 核時 1635 MFLOP/s,約為單核的 3.5 倍
  • 到 8 核時 2263 MFLOP/s,僅約為單核的 4.8 倍
  1. 內存帶寬隨并行數擴展而增長,但4進程開始出現倍數降低,說明共享內存總線成為瓶頸,出現多個進程同時爭奪內存帶寬。

3. 通信瓶頸

  為了測試AMGCL的AMG方法通信瓶頸,需要了解其真實的并行效率。上節中使用Release版本的AMGCL測試發現并行效率受內存帶寬的影響很大,為了避免該因素造成的干擾,我們使用Debug版本進行測試,其計算性能參數結果如下:

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖6

  可見Debug版本計算能力隨進程數增加接近線性提升,帶寬利用率遠低于Release版本,并行效率更多受MPI通信影響。

  此時的400w網格壓力方程使用AMGCL的AMG預條件+BiCGStab并行求解測試結果如下:

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖7

  結果顯示AMGCL以單核為基準計算10倍擴展時的并行效率為50%左右,通過Score-P等性能分析工具分析發現,其在AMG方法中的通信量極大,這和并行稀疏矩陣-矩陣乘(SpGEMM)的實現有關:

以**行式乘法(row-wise-product)數據流**為例,下式描述了通過行式乘法數據流計算結果矩陣一行的方法。行式乘法數據流首先讀取 $A$ 矩陣的一行,然后根據這一行中非零元素的列索引,來讀取 $B$ 矩陣對應的多行。與此同時,從 $B$ 矩陣中讀出的多行需要乘上 $A$ 矩陣中對應的非零元素,然后按照矩陣乘法的要求把得到的中間列表聚合成結果矩陣中的一行。

$$
C_{i*} = \sum_k A_{ik} \cdot B_{k*}
$$

 

  在并行環境中,通常需要通過MPI通信得到矩陣第行一整行的數據,且每個迭代步AMG方法都要重構一次,從而產生了大量的通信。

03 破局:UNAP三大優化策略

 改造使用非結構代數求解庫UNAP(UNstructured Algebra Package),它是一款最初面向異構平臺的、采用 C++ 語言的輕量級代數求解庫,目前支持申威眾核異構平臺和普通 x86 平臺,支持Windows和Linux操作系統。

3.1 稀疏矩陣存儲格式優化

  UNAP的迭代求解部分采用LDU格式。LDU格式矩陣分為三個部分,上三角塊,對角塊和下三角塊(如下圖中所示),由于上三角塊和下三角塊的結構是對稱的,因此可以共用同一套索引,如果以上三角塊為基準,則在上三角塊中的行標、列標剛好變成下三角塊中的列標、行標,這樣就只需要存儲一半的拓撲數據。如果矩陣結構也是對稱的,那么上、下三角的數據就可以只存一塊。

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖8


  LDU矩陣中幾個主要數組的說明如下:

  •  diag  、 lower 、 upper 存具體數值,說明如上圖
  •  upperAddr : 存上三角列索引,相當于圖中的 col 
  •  lowerAddr : 存下三角列索引,相當于圖中的 row 
  • 為了能讓上下三角對應起來,unap的LDU格式要求數據按照一定的順序排列,即上三角按行優先的順序存儲,下三角按列優先的順序存儲 ,上下三角(對角塊)只儲存列標(局部編號)。
  • 進程邊界部分(非對角塊),按照對應進程號的不同分成了若干個部分,實際上可以說是COO格式,儲存行標(局部編號)和列標(全局編號)。

3.2 迭代算法優化

基于LDU格式可以進行多種預條件子/光滑器的優化,以不完全LU分解預條件為例(Incomplete LU Factorization with zero fill-in,ILU(0)),它的核心思想是:
在進行 LU 分解時,只保留原稀疏矩陣中已有的非零位置(pattern),不引入新的填充元素(fill-in)。  
給定稀疏矩陣 $A \in \mathbb{R}^{n \times n}$,目標是找到下三角矩陣 $L$ 和上三角矩陣 $U$,使得:

$$
A≈LU
$$

其中 $L$ 的對角元為 1。

基于CSR的ILU預條件通常需要構建并儲存$L$和$U$矩陣,并通過相應的上下迭代過程完成預條件。而LDU矩陣格式由于天然有上下三角的劃分,只需要通過簡單地循環便可完成與條件過程,減少了構建過程的開銷,其訪存效率也比逐行遍歷時要高。  

3.3 代數多重網格方法優化

AMGCL中使用一種進程邊界優先的PMIS聚集方法,同時使用smooth aggregation以達到更快的收斂效率,但代價是運算量和通信量更大

$$
P = \left( I - \omega D^{-1} A^F \right) \tilde P.
$$

 其中矩陣$A^F = (a_{ij}^F)$ 元素計算如下: 

$$
\begin{aligned}
a_{ij}^F &= \begin{cases}
a_{ij} \quad \text{if} \; j \in N_i\\
0 \quad \text{otherwise}
\end{cases}, \quad \text{if}\; i \neq j, \\
\quad a_{ii}^F &= a_{ii} - \sum\limits_{j=1,j\neq i}^n
\left(a_{ij} - a_{ij}^F \right),
\end{aligned}
$$

式中$N_i$ 是$j$的強連接集合,$D$ 代表$A^F$的對角部分。

    代數多重網格方法優化路線從以下三方面開展

  1. 使用混合的聚集策略以在代價增長一定的前提下提升收斂的可控性與穩健性;

  在細層基于強連接進行圖聚集,得到聚集簇并構建延拓矩陣$P_t$,優化聚集過程以減少粗層矩陣通信界面大??;在粗層適當使用光滑聚集以提高收斂速度。  

   2. SpGEMM過程優化;

  粗層矩陣通過$A_c=RAP$三重乘積構建,需要兩次SpGEMM,其算法過程分為容量估計和數值寫入兩個階段,同時需要整理對應行數據進行通信,可通過一定手段使用非阻塞通信將計算過程與之覆蓋,減少實際耗時開銷。  

   3. 重復使用結構信息的AMG;

  稀疏矩陣的系數值是方程物理特性和幾何信息的結合體,在系統矩陣$A$的稀疏模式不變的場景下,可以通過一定的聚集方法來構造可以重復使用的插值/限制算子($P_l$、$R_l$),從而避免每次都需要從零開始進行AMG的SetUp過程。也可以通過迭代數上限等閾值觸發策略來控制重構插值/限制算子的時機。  

04 驗證:多平臺實測性能對比

  性能評估主要從兩個方面進行:

  • 特定CFD軟件之外、矩陣單獨求解效率對比;
  • 特定CFD軟件之中,整體求解情況對比;

4.1 測試環境與配置

  • x86超算平臺:CentOS7.9系統,64核CPU AMD EPYC 7H12 64-Core Processor 
  • Windows單機: Windows10系統,12核CPU(8個Performance-core) Intel? Core? i7-13700F 

4.2 性能測試方法

  • 主要測試案例:CHT算例矩陣階數為4254912、非零元個數27580140; AHM算例矩陣階數為4639609、非零元個數31434615。
  • 矩陣單獨求解效率對比
  • AMGCL: 使用gcc7.3.1+OpenMPI-4.1.0+cmake編譯
  • UNAP: 使用gcc7.3.1+OpenMPI-4.1.0+cmake編譯
  • 主要對比速度方程/壓力方程的單次求解效率
  • 整體求解情況對比
  • AMGCL: 使用Windows單機默認編譯環境編譯
  • UNAP:使用MSMPI+Visual Studio2019+cmake編譯
  • 對比所有方程的多次求解效率,迭代步數設為100
  • 主要性能指標計算公式
  
    $$\text{加速比} = \frac{T_{\text{串行}}}{T_{\text{并行}}}$$
  
    $T_{\text{串行}}$:程序在單個處理器上運行的執行時間。  
    $T_{\text{并行}}$:程序在 $P$ 個處理器上運行的執行時間。
  
    $$ \text{并行效率} = \frac{\text{加速比}}{P/1} = \frac{T_{\text{串行}} \times 1}{T_{\text{并行}} \times P} $$
  
    $P$:使用的處理器數量。  
    $1$:表示以串行為基準。  
    $T_{\text{串行}}$:程序在單個處理器上運行的執行時間。  
    $T_{\text{并行}}$:程序在 $P$ 個處理器上運行的執行時間。
  • 擴展性:強擴展性(Strong Scaling)和弱擴展性(Weak Scaling)是評估并行程序在不同并行資源(如處理器或節點)配置下性能表現的兩種方法。
  •   強擴展性是指在固定總工作量的情況下,通過增加處理器數量來減少任務完成時間的能力。具體來說,就是在處理器數量增加時,觀察程序執行時間的變化情況。
  •   弱擴展性是指在增加處理器數量的同時,也按比例增加問題規模,以考察程序在不同規模問題上的性能表現。通常,弱擴展性測試的目標是保持每個處理器的計算負載不變,觀察程序的效率。
  •   對于內存受限的程序,強擴展性測試需要特別注意內存帶寬和容量的限制,常見的測試方式有:a. 在內存帶寬不受限的平臺上進行測試; b. 選擇足夠大的問題和并行規模,使得主要矛盾從內存帶寬轉移到通信、同步和負載均衡。

4.3 測試結果與分析

  • 速度方程求解都采用: ILU(0)預條件+ BiCGStab
  • 壓力方程求解都采用: AMG預條件+ BiCGStab
  • $Ax=b$求解結果均通過解向量的輸出對比以保證正確性

1. x86超算平臺

  超算平臺上測試不受內存帶寬影響,并行效率僅受程序的通信、同步和負載均衡影響,其上限由代碼的串行部分占比決定,能夠真實的反映出程序的強擴展性。

  (1)CHT算例

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖9

  從 CHT 速度方程的求解結果來看,UNAP 在并行加速性能上優于 AMGCL。在進程數增加時,UNAP 不僅保持了較低的求解時間,而且在 2-8 進程范圍內出現了超線性加速現象,并行效率一度超過 100%,說明其在緩存利用(cache命中率)和內存帶寬利用(單位訪存帶寬能完成更多的浮點運算)方面具有優勢。

  

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖10

  在壓力方程的求解中,UNAP 與 AMGCL 在單進程下性能相近,隨進程數增加兩者用時均明顯減少,UNAP曲線更偏下一些。在并行擴展性上,UNAP 的效率始終保持在較高水平,而 AMGCL 的效率則隨規模增大快速下降,在大規模并行時僅為前者的約七成,10倍擴展并行效率為 56% ,不及UNAP的 73% ,顯示出AMGCL在通信與負載均衡方面存在瓶頸更明顯。

  (2)AHM算例

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖11

  在 AHM 速度方程算例中,UNAP 在所有進程數下的求解用時均明顯低于 AMGCL,體現出更優的單進程性能與整體效率;而在并行擴展性方面,兩者表現接近,AMGCL 的并行效率略高于 UNAP,但并未以抵消其速度上的性能差距。

   

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖12

  在AHM壓力方程的求解中,UNAP 與 AMGCL 求解用時相近,串行時 AMGCL 略快,并行后 UNAP耗時更短; UNAP的并行效率始終在 AMGCL 的上方,10倍擴展并行效率為 63% 多,而 AMGCL 為 43% 。

2. Windows單機

  (1)正確性驗證

  • AHM算例:通過對比100步迭代過程力系數的變化驗證結果正確性

  結果對比如下,可見當50步之后力系數逐漸趨于穩定后、不同進程數  AMGCL/UNAP并行計算得到的力系數值沒有差別。  

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖13

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖14

  • CHT算例:通過對比800步迭代過程出口熱流的變化驗證結果正確性

  結果對比如下圖,可見AMGCL/UNAP并行計算得到的熱流量變化曲線幾乎沒有差別。

  該算例中線性方程組求解耗時統計如下:

  AMGCL:

Max iluBiCGStab Time: 978.763s, count: 4800,  average: 0.203909s
Max amgCG Time: 976.65s, count: 813,  average: 1.20129s

 UNAP:

Max iluBiCGStab Time: 685.757s, count: 4800,  average: 0.142866s
Max amgCG Time: 616.886s, count: 813,  average: 0.758777s

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖15

  (2)CHT算例

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖16

  CHT算例測試結果如上圖所示,可以發現ILU-BiCGStab方法UNAP始終在AMGCL曲線下方,串行時耗時減少35.3%,10進程并行時耗時減少27.8% 。同樣AMGCG方法曲線UNAP始終在AMGCL下方,串行時耗時減少31.6%、10進程并行時耗時減少17.4% 。

  (3)AHM算例

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖17

  AHM算例測試結果如上圖所示,可以發現ILU-BiCGStab方法UNAP始終在AMGCL曲線下方,串行時耗時減少36.2%,10進程并行時耗時減少23.2% 。同樣AMGCG方法曲線UNAP始終在AMGCL下方, 串行時耗時減少29.3%、10進程并行時耗時減少24.3% 。

05 結語

5.1 成果總結:36%性能突破

  神工坊?技術團隊針對CFD商業軟件在并行場景下代數求解器存在的性能瓶頸,開展了評估與高性能改造工作。通過替換 AMGCL 為 UNAP,并在稀疏矩陣存儲格式、迭代算法及代數多重網格方法等方面開展優化,取得了以下主要成果:

1. 整體性能提升

  在典型算例中,通過混合聚集策略、SpGEMM優化與重用AMG結構信息等方法,UNAP 相較原方法在單次迭代求解時間上實現了加速,100步迭代ILU-BiCGStab耗時減少23.2%-36.2%、AMG-CG減少17.4%-31.6%,提升了矩陣求解效率。

2. 并行效率改善

  在 400 萬以上網格的強擴展性測試中,UNAP 的AMG并行效率相比 AMGCL 明顯提高,10 倍擴展規模下并行效率保持在60-70% 左右,驗證了該方法在并行場景下的有效性。

3. 矩陣數據結構替換

  LDU矩陣格式在減少存儲與訪存冗余的同時,提高了部分預條件子的應用效率;同時LDU格式與CFD軟件內部的網格相關數據結構更加契合,為矩陣構建、組裝等其他模塊的進一步優化奠定了基礎。

5.2 未來展望:邁向AI求解新階段

  本次高性能改造在提升CFD軟件代數求解器的性能與并行效率方面取得了可觀成效,未來神工坊?技術團隊的研究與開發工作將主要在AI 驅動的代數求解加速方向攻堅克難。神工坊?將引入人工智能方法對代數求解過程進行智能化優化,包括:自動選擇最優求解器與預條件子、智能調節參數、以及基于深度學習的迭代初值預測,以進一步提升收斂速度和魯棒性。

技術分享︱突破大規模CFD仿真瓶頸:UNAP代數求解庫性能實測與優化解析的圖18

登錄后免費查看全文
立即登錄
App下載
技術鄰APP
工程師必備
  • 項目客服
  • 培訓客服
  • 平臺客服

TOP