基于流-固耦合的測(cè)量耙設(shè)計(jì)

賈文杰， 孫科， 宋江濤， 牛文敬

(中航工業(yè)中國(guó)飛行試驗(yàn)研究院， 西安 710089)

在航空發(fā)動(dòng)機(jī)飛行試驗(yàn)中，通常需要使用壓力和溫度測(cè)量耙測(cè)量發(fā)動(dòng)機(jī)流道內(nèi)的各種壓力和溫度信號(hào)，通過(guò)當(dāng)?shù)氐目倝?、靜壓、流動(dòng)方向和溫度等數(shù)據(jù)來(lái)確定發(fā)動(dòng)機(jī)的各個(gè)性能參數(shù)，因此測(cè)量耙的設(shè)計(jì)具有重要的意義。

航空發(fā)動(dòng)機(jī)流道參數(shù)測(cè)量主要通過(guò)測(cè)量耙來(lái)實(shí)現(xiàn)[1-2]。目前，針對(duì)測(cè)量耙的結(jié)構(gòu)強(qiáng)度設(shè)計(jì)與數(shù)值仿真方法已經(jīng)有了較多的研究。陶冶等[3-4]、田琳等[5]已經(jīng)針對(duì)測(cè)量耙的模態(tài)特性、動(dòng)態(tài)特性和應(yīng)力分布特性展開了相應(yīng)的研究；李秋鋒等[6]研究了不同結(jié)構(gòu)的測(cè)量耙對(duì)發(fā)動(dòng)機(jī)流場(chǎng)參數(shù)的影響。目前國(guó)內(nèi)外針對(duì)測(cè)量耙的研究已經(jīng)較為成熟，并制定了相應(yīng)的規(guī)范和標(biāo)準(zhǔn)[7-9]。國(guó)內(nèi)測(cè)量耙的設(shè)計(jì)流程為：首先設(shè)計(jì)出測(cè)量耙的基本結(jié)構(gòu)，然后通過(guò)有限元方法進(jìn)行靜強(qiáng)度計(jì)算、模態(tài)分析和根據(jù)《軍用裝備實(shí)驗(yàn)室環(huán)境試驗(yàn)方法 振動(dòng)試驗(yàn)》(GJB 150.16A)中渦輪發(fā)動(dòng)機(jī)振動(dòng)環(huán)境功率譜密度進(jìn)行隨機(jī)振動(dòng)分析，考核測(cè)量耙的強(qiáng)度能否滿足要求。但在某型號(hào)發(fā)動(dòng)機(jī)外涵出口測(cè)量耙設(shè)計(jì)過(guò)程中，發(fā)現(xiàn)通過(guò)有限元計(jì)算得到的靜應(yīng)力與隨機(jī)振動(dòng)應(yīng)力遠(yuǎn)遠(yuǎn)超過(guò)該溫度下的材料許用應(yīng)力，經(jīng)過(guò)理論分析，并結(jié)合以往類似工況下的測(cè)量耙設(shè)計(jì)方案，發(fā)現(xiàn)這是由于現(xiàn)有測(cè)量耙設(shè)計(jì)方法及計(jì)算方法存在的不足所導(dǎo)致的，經(jīng)過(guò)分析，現(xiàn)有的設(shè)計(jì)方法存在以下幾點(diǎn)問(wèn)題：一是靜強(qiáng)度計(jì)算中，使用測(cè)量耙安裝截面的動(dòng)態(tài)壓力作為均布面壓力載荷，模擬測(cè)量耙迎風(fēng)面受到的氣動(dòng)載荷，但實(shí)際上測(cè)量耙在流道的受力可以近似為圓柱繞流問(wèn)題[10-11]，測(cè)量耙無(wú)法將氣流完全滯止在耙體前方，現(xiàn)有的加載方式與測(cè)量耙實(shí)際受力條件相差較大，導(dǎo)致靜強(qiáng)度計(jì)算結(jié)果可靠性較低。二是模態(tài)分析時(shí)，未考慮測(cè)量耙所受氣動(dòng)力、溫度載荷等預(yù)應(yīng)力對(duì)測(cè)量耙模態(tài)特性的影響?，F(xiàn)有研究表明，結(jié)構(gòu)受到不可忽略的應(yīng)力載荷時(shí)，模態(tài)會(huì)發(fā)生變化，即“模態(tài)預(yù)應(yīng)力效應(yīng)”[12]。測(cè)量耙安裝在發(fā)動(dòng)機(jī)流道內(nèi)部，受到沿來(lái)流方向的氣動(dòng)力和流場(chǎng)中內(nèi)外涵溫度不均勻帶來(lái)的熱應(yīng)力，模態(tài)特性與靜止大氣環(huán)境中的不同。三是在實(shí)際的環(huán)境試驗(yàn)中，常規(guī)的掃頻試驗(yàn)無(wú)法準(zhǔn)確模擬測(cè)量耙在真實(shí)流道內(nèi)所受的氣動(dòng)力與熱載荷下的頻率，現(xiàn)有的試驗(yàn)方法也很難模擬真實(shí)裝機(jī)狀態(tài)下的測(cè)量耙模態(tài)特性和所受應(yīng)力載荷。

為解決現(xiàn)有設(shè)計(jì)方法中存在的不足，現(xiàn)綜合考慮測(cè)量耙所處的流場(chǎng)參數(shù)對(duì)測(cè)量耙強(qiáng)度和模態(tài)特性的影響，模擬測(cè)量耙在流場(chǎng)中所受到的真實(shí)載荷，選擇采用單向流-固耦合方法進(jìn)行測(cè)量耙強(qiáng)度計(jì)算，以提高現(xiàn)有測(cè)量耙設(shè)計(jì)水平，并應(yīng)用于后期測(cè)量耙設(shè)計(jì)中。

1 理論分析

流固耦合方法主要用于研究固體在流場(chǎng)作用下的受力以及固體位移對(duì)流場(chǎng)產(chǎn)生的影響[13-14]，主要分為單向流固耦合和雙向流固耦合，其中雙向流固耦合適用于流場(chǎng)產(chǎn)生的載荷使固體結(jié)構(gòu)產(chǎn)生較大的變形，從而對(duì)流場(chǎng)分布產(chǎn)生影響的算例中，雙向流固耦合計(jì)算精度較高，但占用的計(jì)算資源較大且收斂緩慢；單向流固耦合適用于固體域形變較小，且對(duì)流場(chǎng)帶來(lái)的影響較小的工況，所用計(jì)算資源小且容易收斂，在小變形條件下有著較好的計(jì)算精度，考慮到測(cè)量耙通常由高強(qiáng)度鋼或高溫合金制成，安裝在發(fā)動(dòng)機(jī)流道中變形較小，因此選用單向流固耦合方法進(jìn)行研究，以下流固耦合均指單向流固耦合。

在流固耦合計(jì)算過(guò)程中，流體計(jì)算域通過(guò)流固耦合交界面實(shí)現(xiàn)與固體計(jì)算域之間的數(shù)據(jù)傳遞，因此在流固耦合交界面上需要保持位移、應(yīng)力、熱流量以及溫度等變量之間的協(xié)調(diào)性，以保證數(shù)據(jù)傳遞的可靠性，需滿足以下4個(gè)守恒方程[15]，即

(1)

式(1)中：τ為應(yīng)力；n為法向量；d為位移；q為熱流量；T為溫度；下標(biāo)f為流體計(jì)算域，s為固體計(jì)算域。

考慮到熱應(yīng)力和氣動(dòng)力對(duì)測(cè)量耙模態(tài)特性的影響，在結(jié)構(gòu)面上附加一個(gè)預(yù)應(yīng)力剛度矩陣，結(jié)合幾何非線性問(wèn)題中的有限變形理論，采用拉格朗日坐標(biāo)系中的增量位能原理，建立有限元方程，將這些影響因素轉(zhuǎn)換為測(cè)量耙結(jié)構(gòu)的預(yù)應(yīng)力剛度矩陣，并將其疊加到結(jié)構(gòu)的剛度矩陣上，最后生成考慮預(yù)應(yīng)力和流固耦合的方程[16]，具體為

(2)

式(2)中：M為結(jié)構(gòu)質(zhì)量矩陣；Mf為流體等效質(zhì)量矩陣；ρf為流體密度；R為流固耦合矩陣；C為阻尼；K為結(jié)構(gòu)剛度矩陣；Kr和Kf分別為預(yù)應(yīng)力剛度矩陣及流體等效剛度矩陣；F為結(jié)構(gòu)外載荷；p為流體壓力；δ為結(jié)構(gòu)位移，上標(biāo)一個(gè)點(diǎn)和兩個(gè)點(diǎn)分別表示該參數(shù)的一階導(dǎo)和二階導(dǎo)。

求解上述問(wèn)題時(shí)，首先進(jìn)行靜態(tài)條件下的線性應(yīng)力分析，將外部載荷轉(zhuǎn)化為結(jié)構(gòu)的額外剛度。靜態(tài)分析時(shí)不考慮阻尼，且不考慮自由流場(chǎng)中流體的位移帶來(lái)的壓力變化，即

Kσθ=F

(3)

然后，基于靜態(tài)分析的應(yīng)力位移σθ獲得預(yù)應(yīng)力等效矩陣S，結(jié)合預(yù)應(yīng)力等效矩陣并忽略阻尼，得到測(cè)量耙的模態(tài)方程為

(4)

2 仿真計(jì)算

2.1 模型介紹

圖1(a)所示為某型發(fā)動(dòng)機(jī)外涵出口總溫總壓測(cè)量耙在發(fā)動(dòng)機(jī)流道中的裝配圖，測(cè)量耙使用外裝方式固定在發(fā)動(dòng)機(jī)混合器機(jī)匣壁面上，在內(nèi)涵機(jī)匣波瓣混合器上留有跑道形的安裝孔，測(cè)量耙通過(guò)該跑道形的小孔伸入內(nèi)涵流道中，波瓣混合器幾何結(jié)構(gòu)如圖1(b)所示。

對(duì)上述模型進(jìn)行單向流固耦合計(jì)算，考慮到計(jì)算資源的合理分配，對(duì)計(jì)算域進(jìn)行如下簡(jiǎn)化。

(1)測(cè)量耙在加工過(guò)程中內(nèi)部鋪設(shè)有引氣管路及相關(guān)附件，在測(cè)量耙模型中對(duì)其管路結(jié)構(gòu)進(jìn)行了相應(yīng)的簡(jiǎn)化。

(2)發(fā)動(dòng)機(jī)機(jī)匣和波瓣混合器均為圓柱形中心對(duì)稱結(jié)構(gòu)，為節(jié)約計(jì)算資源，取發(fā)動(dòng)機(jī)流道1/8扇段進(jìn)行研究。采用的計(jì)算流程如圖2所示。

2.2 流場(chǎng)分析

首先對(duì)測(cè)量耙進(jìn)行流場(chǎng)分析，使用的計(jì)算流體動(dòng)力學(xué)(computational fluid dynamics，CFD)流場(chǎng)計(jì)算軟件為Fluent，使用壓力求解器，湍流模型為Standardk-ε，近壁面處使用標(biāo)準(zhǔn)壁面函數(shù)，求解采用SIMPLEC方法，使用壓力進(jìn)口和壓力出口邊界條件，進(jìn)出口條件參數(shù)如表1所示。選用四面體網(wǎng)格進(jìn)行網(wǎng)格劃分，經(jīng)過(guò)網(wǎng)格無(wú)關(guān)性驗(yàn)證，最終選擇的流場(chǎng)計(jì)算域網(wǎng)格數(shù)量為320萬(wàn)。

圖1 計(jì)算模型介紹Fig.1 Introduction to calculation model

如圖3(a)所示為流場(chǎng)計(jì)算得到的測(cè)量耙表面壓力分布云圖，圖3(b)為采用傳統(tǒng)方法進(jìn)行靜強(qiáng)度計(jì)算時(shí)對(duì)測(cè)量耙迎風(fēng)面施加的均布面壓力載荷?？?/p>

表1 進(jìn)出口參數(shù)

圖2 流固耦合計(jì)算流程示意圖Fig.2 Schematic diagram of coupling calculation process

圖3 測(cè)量耙表面壓力載荷云圖Fig.3 Pressure contour of the rake surface

以看出，這與流場(chǎng)計(jì)算得到的測(cè)量耙表面真實(shí)受力有很大不同，在實(shí)際的發(fā)動(dòng)機(jī)流道中，測(cè)量耙迎風(fēng)面的壓力載荷并非均布載荷，靠近耙體中心位置可以完全滯止氣流，因此壓力最大；氣流從耙體兩側(cè)流過(guò)，因此迎風(fēng)面兩側(cè)壓力較低，這與現(xiàn)有圓柱繞流中的研究結(jié)果相符。

2.3 靜強(qiáng)度分析

測(cè)量耙耙體暴露在外涵和內(nèi)涵流道中，內(nèi)外涵溫度和壓力差別較大，因此受到熱載荷與氣動(dòng)壓力載荷，將流場(chǎng)計(jì)算得到的壓力載荷與溫度載荷作為輸入?yún)?shù)，進(jìn)行氣動(dòng)載荷分析。該測(cè)量耙由GH4169材料制成，材料密度為8 240 kg/m3，熔化溫度范圍為1 260～1 320 ℃，各溫度下彈性模量E及泊松比μ如表2所示，雙真空熔煉的圓餅鍛件持久強(qiáng)度如表3所示。

在ANSYS中進(jìn)行分析，采用以四面體為主的網(wǎng)格劃分方法，共劃分402 446個(gè)網(wǎng)格節(jié)點(diǎn)，確保有限元網(wǎng)格尺寸與流場(chǎng)計(jì)算中流固耦合面的網(wǎng)格尺寸完全相同。圖4為耙體溫度分布，可以看出，測(cè)量耙內(nèi)涵部分和外涵部分溫度相差較大，考慮到發(fā)動(dòng)機(jī)內(nèi)涵熱變形，波瓣混合器上預(yù)留的跑道型孔尺寸大于測(cè)量耙內(nèi)涵部分直徑，導(dǎo)致有一部分內(nèi)涵高溫氣體流入外涵，因此測(cè)量耙過(guò)渡段背風(fēng)面存在局部高溫區(qū)。結(jié)合表3中不同溫度下材料的持久強(qiáng)度變化可以看出，溫度對(duì)測(cè)量耙的影響無(wú)法忽略，原有靜強(qiáng)度計(jì)算方法存在顯著不足，流固耦合方法更能準(zhǔn)確模擬測(cè)量耙真實(shí)的受力情況。

表2 GH4169材料不同溫度下的彈性模量和泊松比

表3 GH4169材料持久強(qiáng)度

從表4中可以看出，兩種方法靜強(qiáng)度計(jì)算結(jié)果相差較大，靜應(yīng)力由248.69 MPa減小為51.49 MPa，變形由0.81 mm減小為0.16 mm，分析原因如下。

(1)原有計(jì)算方法中假設(shè)氣體被測(cè)量耙完全滯止，因此在靜強(qiáng)度計(jì)算時(shí)施加的壓力載荷過(guò)大。

(2)原始方法不考慮熱應(yīng)力的影響，考慮到測(cè)量耙的安全裕度，計(jì)算時(shí)一般通過(guò)插值得到測(cè)量耙在流場(chǎng)最高溫度下的彈性模型、泊松比等材料屬性，過(guò)于保守導(dǎo)致材料本身強(qiáng)度較低。

經(jīng)過(guò)對(duì)比，可以看出流固耦合方法在測(cè)量耙靜強(qiáng)度計(jì)算方面效果更優(yōu)。

表4 靜強(qiáng)度及變形分布

2.4 模態(tài)分析

測(cè)量耙安裝于發(fā)動(dòng)機(jī)流道內(nèi)，受到發(fā)動(dòng)機(jī)轉(zhuǎn)子所產(chǎn)生的振動(dòng)載荷，當(dāng)測(cè)量耙的固有頻率和氣流的激振頻率或轉(zhuǎn)子的轉(zhuǎn)動(dòng)頻率重合時(shí)，就會(huì)產(chǎn)生諧共振，造成耙體的損壞，因此需要進(jìn)行模態(tài)分析，確保測(cè)量耙頻率可以避開發(fā)動(dòng)機(jī)轉(zhuǎn)子固有頻率。

將前文計(jì)算得到的靜應(yīng)力作為預(yù)應(yīng)力，開展預(yù)應(yīng)力載荷下的測(cè)量耙模態(tài)分析，表5為采用原有設(shè)計(jì)方法和流固耦合方法計(jì)算得到的模態(tài)對(duì)比，可以看出，采用流固耦合方法計(jì)算得到的測(cè)量耙前兩階固有頻率分別提高了13.96%和11.9%，分析原因可能是由于原有方法未考慮測(cè)量耙在實(shí)際流場(chǎng)中受到的氣動(dòng)力對(duì)模態(tài)的影響，以及兩種計(jì)算方法中材料參數(shù)設(shè)置的不同導(dǎo)致。采用流固耦合的方法可以更加真實(shí)的模擬測(cè)量耙在發(fā)動(dòng)機(jī)流道內(nèi)所受的熱載荷和氣動(dòng)載荷。

圖4 耙體溫度分布云圖Fig.4 Temperature contour of the rake

表5 兩種方法計(jì)算模態(tài)結(jié)果對(duì)比

3 結(jié)論

針對(duì)現(xiàn)有測(cè)量耙設(shè)計(jì)方法中存在的不足，提出采用流固耦合的方法開展測(cè)量耙在流場(chǎng)中的氣動(dòng)載荷分析，然后以流場(chǎng)分析的結(jié)果作為輸入?yún)?shù)進(jìn)行測(cè)量耙的模態(tài)特性分析，結(jié)論如下。

(1)原有方法進(jìn)行靜強(qiáng)度計(jì)算時(shí)，將測(cè)量耙安裝截面動(dòng)壓作為均布載荷施加在耙體迎風(fēng)面上，這與實(shí)際裝機(jī)狀態(tài)相差較大，原有模態(tài)計(jì)算方法未考慮氣動(dòng)力等預(yù)應(yīng)力對(duì)測(cè)量耙模態(tài)特性的影響，導(dǎo)致計(jì)算結(jié)果誤差較大。

(2)采用流固耦合方法進(jìn)行計(jì)算，將CFD流場(chǎng)計(jì)算的溫度和壓力載荷作為靜強(qiáng)度輸入?yún)?shù)，與傳統(tǒng)方法對(duì)比發(fā)現(xiàn)，靜應(yīng)力由248.69 MPa減小為51.49 MPa，變形由0.81 mm減小為0.16 mm。

(3)將計(jì)算得到的靜應(yīng)力作為預(yù)應(yīng)力進(jìn)行模態(tài)分析，與原有方法相比，前六階模態(tài)均增長(zhǎng)超過(guò)7%，分析是由于傳統(tǒng)測(cè)量耙設(shè)計(jì)方法不能準(zhǔn)確模擬測(cè)量耙真實(shí)的裝機(jī)狀態(tài)，且設(shè)計(jì)方法過(guò)于保守，導(dǎo)致計(jì)算得到的靜應(yīng)力偏大，模態(tài)值偏低。

(4)考慮到流固耦合方法在渦輪葉片設(shè)計(jì)等領(lǐng)域已經(jīng)有了越來(lái)越廣泛的應(yīng)用，其可靠性已經(jīng)得到充分驗(yàn)證，因此流固耦合設(shè)計(jì)方法可以為提高現(xiàn)有的測(cè)量耙設(shè)計(jì)能力提供參考，同時(shí)應(yīng)用于后續(xù)測(cè)量耙設(shè)計(jì)及校核計(jì)算中。