壓氣機(jī)數(shù)字孿生模型中的仿真技術(shù)探索

2023-07-07 09:46:46高麗敏俞一波李瑞宇王可鑫

航空發(fā)動(dòng)機(jī) 2023年3期

趙磊，高麗敏，俞一波，李瑞宇，王可鑫

（西北工業(yè)大學(xué)動(dòng)力與能源學(xué)院，西安 710021）

0 引言

數(shù)字孿生技術(shù)是未來(lái)工業(yè)極具潛力和挑戰(zhàn)的發(fā)展方向，依托該技術(shù)所構(gòu)建形成的航空發(fā)動(dòng)機(jī)數(shù)字孿生體[1]，對(duì)于解決航空發(fā)動(dòng)機(jī)研制中面臨的多系統(tǒng)、多部件、多物理場(chǎng)耦合等挑戰(zhàn)性問(wèn)題提供了新的思路。壓氣機(jī)作為航空發(fā)動(dòng)機(jī)的3 大核心部件，其數(shù)字孿生體的研制對(duì)于形成航發(fā)整機(jī)的孿生體具有重要意義。

以壓氣機(jī)為物理實(shí)體構(gòu)建孿生模型，首先要保證壓氣機(jī)的數(shù)值模型能夠與現(xiàn)實(shí)環(huán)境中的物理實(shí)體一一對(duì)應(yīng)，即壓氣機(jī)數(shù)值模型具備完備性要求。“孿生”技術(shù)最早由美國(guó)國(guó)家航空航天局（NASA）在阿波羅項(xiàng)目[2]中提出；Grieves 等[3]指出孿生體是物理實(shí)體在虛擬空間的“數(shù)字復(fù)制品”；NASA 后期[4]發(fā)布了數(shù)字孿生的路線(xiàn)圖，其中數(shù)值仿真模擬為該技術(shù)的主要部分；Tuegel 等[5]、Stargel 等[6]、Cerrone 等[7]也指出數(shù)字孿生體需要對(duì)復(fù)雜系統(tǒng)中不同空間和時(shí)間的事件進(jìn)行準(zhǔn)確數(shù)值仿真預(yù)測(cè)；孫明霞等[8]、張志博等[9]、劉永泉等[10]以及王樂(lè)等[11]通過(guò)分析現(xiàn)代數(shù)字孿生技術(shù)，指出航空發(fā)動(dòng)機(jī)的數(shù)字孿生技術(shù)中最核心的部分為孿生模型相關(guān)的數(shù)值模型的構(gòu)建。可見(jiàn)，數(shù)字孿生技術(shù)以數(shù)值仿真作為其核心技術(shù)手段。然而，在現(xiàn)有航空發(fā)動(dòng)機(jī)研制過(guò)程中，壓氣機(jī)的數(shù)值仿真多以定常流動(dòng)、均勻進(jìn)氣、單通道周期性假設(shè)等建立數(shù)值模型，所建立的數(shù)字孿生體不能完全反映壓氣機(jī)在真實(shí)環(huán)境下的工作性能等參數(shù)變化。因此，從計(jì)算流體力學(xué)角度出發(fā)構(gòu)建壓氣機(jī)數(shù)字孿生體，其孿生模型需要具備開(kāi)展真實(shí)邊界條件下多級(jí)全環(huán)的非定常數(shù)值仿真能力。然而，壓氣機(jī)數(shù)字孿生模型中的全環(huán)非定常的計(jì)算需求巨大。Garbor 等[12]提出數(shù)字孿生的計(jì)算成本在于仿真所需的大量計(jì)算資源；Yamada 等[13]對(duì)某燃?xì)廨啓C(jī)7 級(jí)壓氣機(jī)開(kāi)展全環(huán)計(jì)算，使用了近7300 多CPU核心的計(jì)算資源；Gourdain 等[14]指出未來(lái)采用大渦模擬開(kāi)展多級(jí)壓氣機(jī)全環(huán)計(jì)算需要具有E級(jí)（運(yùn)算峰值性能達(dá)到百億億次，即1018flop/s）計(jì)算能力的高性能計(jì)算機(jī)；Slotnick 等[15]指出航空發(fā)動(dòng)機(jī)整機(jī)在非設(shè)計(jì)工況流動(dòng)的數(shù)值模擬（其中包含壓氣機(jī)多級(jí)全環(huán)非定常數(shù)值模擬）將是未來(lái)高性能并行計(jì)算領(lǐng)域的重大挑戰(zhàn)課題；Belzner等[16]曾提出數(shù)字孿生體中事件的模擬速度要盡可能比現(xiàn)實(shí)中的事件發(fā)生速度更快。可見(jiàn)，針對(duì)未來(lái)壓氣機(jī)數(shù)字孿生技術(shù)的發(fā)展必然要借助高性能并行計(jì)算技術(shù)，以提升數(shù)值仿真手段的計(jì)算能力和速度。

本文結(jié)合數(shù)值仿真方法給出了航發(fā)數(shù)字孿生技術(shù)中數(shù)值模型完備性的定義；針對(duì)某S 彎進(jìn)氣道中涵道風(fēng)扇的數(shù)字孿生體，對(duì)比分析了不同數(shù)值模型對(duì)孿生結(jié)果的影響規(guī)律；針對(duì)現(xiàn)代高性能計(jì)算平臺(tái)，給出了壓氣機(jī)數(shù)字孿生模型中的2 層次并行計(jì)算模型，并開(kāi)展了多級(jí)壓氣機(jī)全環(huán)非定常計(jì)算驗(yàn)證。

1 航發(fā)孿生技術(shù)中數(shù)值模型完備性定義

航空發(fā)動(dòng)機(jī)的數(shù)字孿生技術(shù)是集成多學(xué)科、多物理場(chǎng)、多尺度、概率性的仿真系統(tǒng)。航發(fā)孿生體是航空發(fā)動(dòng)機(jī)的物理過(guò)程在數(shù)字虛擬世界中的完美映射，該孿生體由許多相互關(guān)聯(lián)的數(shù)值模型組成，各模型綜合起來(lái)用以描述發(fā)動(dòng)機(jī)運(yùn)行過(guò)程中的氣動(dòng)、傳熱、燃燒、強(qiáng)度等物理過(guò)程及其性能。要保證數(shù)字孿生體的正確性，首先需要滿(mǎn)足組成孿生體相關(guān)不同數(shù)值模型的完備性要求。由此，定義針對(duì)任意給定誤差δ≥0 和可描述航空發(fā)動(dòng)機(jī)實(shí)時(shí)物理過(guò)程的控制方程Φ(x?)及其定解條件，若存在一個(gè)數(shù)值可解的計(jì)算模型φ(x?)使得

在給定解條件下的所有合法輸入x?均成立，那么認(rèn)為該數(shù)值可解的計(jì)算模型及其定解條件對(duì)于航發(fā)數(shù)字孿生體是計(jì)算完備的，即數(shù)值模型具備完備性要求。工程實(shí)踐中多采用計(jì)算模型的結(jié)果與真實(shí)環(huán)境下試驗(yàn)結(jié)果之間的誤差來(lái)定義δ。

壓氣機(jī)中氣動(dòng)物理過(guò)程可用N-S方程及定解條件對(duì)其進(jìn)行描述，相應(yīng)的數(shù)值模型可結(jié)合現(xiàn)代計(jì)算流體力學(xué)方法（Computational Fluid Dynamics，CFD）給出?，F(xiàn)有的CFD 中以有限體積方法居多，該方法將連續(xù)N-S方程及其定解條件離散化為物理空間中多個(gè)網(wǎng)格計(jì)算單元上的代數(shù)方程組，通過(guò)數(shù)值求解該方程組可近似獲得壓氣機(jī)中的氣動(dòng)物理過(guò)程。因此，從計(jì)算流體力學(xué)角度出發(fā)，滿(mǎn)足壓氣機(jī)數(shù)值模型完備性的要求不僅要保證離散化連續(xù)方程過(guò)程中引入的離散誤差、舍入誤差等盡可能小，同時(shí)描述氣體流動(dòng)湍流狀態(tài)等物理模型的誤差也要盡可能小。此外，幾何模型、邊界條件、初始條件等也要滿(mǎn)足一定的誤差要求。綜上，在滿(mǎn)足給定誤差δ的前提下，采用CFD 方法的壓氣機(jī)數(shù)值模型具備對(duì)壓氣機(jī)氣動(dòng)物理過(guò)程的計(jì)算完備性要求。

2 某涵道風(fēng)扇/壓氣機(jī)數(shù)字孿生中數(shù)值模型完備性分析

選取某帶S 彎進(jìn)氣道的涵道風(fēng)扇/壓氣機(jī)（如圖1所示）作為研究對(duì)象，氣流過(guò)S 彎進(jìn)氣道后，在進(jìn)氣道與涵道風(fēng)扇/壓氣機(jī)的氣動(dòng)交界面（Aerodynamic Inter?face Plane，AIP）上形成總壓畸變圖譜，如圖2 所示。對(duì)該涵道風(fēng)扇/壓氣機(jī)建立數(shù)字孿生模型可以采用3種不同數(shù)值模型方案（暫不考慮進(jìn)氣道和壓氣機(jī)的耦合效應(yīng)）：（1）將AIP 界面總壓進(jìn)行平均后取平均總壓作為進(jìn)氣邊界，即可采用現(xiàn)有的均勻進(jìn)氣周期性假設(shè)的單通道計(jì)算模型，如圖3 所示；（2）由于AIP 界面主要以總壓畸變?yōu)橹?，且以周向總壓畸變?yōu)橹鳎虼丝蓮脑糀IP 界面總壓畸變圖譜上將其簡(jiǎn)化為周向總壓畸變，如圖4 所示。其中，從AIP 界面測(cè)得周向畸變角度為130.9°，總壓畸變強(qiáng)度為8%。選取簡(jiǎn)化后的周向總壓畸變圖譜作為壓氣機(jī)進(jìn)氣邊界條件，并采用全通道幾何模型，如圖5 所示；（3）直接將AIP 界面的總壓畸變圖譜（圖2）給定為壓氣機(jī)進(jìn)口邊界條件，同時(shí)采用全通道計(jì)算模型（圖5）。以上3 種數(shù)值模型的邊界條件和幾何模型見(jiàn)表1。

表1 涵道風(fēng)扇不同數(shù)字孿生模型對(duì)比

圖1 帶S彎進(jìn)氣道涵道風(fēng)扇/壓氣機(jī)

圖2 AIP界面總壓畸變圖譜

圖3 均勻進(jìn)氣與周期性假設(shè)的單通道數(shù)值模型

圖4 通過(guò)AIP界面獲得的簡(jiǎn)化周向總壓畸變圖譜

圖5 涵道風(fēng)扇/壓氣機(jī)全通道計(jì)算模型

3 種孿生模型獲得的壓氣機(jī)特性如圖6 所示，其中數(shù)值方法以及網(wǎng)格無(wú)關(guān)性驗(yàn)證見(jiàn)文獻(xiàn)[17]。通過(guò)對(duì)比3 種不同數(shù)值模型可見(jiàn)，在相同轉(zhuǎn)速工況下，采用均勻進(jìn)氣和單通道假設(shè)的數(shù)值模型1 的壓比和效率要高于模型2、3 的，其中模型1 最高效率比模型2的高3.8%，比模型3 的高5.5%，同時(shí)相應(yīng)壓比也有不同程度的差異。此外從所獲得的壓氣機(jī)穩(wěn)定裕度可以看出，采用模型1 所得到的該涵道風(fēng)扇的穩(wěn)定裕度要遠(yuǎn)大于采用模型2、3 所獲得的穩(wěn)定裕度，其中采用模型2 所獲得的穩(wěn)定裕度僅為采用模型1 所獲得的50.5%，裕度范圍幾乎縮小一半。同樣的，對(duì)比模型2、3 可見(jiàn)，數(shù)值模型都采用全通道計(jì)算模型，僅是進(jìn)氣邊界的不同，采用模型3 所獲得的穩(wěn)定裕度范圍已經(jīng)縮小到采用模型2所獲得的29.4%。如果假設(shè)模型3（AIP 界面畸變圖譜為進(jìn)氣邊界條件，采用全通道計(jì)算模型）是該涵道風(fēng)扇的完備數(shù)值模型，通過(guò)對(duì)比可見(jiàn)，采用傳統(tǒng)單通道周期性和均勻進(jìn)氣假設(shè)所獲得的壓氣機(jī)性能整體參數(shù)將完全不適合作為該壓氣機(jī)的數(shù)字孿生模型。

圖6 3種孿生模型獲得的壓氣機(jī)特性

全通道模型不同進(jìn)氣邊界條件下涵道風(fēng)扇出口總壓分布如圖7 所示。二者均采用全通道計(jì)算模型，僅是進(jìn)口邊界條件不同。對(duì)比出口總壓分布云圖可見(jiàn)，不同進(jìn)口總壓畸變分布（兩者總壓畸變強(qiáng)度均相同），所得到的涵道風(fēng)扇出口的總壓分布存在較大差異，其中采用模型3 出口的高壓區(qū)的分布不均勻性明顯高于采用模型2 的，其壓力分布的不均勻更加明顯。由此可見(jiàn)，在數(shù)字孿生模型中，即使幾何模型一致，進(jìn)、出口邊界條件的差異（保持總壓畸變強(qiáng)度一定）也會(huì)使得進(jìn)氣畸變經(jīng)過(guò)風(fēng)扇轉(zhuǎn)子和靜子傳遞后，內(nèi)部流場(chǎng)特征存在較大差異。

圖7 全通道模型不同進(jìn)氣邊界條件下涵道風(fēng)扇出口總壓分布

綜合以上分析可知，涵道風(fēng)扇孿生模型的不同對(duì)“現(xiàn)實(shí)環(huán)境”中的風(fēng)扇壓氣機(jī)總體性能參數(shù)影響非常大，孿生模型即壓氣機(jī)數(shù)值計(jì)算模型的不同會(huì)對(duì)壓氣機(jī)本身數(shù)字孿生體的內(nèi)部流場(chǎng)以及宏觀性能參數(shù)會(huì)產(chǎn)生巨大影響。針對(duì)進(jìn)氣畸變的孿生場(chǎng)景研究中發(fā)現(xiàn)，幾何模型（單通道周期性假設(shè)）的不同對(duì)于所獲得的性能總體參數(shù)的影響最大。

3 數(shù)字孿生模型中的高效計(jì)算技術(shù)

航發(fā)孿生技術(shù)對(duì)高性能計(jì)算技術(shù)有著較高要求，如何充分挖掘現(xiàn)有高性能計(jì)算潛能以推動(dòng)數(shù)值孿生技術(shù)的發(fā)展，是未來(lái)航發(fā)數(shù)字孿生需要突破的瓶頸之一。一方面，計(jì)算機(jī)技術(shù)的發(fā)展正在促使高性能的浮點(diǎn)計(jì)算性能朝著E 級(jí)（百億億次浮點(diǎn)計(jì)算能力）計(jì)算發(fā)展，甚或是未來(lái)量子計(jì)算成為可能，都有望促使基于3 維的物理過(guò)程仿真與現(xiàn)實(shí)世界同步；另一方面，對(duì)諸如壓氣機(jī)等航發(fā)部件實(shí)現(xiàn)部件級(jí)物理過(guò)程仿真（例如實(shí)現(xiàn)壓氣機(jī)整機(jī)全環(huán)非定常數(shù)值模擬），對(duì)于航發(fā)孿生體不同部件在不同工作環(huán)境下數(shù)字孿生體的準(zhǔn)確構(gòu)建與驗(yàn)證具有重要意義。因此，針對(duì)壓氣機(jī)，未來(lái)壓氣機(jī)數(shù)字孿生模型不僅需要具備多級(jí)全環(huán)非定常的并行計(jì)算能力，同時(shí)相應(yīng)仿真計(jì)算能力的構(gòu)建需要能夠充分挖掘高性能計(jì)算的浮點(diǎn)計(jì)算潛力。

現(xiàn)有高性能計(jì)算平臺(tái)均采用分布式計(jì)算系統(tǒng)，大都具有機(jī)柜、節(jié)點(diǎn)、多CPU、多計(jì)算核心等不同層次的硬件架構(gòu)特征。本文結(jié)合該計(jì)算架構(gòu)特征，發(fā)展了針對(duì)壓氣機(jī)多級(jí)計(jì)算的兩層次并行計(jì)算模型，2 層次并行計(jì)算模型如圖8 所示。結(jié)合計(jì)算系統(tǒng)中計(jì)算進(jìn)程、計(jì)算線(xiàn)程在不同層次架構(gòu)中的計(jì)算特征，可將上述硬件結(jié)構(gòu)特征抽象為2 層計(jì)算模型，上層主要以計(jì)算進(jìn)程為主，下層主要以計(jì)算線(xiàn)程為主。針對(duì)該2 層計(jì)算特征，可分別采用不同的并行編程模型。例如針對(duì)計(jì)算進(jìn)程特征層，一般采用MPI 編程模型；針對(duì)計(jì)算線(xiàn)程特征層，可結(jié)合具體的硬件架構(gòu)特征，采用OPENMP、CUDA 等編程模型。計(jì)算進(jìn)程可采用區(qū)域分解算法對(duì)計(jì)算對(duì)象進(jìn)行任務(wù)級(jí)分解，可充分發(fā)揮不同節(jié)點(diǎn)、不同CPU 上的計(jì)算潛力。計(jì)算線(xiàn)程層次主要針對(duì)現(xiàn)有的異構(gòu)加速硬件，針對(duì)任務(wù)級(jí)中的計(jì)算密集型任務(wù)實(shí)現(xiàn)高效加速。

圖8 2層次并行計(jì)算模型

結(jié)合上述2 層次高性能計(jì)算模型，本文自主發(fā)展了相應(yīng)的高性能計(jì)算軟件，對(duì)某15 級(jí)壓氣機(jī)開(kāi)展了全環(huán)非定常計(jì)算驗(yàn)證。15 級(jí)壓氣機(jī)全通道計(jì)算網(wǎng)格如圖9所示。針對(duì)該15級(jí)壓氣機(jī)的全環(huán)模型，轉(zhuǎn)靜交界面的處理采用滑移網(wǎng)格法，計(jì)算時(shí)網(wǎng)格量約為1.5×109，每排葉排網(wǎng)格數(shù)約為4.6×107，單個(gè)葉片通道的網(wǎng)格量大于106，y+為1。共使用256 個(gè)計(jì)算節(jié)點(diǎn)進(jìn)行并行計(jì)算，每個(gè)計(jì)算節(jié)點(diǎn)為2 路CPU，每顆CPU 均具有48個(gè)計(jì)算核心。該多級(jí)壓氣機(jī)全環(huán)非定常計(jì)算采用3072 個(gè)計(jì)算進(jìn)程，每個(gè)計(jì)算進(jìn)程采用4 個(gè)計(jì)算線(xiàn)程實(shí)現(xiàn)加速計(jì)算。

圖9 15級(jí)壓氣機(jī)全通道計(jì)算網(wǎng)格

15 級(jí)壓氣機(jī)軸向速度瞬時(shí)流場(chǎng)分布如圖10 所示。為了深入對(duì)比不同時(shí)刻多級(jí)壓氣機(jī)內(nèi)部流場(chǎng)的非定常變化，前4級(jí)壓氣機(jī)50%葉不同時(shí)刻的壓力分布如圖11 所示。從圖中可見(jiàn)，受到上下游葉片排的相互影響，在葉片排軸向間隙中的流場(chǎng)變化隨著動(dòng)葉的旋轉(zhuǎn)都在變化。第3級(jí)動(dòng)葉排出口靜壓分布隨時(shí)間的變化如圖12 所示。從圖中第3 級(jí)動(dòng)葉排出口的壓力分布可見(jiàn)，不同時(shí)刻第3 級(jí)動(dòng)葉排出口的壓力分布的周向不均勻性變化明顯，且在不同時(shí)刻的壓力分布的不均勻性均存在。第10~12 級(jí)壓氣機(jī)50%葉高在不同時(shí)刻的壓力分布和第11 級(jí)動(dòng)葉排出口靜壓分布隨時(shí)間的變化如圖13、14 所示。對(duì)比前面級(jí)的仿真結(jié)果可見(jiàn)，多級(jí)計(jì)算中不同葉排間、單葉排不同葉道在不同時(shí)刻的流場(chǎng)結(jié)果均呈現(xiàn)出差異，若僅是采用單通道定常假設(shè)的數(shù)值模型是無(wú)法獲得該流場(chǎng)的細(xì)節(jié)差異。

圖10 15級(jí)壓氣機(jī)軸向速度瞬時(shí)流場(chǎng)分布

圖11 前4級(jí)壓氣機(jī)50%葉不同時(shí)刻的壓力分布

圖12 第3級(jí)動(dòng)葉排出口靜壓分布隨時(shí)間的變化

圖13 第10~12級(jí)壓氣機(jī)50%葉高在不同時(shí)刻的壓力分布

圖14 第11級(jí)動(dòng)葉排出口靜壓分布隨時(shí)間的變化

綜合上述分析可見(jiàn)，壓氣機(jī)數(shù)字孿生模型要具備多級(jí)全環(huán)非定常計(jì)算能力離不開(kāi)高性能計(jì)算技術(shù)支撐。同時(shí)結(jié)合現(xiàn)有高性能計(jì)算平臺(tái)發(fā)展適合于壓氣機(jī)整機(jī)非定常計(jì)算的并行計(jì)算模型，對(duì)于發(fā)掘壓氣機(jī)內(nèi)部流動(dòng)的非定常流動(dòng)特性極具必要性。

4 總結(jié)

（1）針對(duì)S 彎進(jìn)氣道中涵道風(fēng)扇建立了數(shù)字孿生的不同數(shù)值模型，主要從進(jìn)氣邊界模型和幾何模型2方面進(jìn)行了分析。通過(guò)研究發(fā)現(xiàn)，采用均勻進(jìn)氣邊界和單通道周期性假設(shè)計(jì)算模型與采用AIP 界面原始畸變圖譜、全通道計(jì)算模型所獲得的壓氣機(jī)性能參數(shù)差別巨大，尤其是壓氣機(jī)穩(wěn)定裕度范圍完全不同?？梢?jiàn)，從數(shù)值模型完備性角度出發(fā)，壓氣機(jī)數(shù)字孿生模型需要同時(shí)在均勻進(jìn)氣和畸變進(jìn)氣環(huán)境下獲得準(zhǔn)確結(jié)果時(shí)，采用全環(huán)計(jì)算模型構(gòu)建壓氣機(jī)數(shù)字孿生體較為合適。