渦輪葉片累積損傷指數(shù)模型及服役可靠性評(píng)估

雷世英，孫見忠，劉赫

南京航空航天大學(xué) 民航學(xué)院，南京 211106

民航發(fā)動(dòng)機(jī)熱端部件的故障和失效是造成發(fā)動(dòng)機(jī)性能退化及非計(jì)劃換發(fā)的一個(gè)重要因素，其很大程度上制約了發(fā)動(dòng)機(jī)整機(jī)的使用壽命。按照失效造成的影響嚴(yán)重性，可將熱端部件分為A、B兩大類。A類件又稱為限壽件(Life Limited Part, LLP)，包括渦輪軸、渦輪盤等主要旋轉(zhuǎn)件，其失效后會(huì)造成機(jī)毀人亡的嚴(yán)重后果，工程上一般通過安全性分析給定安全壽命，當(dāng)達(dá)到規(guī)定的使用循環(huán)后進(jìn)行強(qiáng)制報(bào)廢，這種方法由于沒有考慮個(gè)體發(fā)動(dòng)機(jī)實(shí)際的載荷差異，壽命估計(jì)一般較為保守；B類件又稱為關(guān)鍵件，主要是對(duì)發(fā)動(dòng)機(jī)性能影響較大的結(jié)構(gòu)件，如壓氣機(jī)/渦輪葉片、軸承等，其失效的后果沒有A類件嚴(yán)重，但維修成本較高，發(fā)生故障時(shí)還可能造成較大的二次破壞，工程上一般采取視情維修的策略，通過評(píng)估剩余壽命來制訂相應(yīng)的維修計(jì)劃，因此其壽命預(yù)測(cè)的準(zhǔn)確與否對(duì)發(fā)動(dòng)機(jī)的安全性和經(jīng)濟(jì)性有著重要影響，也是發(fā)動(dòng)機(jī)壽命預(yù)測(cè)領(lǐng)域廣泛關(guān)注的問題。如今發(fā)動(dòng)機(jī)制造商和運(yùn)營(yíng)商對(duì)安全性和經(jīng)濟(jì)性要求越來越高，越來越趨向于在確保安全的前提下，盡可能地提高發(fā)動(dòng)機(jī)的可用性，最大限度降低壽命周期運(yùn)行成本，實(shí)現(xiàn)安全性和經(jīng)濟(jì)性的最大化。然而實(shí)際中要想實(shí)現(xiàn)壽命的準(zhǔn)確預(yù)測(cè)是相當(dāng)困難的，熱端部件的使用壽命不僅受設(shè)計(jì)制造、材料工藝水平的影響，還與發(fā)動(dòng)機(jī)的實(shí)際運(yùn)行環(huán)境、使用方式以及維修模式等有著密切關(guān)聯(lián)，這些多方面的因素使得預(yù)測(cè)的不確定性大大增加，預(yù)測(cè)精度很難滿足實(shí)際需要。

目前對(duì)于關(guān)鍵件服役壽命的預(yù)測(cè)，總體方法上可以分為統(tǒng)計(jì)模型法、物理模型法和多模型、多數(shù)據(jù)融合的方法。在不具備零部件設(shè)計(jì)信息的情況下，一般使用可靠性分析的方法評(píng)估葉片的使用可靠性及剩余壽命，如An等采用貝葉斯方法融合外場(chǎng)可靠性數(shù)據(jù)，得到更新后的壽命分布來確定葉片壽命；Zaretsky等根據(jù)發(fā)動(dòng)機(jī)外場(chǎng)使用數(shù)據(jù)分析了葉片的使用壽命，最終推導(dǎo)出一個(gè)估算葉片剩余壽命的計(jì)算公式；Li等采用模糊理論來表示渦輪葉片壽命預(yù)測(cè)中涉及的不確定性，利用線性融合算法對(duì)渦輪葉片的可靠性壽命進(jìn)行區(qū)間評(píng)估。上述統(tǒng)計(jì)模型法需要提供大量的輸入數(shù)據(jù)，但搜集數(shù)據(jù)的時(shí)間成本較高，具體的實(shí)驗(yàn)數(shù)據(jù)和現(xiàn)場(chǎng)數(shù)據(jù)也較為缺乏；除此之外，壽命可靠性預(yù)測(cè)涉及小樣本問題，得到的結(jié)果反映的是相似使用條件下零部件使用可靠性的平均屬性，難以體現(xiàn)單機(jī)狀態(tài)下發(fā)動(dòng)機(jī)自身及運(yùn)行環(huán)境、運(yùn)行載荷等方面的差異性。物理模型法則是建立在掌握了詳細(xì)的材料參數(shù)與部件設(shè)計(jì)參數(shù)的基礎(chǔ)之上，通過微觀失效機(jī)理研究、計(jì)算流體力學(xué)與結(jié)構(gòu)有限元法等分析手段，得到關(guān)鍵部位的應(yīng)力、溫度載荷，進(jìn)而憑借失效機(jī)理對(duì)應(yīng)的壽命計(jì)算模型，對(duì)壽命損耗進(jìn)行評(píng)估，如Staroselsky等建立了基于晶體形態(tài)的彈粘塑性模型，考慮了材料滑移層的本構(gòu)行為，并將其耦合到有限元模型中，使用半經(jīng)驗(yàn)壽命預(yù)測(cè)模型估算葉片的剩余壽命；荷蘭NLR中心和德國(guó)DLR中心等分別針對(duì)發(fā)動(dòng)機(jī)渦輪葉片，通過實(shí)驗(yàn)設(shè)計(jì)(Design of Experience, DOE)方法構(gòu)建了葉片溫度和應(yīng)力計(jì)算的代理模型，并考慮了蠕變和熱機(jī)疲勞兩種損傷模式來預(yù)測(cè)葉片剩余壽命；Zhu等通過引入種載荷參數(shù)來考慮渦輪葉片高低周疲勞損傷的耦合效應(yīng)，使得預(yù)測(cè)誤差的均值和標(biāo)準(zhǔn)差都達(dá)到了較高的精度；Brand?o等等通過逆向工程得到了渦輪葉片的材料成分，借助有限元仿真比較了等溫模型和非等溫模型的蠕變壽命預(yù)測(cè)結(jié)果。上述物理模型的預(yù)測(cè)方法依托于材料力學(xué)、流體力學(xué)、熱力學(xué)、有限單元法等基本理論，利用發(fā)動(dòng)機(jī)的實(shí)際運(yùn)行參數(shù)，獲得葉片關(guān)鍵部位的應(yīng)力和溫度載荷譜，進(jìn)而借助壽命損耗模型來評(píng)估渦輪葉片的剩余壽命。此類方法需要掌握發(fā)動(dòng)機(jī)詳細(xì)的設(shè)計(jì)數(shù)據(jù)，葉片的應(yīng)力與溫度載荷分析需要進(jìn)行復(fù)雜的計(jì)算流體力學(xué)與結(jié)構(gòu)有限元分析，但在模型計(jì)算時(shí)做了大量假設(shè)，會(huì)造成較大的預(yù)測(cè)不確定性。而多模型、多數(shù)據(jù)融合的方法則可以充分利用已有的多源異構(gòu)數(shù)據(jù)，在葉片失效機(jī)理的基礎(chǔ)上綜合考量使用條件、運(yùn)行環(huán)境等差異，相較于單模型方法能夠提供更高的預(yù)測(cè)精度。Pillai等提出了基于有限元代理模型的多源信息融合方法，利用機(jī)器學(xué)習(xí)方法融合外場(chǎng)失效等數(shù)據(jù)，進(jìn)一步降低預(yù)測(cè)結(jié)果的不確定性；Giesecke等結(jié)合發(fā)動(dòng)機(jī)使用參數(shù)及運(yùn)行區(qū)域大氣環(huán)境條件建立貝葉斯置信網(wǎng)絡(luò)用來預(yù)測(cè)葉片壽命損耗，以輔助決策葉片修理級(jí)別；Fu等考慮了微觀結(jié)構(gòu)準(zhǔn)則和蠕變應(yīng)變準(zhǔn)則，并結(jié)合人工神經(jīng)網(wǎng)絡(luò)和改進(jìn)的投影模型來評(píng)估使用中的渦輪葉片的剩余蠕變壽命，開發(fā)了一種基于材料基因工程的渦輪葉片在役蠕變剩余壽命的集成計(jì)算方法。

實(shí)際上，反映發(fā)動(dòng)機(jī)關(guān)鍵件服役可靠性降級(jí)等行為特性的信息存在于多模態(tài)數(shù)據(jù)中，包括產(chǎn)品設(shè)計(jì)數(shù)據(jù)與使用維護(hù)數(shù)據(jù)，如產(chǎn)品幾何參數(shù)、材料數(shù)據(jù)、狀態(tài)監(jiān)測(cè)數(shù)據(jù)、離線檢測(cè)/檢查報(bào)告、故障記錄和大修報(bào)告等，貫穿產(chǎn)品壽命周期的不同階段，如何充分利用這些多模態(tài)數(shù)據(jù)，實(shí)現(xiàn)多領(lǐng)域物理模型與數(shù)據(jù)融合，是提高服役發(fā)動(dòng)機(jī)關(guān)鍵件壽命預(yù)測(cè)分析準(zhǔn)確性的重要途徑。本文以民航發(fā)動(dòng)機(jī)高壓渦輪(HPT)葉片蠕變失效為例，在蠕變失效機(jī)理分析基礎(chǔ)上，同時(shí)考慮服役條件下的使用信息、狀態(tài)參數(shù)及截尾失效時(shí)間等數(shù)據(jù)，提出的蠕變累積損傷指數(shù)模型實(shí)現(xiàn)了多模態(tài)信息的融合，對(duì)渦輪葉片服役可靠性進(jìn)行評(píng)估，并預(yù)測(cè)其剩余壽命?，F(xiàn)階段外場(chǎng)數(shù)據(jù)獲取有限，因此論文通過仿真數(shù)據(jù)驗(yàn)證提出方法的可行性和有效性，結(jié)果表明，本文方法比傳統(tǒng)可靠性分析方法預(yù)測(cè)精度有顯著提高。

文章的主要內(nèi)容安排如下：第1部分主要介紹渦輪葉片的蠕變失效機(jī)理，闡述了蠕變的發(fā)展階段、微觀失效原理以及影響蠕變損傷的常見因素，并引入蠕變損傷評(píng)估模型；第2部分主要介紹葉片蠕變壽命計(jì)算的整個(gè)流程，包括發(fā)動(dòng)機(jī)性能仿真、葉片應(yīng)力計(jì)算、葉片溫度計(jì)算及蠕變損傷評(píng)估等過程；第3部分介紹葉片蠕變累積損傷指數(shù)模型，給出了損傷指數(shù)模型的定義方法以及使用該模型進(jìn)行剩余壽命預(yù)測(cè)的流程；第4部分為模型的仿真驗(yàn)證部分，使用實(shí)際發(fā)動(dòng)機(jī)的狀態(tài)監(jiān)測(cè)數(shù)據(jù)進(jìn)行了全壽命仿真，基于仿真數(shù)據(jù)驗(yàn)證蠕變累積損傷指數(shù)模型，預(yù)測(cè)HPT葉片的可靠性壽命，并與基準(zhǔn)模型進(jìn)行了結(jié)果對(duì)比分析；第5部分是對(duì)全文的總結(jié)以及后續(xù)工作的展望。

1 渦輪葉片蠕變失效機(jī)理

渦輪葉片作為典型的熱端部件，長(zhǎng)期工作在高溫、高壓、高轉(zhuǎn)速的極限環(huán)境下，承受著離心載荷、氣動(dòng)載荷、溫度載荷、振動(dòng)載荷及高溫氧化腐蝕的綜合作用，其主要的失效模式包括：熱機(jī)疲勞在內(nèi)的低循環(huán)疲勞、振動(dòng)引起的高循環(huán)疲勞、高溫長(zhǎng)時(shí)間載荷作用下的蠕變變形和蠕變應(yīng)力斷裂；高溫燃?xì)鉀_刷腐蝕和氧化以及外物損傷等。對(duì)于現(xiàn)代高涵道比民航發(fā)動(dòng)機(jī)而言，其飛行過程平穩(wěn)、飛行任務(wù)較為單一、飛行時(shí)間較長(zhǎng)，這種工況下HPT葉片的蠕變損傷是決定其使用壽命的主要因素之一。同時(shí)，對(duì)于蠕變壽命的數(shù)值方法研究已經(jīng)有了較為成熟的模型和方法，在工程實(shí)踐中得到了廣泛的應(yīng)用。因此以葉片的典型失效模式—蠕變失效為例，來研究累積損傷指數(shù)建模以及外場(chǎng)數(shù)據(jù)融合的方法。

1.1 蠕變失效機(jī)理

蠕變是指工作在高溫、高應(yīng)力水平條件下的部件，隨著加載時(shí)間的延長(zhǎng)緩慢地發(fā)生塑性變形的現(xiàn)象。一般情況下，當(dāng)材料的工作溫度大于或等于其0.4～0.6倍熔點(diǎn)溫度，且工作應(yīng)力遠(yuǎn)低于材料的屈服強(qiáng)度時(shí)就會(huì)產(chǎn)生較為顯著的蠕變變形。如圖1所示，以蠕變時(shí)間為橫軸，蠕變導(dǎo)致的材料變形量為縱軸，可以發(fā)現(xiàn)蠕變的發(fā)展一般分為3個(gè)階段：① 初始蠕變階段(→)，材料已發(fā)生一定形變，初始應(yīng)變?yōu)?，蠕變速率隨時(shí)間逐漸降低；② 二次蠕變階段(→)，蠕變速率基本保持恒定；③ 三次蠕變階段(→)，蠕變速率隨時(shí)間逐漸增加，直到蠕變時(shí)間發(fā)生蠕變斷裂,此時(shí)變形量為。

圖1 蠕變發(fā)展階段Fig.1 Development stages of creep

渦輪葉片的蠕變壽命對(duì)溫度的升高十分敏感，研究表明，葉片的溫度每升高10～15 ℃，其蠕變使用壽命將減少一半。同時(shí)，蠕變是一個(gè)依賴時(shí)間的熱變形過程，葉片暴露在高溫下的時(shí)間越長(zhǎng)，其累積變形量越大。對(duì)于發(fā)動(dòng)機(jī)渦輪葉片而言，工程上一般將蠕變損傷作為葉片應(yīng)力、葉片溫度和持續(xù)時(shí)間的函數(shù)。

1.2 蠕變損傷評(píng)估模型

葉片設(shè)計(jì)時(shí)所要求的持久強(qiáng)度一般在數(shù)萬小時(shí)，通過長(zhǎng)時(shí)間的試驗(yàn)獲得相關(guān)的蠕變數(shù)據(jù)是不現(xiàn)實(shí)的，而想要精確地建立起預(yù)測(cè)蠕變損傷和剩余壽命的模型就越發(fā)困難。工程上一般使用由實(shí)驗(yàn)驗(yàn)證的蠕變持久方程來預(yù)測(cè)葉片的蠕變持久強(qiáng)度。目前，國(guó)內(nèi)外較為成熟的葉片設(shè)計(jì)規(guī)范和手冊(cè)普遍使用拉森-米勒(L-M)參數(shù)方程、葛-唐吾(G-D)方程、曼森-索柯普(M-S)方程、曼森-哈弗特(M-H)方程等方法進(jìn)行蠕變壽命的計(jì)算和實(shí)驗(yàn)驗(yàn)證。選用工程上常用的L-M參數(shù)方程預(yù)測(cè)渦輪葉片的蠕變持久強(qiáng)度，該方法基于材料的蠕變?cè)囼?yàn)數(shù)據(jù)，使用拉森米勒參數(shù)來表征材料在不同溫度和應(yīng)力組合下的蠕變斷裂行為。針對(duì)葉片特定的某點(diǎn)，拉森米勒參數(shù)滿足：

(1)

式中：為該點(diǎn)的溫度載荷，K；為該點(diǎn)的應(yīng)力水平，MPa；為蠕變斷裂時(shí)間，h；為介于17～23之間的常數(shù)，一般情況取=20。同時(shí)，應(yīng)力水平與拉森米勒參數(shù)之間的關(guān)系還可以用一個(gè)玻爾茲曼函數(shù)來表示：

(2)

式中：、、、是玻爾茲曼函數(shù)的參數(shù)。式(2) 可以看作是拉森米勒方程在全應(yīng)力和溫度范圍內(nèi)的擴(kuò)展，模型的參數(shù)可以由材料在不同溫度下的蠕變斷裂實(shí)驗(yàn)數(shù)據(jù)和抗拉強(qiáng)度數(shù)據(jù)通過插值得到。由式(1)可以推導(dǎo)出在一定的應(yīng)力、溫度載荷作用下，持續(xù)時(shí)間內(nèi)的蠕變損傷為

(3)

2 渦輪葉片蠕變壽命仿真計(jì)算

由于目前航空公司實(shí)際數(shù)據(jù)樣本比較難以獲取，基于某機(jī)隊(duì)歷史飛行數(shù)據(jù)，借助發(fā)動(dòng)機(jī)性能仿真模型仿真得到整機(jī)關(guān)鍵站位的溫度、壓力、流量等性能參數(shù)，進(jìn)一步基于HPT葉片的幾何參數(shù)、材料參數(shù)，構(gòu)建起葉片截面平均溫度和參考點(diǎn)應(yīng)力的數(shù)學(xué)計(jì)算經(jīng)驗(yàn)?zāi)Ｐ?，得到每個(gè)航班飛行任務(wù)中葉片關(guān)鍵部位的應(yīng)力譜和溫度譜，并基于L-M參數(shù)方程計(jì)算每個(gè)航班的蠕變損傷量。

2.1 渦輪葉片蠕變壽命計(jì)算流程

蠕變損傷計(jì)算涉及到溫度和應(yīng)力，需要進(jìn)行相應(yīng)的熱分析和應(yīng)力分析。按照分析維度的不同可劃分為：考慮面平均或體平均的0D分析、考慮直線分布的1D分析、考慮平面分布的2D分析、考慮體分布的3D分析。隨著分析維度的提高，模型中包含的單元類型和數(shù)量越來越多，建模的復(fù)雜程度逐漸提高，求解結(jié)果也會(huì)越來越趨近于真實(shí)情況?？梢姺治龅木S度從根本上決定了分析結(jié)果的準(zhǔn)確性，但單純?cè)黾臃治鼍S度會(huì)迅速延長(zhǎng)模型的求解時(shí)間，模型的輸入也越來越趨向于精細(xì)化，便利性大大降低。參考發(fā)動(dòng)機(jī)部分設(shè)計(jì)公式，葉片關(guān)鍵部位溫度和應(yīng)力載荷歷程計(jì)算用0D和1D的計(jì)算方法進(jìn)行簡(jiǎn)化，在提高計(jì)算速度的同時(shí)，保證了一定的計(jì)算精度。所使用的蠕變損傷集成計(jì)算流程如圖2所示，其中包含了用于中間環(huán)節(jié)計(jì)算和參數(shù)傳遞的子模型，同時(shí)也顯示了不同子模型所需要的輸入數(shù)據(jù)。

從圖2中可以看出，蠕變損傷綜合評(píng)估模型由4個(gè)子模型組成：

1) 性能仿真模型，仿真發(fā)動(dòng)機(jī)HPT葉片進(jìn)出口相關(guān)的氣動(dòng)參數(shù)和熱參數(shù)，為葉片的溫度和應(yīng)力計(jì)算模型提供輸入條件。

2) 葉片溫度計(jì)算模型，對(duì)葉片進(jìn)行溫度分析，并為葉片應(yīng)力計(jì)算模型提供輸入條件。

3) 葉片應(yīng)力計(jì)算模型，計(jì)算葉片參考點(diǎn)的應(yīng)力。

4) 蠕變損傷計(jì)算模型，根據(jù)給定的飛參數(shù)據(jù)計(jì)算葉片蠕變累積損傷。

圖2 蠕變損傷集成計(jì)算流程Fig.2 Integrated calculation process of creep damage

2.2 發(fā)動(dòng)機(jī)性能模型

進(jìn)行葉片的溫度和應(yīng)力的計(jì)算時(shí)需要提供葉片進(jìn)出口的速度、溫度、壓力、流量等參數(shù)作為仿真計(jì)算輸入的邊界條件。這些參數(shù)往往不能由傳感器直接測(cè)量，工程中通常使用一個(gè)能夠反映發(fā)動(dòng)機(jī)整機(jī)性能的模型來近似替代?；诎l(fā)動(dòng)機(jī)性能軟件Turbofan Simulator，通過使用動(dòng)態(tài)鏈接庫(kù)的方法編程調(diào)用，建立能夠匹配實(shí)際發(fā)動(dòng)機(jī)的性能參數(shù)的部件級(jí)性能仿真模型。性能仿真實(shí)際上是對(duì)發(fā)動(dòng)機(jī)的氣路進(jìn)行分析，通過調(diào)整相關(guān)的健康參數(shù)，使得發(fā)動(dòng)機(jī)性能模型的輸出與實(shí)測(cè)的氣路參數(shù)相匹配。其中非線性模型由于精度較高的優(yōu)勢(shì)，最先是由Stamatis等引進(jìn)到工業(yè)燃?xì)廨啓C(jī)的氣路分析當(dāng)中,之后這種方法被推廣到了航空發(fā)動(dòng)機(jī)上，隨著應(yīng)用技術(shù)的成熟，最終開發(fā)出了民航發(fā)動(dòng)機(jī)性能模型商業(yè)軟件Turbofan Simulator。Turbofan Simulator能夠?qū)崿F(xiàn)海拔高度、環(huán)境溫度、飛行馬赫數(shù)、控制變量等運(yùn)行參數(shù)同發(fā)動(dòng)機(jī)整機(jī)各部位的溫度、壓力、流量等性能參數(shù)的快速匹配，同時(shí)可以調(diào)節(jié)表征發(fā)動(dòng)機(jī)氣路部件性能退化的健康指數(shù)，模擬發(fā)動(dòng)機(jī)服役條件下的性能退化行為。性能模型包含了風(fēng)扇(Fan)、低壓壓氣機(jī)(Low Pressure Compressor, LPC)、高壓壓氣機(jī)(High Pressure Compressor, HPC)、燃燒室、高壓渦輪(High Pressure Turbine, HPT)和低壓渦輪(Low Pressure Turbine, LPT)5個(gè)單元體，結(jié)構(gòu)如圖3所示。

性能模型輸入?yún)?shù)包括飛行高度、標(biāo)準(zhǔn)大氣溫差、飛行馬赫數(shù)等；蠕變損傷集成計(jì)算模型需要的性能模型輸出參數(shù)包括高壓轉(zhuǎn)子轉(zhuǎn)速、高壓壓氣機(jī)進(jìn)出口壓力、HPT進(jìn)出口溫度等，詳細(xì)數(shù)據(jù)見表1。

圖3 發(fā)動(dòng)機(jī)性能模型結(jié)構(gòu)示意圖Fig.3 Structural diagram of engine performance model

表1 性能仿真模型輸入和部分輸出參數(shù)

2.3 渦輪葉片應(yīng)力計(jì)算模型

應(yīng)力計(jì)算模型主要計(jì)算葉片參考點(diǎn)的應(yīng)力，根據(jù)對(duì)應(yīng)力成分的分析，將兩種主要的應(yīng)力來源納入應(yīng)力計(jì)算模型當(dāng)中：由發(fā)動(dòng)機(jī)轉(zhuǎn)動(dòng)的離心載荷引起的離心應(yīng)力和由氣體彎矩引起的氣動(dòng)應(yīng)力。

為了較為全面地考核葉片的應(yīng)力分布情況，針對(duì)某型民航發(fā)動(dòng)機(jī)HPT渦輪葉片，首先通過逆向工程獲取其計(jì)算機(jī)輔助設(shè)計(jì)模型，如圖4所示。

圖4 逆向工程得到的HPT葉片三維數(shù)模Fig.4 Three-dimensional model of HPT blade obtained by reverse engineering

進(jìn)一步采用了一種偽二維解析的方法，按葉片高度取渦輪葉片的4個(gè)截面：葉根、25%葉高、50%葉高、75%葉高；根據(jù)葉片可能出現(xiàn)蠕變失效的危險(xiǎn)點(diǎn)位置，每個(gè)葉片截面上分別取前緣最遠(yuǎn)點(diǎn)(Leading Edge, LE)、后緣最遠(yuǎn)點(diǎn)(Trailing Edge, TE)、葉背最遠(yuǎn)點(diǎn)(Blade Back, BB)作為計(jì)算參考點(diǎn)，并使用三維建模軟件對(duì)葉片截面的幾何尺寸進(jìn)行測(cè)繪，截面選取的位置和計(jì)算參考點(diǎn)位置示意如圖5 和圖6所示。

圖5 葉片截面示意圖Fig.5 Schematic diagram of blade section

圖6 葉片截面尺寸示意圖Fig.6 Schematic diagram of blade section size

取葉根截面、50%葉高截面和葉尖截面，將葉片分為“葉根-葉中”和“葉中-葉尖”兩部分，由于本文中的葉片型號(hào)沒有葉冠，因此可近似認(rèn)為葉尖截面的離心應(yīng)力為0、25%葉高截面和75%葉高截面的應(yīng)力可由葉根和葉中截面的計(jì)算結(jié)果進(jìn)行插值得到。

葉根和50%葉高截面的平均離心力和離心應(yīng)力分別為

=

(4)

(5)

式中：為葉片材料密度；為計(jì)算截面和葉尖截面的平均面積；為截面到葉尖的高度；為旋轉(zhuǎn)角速度；為“截面—葉尖”部分的重心到旋轉(zhuǎn)中心的距離；為對(duì)應(yīng)的截面面積。葉根截面和葉中截面的壓力、彎矩分別為

(6)

=∑()

(7)

式中：為截面截取的葉片段繞旋轉(zhuǎn)軸一周在垂直于旋轉(zhuǎn)軸方向投影形成的環(huán)形面積；為葉片段的平均靜壓差；為葉片的個(gè)數(shù)；為葉片段的重心到葉根的距離。

流經(jīng)葉片進(jìn)口和出口的氣流速度變化會(huì)導(dǎo)致氣體沿軸向和切向的動(dòng)量變化，葉根和葉中截面軸向和切向的動(dòng)量力以及動(dòng)量彎矩分別為

(8)

(9)

=∑()

(10)

=∑(Tansec)

(11)

式中：為投影的環(huán)狀面積單位面積的質(zhì)量流量；和分別為軸向和切向速度差。葉片截面的方向彎矩和截面參考點(diǎn)的彎矩應(yīng)力分別為

=(+)sin+cos

(12)

=(+)cos-sin

(13)

(14)

式中：為葉片安裝角；和為葉片段產(chǎn)生的彎矩沿主軸方向的投影分量；和分別為截面參考點(diǎn)到截面慣性主軸的距離；和分別為截面的最大和最小慣性矩。葉片截面尺寸的定義如圖6所示，、、分別代表葉片截面前緣、葉背最遠(yuǎn)點(diǎn)、尾緣到軸的距離；、、分別代表葉片截面前緣、葉背最遠(yuǎn)點(diǎn)、尾緣到軸的距離；為葉片安裝角。因此，截面參考點(diǎn)的總應(yīng)力可表示為

=+

(15)

2.4 渦輪葉片溫度計(jì)算模型

溫度計(jì)算模型主要用來計(jì)算葉片的基體溫度，本文中涉及的葉片型號(hào)主要冷卻方式為氣膜冷卻，采用1D的溫度計(jì)算模型不僅能反映蠕變失效的點(diǎn)位，結(jié)合上文中的應(yīng)力計(jì)算結(jié)果和后續(xù)的蠕變損傷評(píng)估結(jié)果，還能反映不同位置處溫度和應(yīng)力的不同組合對(duì)蠕變損傷的影響。類比于應(yīng)力計(jì)算模型，計(jì)算溫度時(shí)將葉片按照葉片長(zhǎng)度的0～25%、25%～50%、50%～75%、75%～100%依次劃分為4段，從葉根到葉尖，把前一段的冷卻氣體出口溫度作為后一段的冷卻氣體入口溫度，計(jì)算得到葉片葉根、25%葉高、50%葉高、75%葉高共4個(gè)截面的平均溫度。

計(jì)算時(shí)首先需要指定徑向溫度分布因子(Radial Temperature Distribution Factor, RTDF)，RTDF的值一般不超過0.2。同時(shí)在進(jìn)行溫度計(jì)算前，需要做如下假設(shè)：

1) 由于葉片的旋轉(zhuǎn)效應(yīng)，葉片最高溫度的位置會(huì)從葉片中部向葉尖區(qū)域轉(zhuǎn)移(圖7)，因此假設(shè)最大燃?xì)鉁囟任挥诰嚯x葉根75%葉高的位置。

2) 葉根和葉尖的燃?xì)鉁囟葹樽畹腿細(xì)鉁囟取?/p>

3) 從75%葉高到葉根的溫度降低是線性的。

4) 從75%葉高到葉尖的溫度降低是線性的。

最大和最小燃?xì)鉁囟葹?/p>

=+

(16)

(17)

式中：為渦輪進(jìn)口溫度；為燃燒室溫升；為徑向溫度分布因子。葉片段的平均金屬溫度為

=-(-)

(18)

(19)

(20)

(21)

圖7 徑向溫度分布Fig.7 Radial temperature distribution

2.5 蠕變損傷評(píng)估模型

使用上述的溫度、應(yīng)力計(jì)算模型，每次計(jì)算可以得到葉片4個(gè)截面、每個(gè)截面3個(gè)參考點(diǎn)共12個(gè)點(diǎn)的溫度和應(yīng)力情況，由于葉片的葉型尺寸和工作條件決定了葉片的蠕變失效位置，因此蠕變失效最先出現(xiàn)在葉片的哪個(gè)位置并不固定，表2給出了某航班爬升工況下的溫度和應(yīng)力計(jì)算結(jié)果。

從表2中可以看出，截面溫度隨著葉片跨度的增加而逐漸升高，截面應(yīng)力隨著葉片跨度的增加而逐漸減??；在同一截面內(nèi)，前緣參考點(diǎn)應(yīng)力最大，尾緣參考點(diǎn)應(yīng)力次之，葉背參考點(diǎn)應(yīng)力最小，但應(yīng)力水平差距較小。綜合式(1)和式(2)可計(jì)算出表2中溫度和應(yīng)力對(duì)應(yīng)的蠕變壽命，見表3。

表2 溫度和應(yīng)力計(jì)算結(jié)果Table 2 Temperature and stress calculation results

表3 蠕變壽命計(jì)算結(jié)果Table 3 Creep life calculation results

從蠕變壽命計(jì)算結(jié)果中可以看出，50%葉高截面的前緣蠕變壽命最短，可以作為整個(gè)葉片的危險(xiǎn)點(diǎn)。由式(3)可求出所有時(shí)刻對(duì)應(yīng)的單位時(shí)間危險(xiǎn)點(diǎn)蠕變損傷，之后遵循Minner線性法則，對(duì)一個(gè)飛行循環(huán)持續(xù)時(shí)間內(nèi)的蠕變損傷進(jìn)行累加，即可得到某個(gè)飛行循環(huán)對(duì)應(yīng)的累積蠕變損傷：

(22)

式中：為在第個(gè)飛行循環(huán)、第個(gè)載荷狀態(tài)下葉片的蠕變持久壽命；為在第個(gè)飛行循環(huán)、第個(gè)載荷狀態(tài)下的持續(xù)時(shí)間；即為第個(gè)飛行循環(huán)對(duì)應(yīng)的累積蠕變損傷。

3 渦輪葉片蠕變累積損傷指數(shù)模型

在飛機(jī)服役過程中，即使有兩臺(tái)發(fā)動(dòng)機(jī)的飛行任務(wù)剖面完全一致，甚至于渦輪葉片表面的載荷歷程完全一致，但由于實(shí)際的運(yùn)行環(huán)境和操作條件等因素不可能完全相同，最終導(dǎo)致兩臺(tái)發(fā)動(dòng)機(jī)的葉片使用壽命存在較大差異。傳統(tǒng)的基于統(tǒng)計(jì)分布的方法雖然可以擬合出葉片失效樣本整體的壽命分布，但很難基于歷史的環(huán)境載荷對(duì)單機(jī)狀態(tài)進(jìn)行有效的壽命評(píng)估；而對(duì)于非設(shè)計(jì)單位而言，葉片材料數(shù)據(jù)、三維模型等資料難以獲取，使得精確的有限元仿真計(jì)算難以進(jìn)行，其計(jì)算結(jié)果很難反映出實(shí)際情況下的不確定性因素對(duì)壽命的影響。

對(duì)于發(fā)動(dòng)機(jī)來說，飛參是最常見的發(fā)動(dòng)機(jī)動(dòng)態(tài)時(shí)序數(shù)據(jù)，可以從機(jī)載設(shè)備上直接獲取，其次它作為傳感器數(shù)據(jù)，本身就集成了運(yùn)行環(huán)境、運(yùn)行載荷等因素的綜合作用，直接反映了發(fā)動(dòng)機(jī)的使用歷程，體現(xiàn)了單機(jī)運(yùn)行條件。將這種動(dòng)態(tài)數(shù)據(jù)納入分析的范圍也可以更好地解釋發(fā)動(dòng)機(jī)個(gè)體剩余壽命的差異性。由第1部分蠕變損傷評(píng)估模型可知，葉片的溫度和應(yīng)力直接決定了葉片的蠕變損傷，而高壓轉(zhuǎn)子轉(zhuǎn)速直接決定了葉片的離心應(yīng)力；排氣溫度通常與葉片溫度有著較好的相關(guān)性，這兩種時(shí)序數(shù)據(jù)作為常見的狀態(tài)監(jiān)測(cè)數(shù)據(jù)記錄在機(jī)載快速存儲(chǔ)設(shè)備中，非常容易獲取。因此基于動(dòng)態(tài)時(shí)序數(shù)據(jù)和,定義了一種HPT葉片的蠕變累積損傷指數(shù)模型，用于表征使用環(huán)境載荷等參數(shù)對(duì)蠕變壽命分布的影響。通過對(duì)Ⅰ型截尾數(shù)據(jù)累積損傷指數(shù)失效分布相關(guān)參數(shù)的最大似然估計(jì)，得到基于累積損傷指數(shù)的失效分布函數(shù)，并在第4節(jié)對(duì)葉片壽命進(jìn)行評(píng)估。

3.1 累積損傷指數(shù)模型

累積損傷指數(shù)(Cumulative damage index, CDI)描述了動(dòng)態(tài)環(huán)境載荷因素對(duì)觀測(cè)樣本失效時(shí)間分布的影響，環(huán)境載荷歷程稱為樣本的動(dòng)態(tài)協(xié)變量過程。假設(shè)隨機(jī)變量表示樣本的失效時(shí)間，隨機(jī)變量表示截尾數(shù)據(jù)的失效標(biāo)志位(=1表示觀察結(jié)束前樣本已失效；=0表示在觀察結(jié)束時(shí)樣本未失效)，()表示樣本空間內(nèi)的個(gè)體在一段時(shí)間間隔下的協(xié)變量集合。=(,,())為可觀測(cè)到的組合隨機(jī)變量，對(duì)于樣本內(nèi)的第個(gè)樣本個(gè)體，具體數(shù)據(jù)可表示為(,,()),=1,2,…,，其中為樣本總量；表示第個(gè)樣本個(gè)體的失效時(shí)間；表示第個(gè)樣本個(gè)體的失效標(biāo)志位；()表示第個(gè)樣本個(gè)體失效時(shí)間內(nèi)的協(xié)變量過程，則時(shí)間內(nèi)的累積損傷指數(shù)定義為

(23)

圖8 葉片溫度與排氣溫度散點(diǎn)圖Fig.8 Scatter diagram of blade temperature and exhaust temperature

圖9 葉片應(yīng)力與高壓轉(zhuǎn)子轉(zhuǎn)速平方散點(diǎn)圖Fig.9 Scatter diagram of blade stress and high pressure rotor speed

(24)

式中：、為線性回歸的回歸系數(shù)；、為常數(shù)項(xiàng)。忽略、并將式(24)代入式(2)和式(3)中可得：

(25)

忽略式(25)中的次要常數(shù)項(xiàng)，并作進(jìn)一步化簡(jiǎn)

(26)

式中：、、為相關(guān)的系數(shù)。綜合上述分析，使用協(xié)變量和的函數(shù)形式來定義蠕變損傷指數(shù)：

(27)

式中：和是損傷指數(shù)的相關(guān)未知參數(shù)；=8為給定的常數(shù)項(xiàng)。進(jìn)而一段時(shí)間內(nèi)的蠕變累積損傷指數(shù)()可表示為

(28)

在失效時(shí)間內(nèi)，每個(gè)樣本依據(jù)隨機(jī)的協(xié)變量過程(∞)產(chǎn)生對(duì)應(yīng)的累積損傷指數(shù)，當(dāng)其達(dá)到一定的閾值時(shí)樣本失效，假定的累積分布函數(shù)(,)服從威布爾分布，轉(zhuǎn)化為極值分布后滿足：

(29)

(30)

式中：為極值分布函數(shù)；為極值分布概率密度函數(shù)；=(,)′為極值分布的特征參數(shù)。與失效時(shí)間之間的關(guān)系為

(31)

故的分布函數(shù)和概率密度函數(shù)分別為

(;,)=((;,());)

(32)

(;,)=′(();)((;,());)

(33)

3.2 模型參數(shù)估計(jì)方法

渦輪葉片的服役壽命普遍在上萬小時(shí)，實(shí)際當(dāng)中一般統(tǒng)計(jì)截止時(shí)間內(nèi)的飛行循環(huán)數(shù)作為葉片的壽命信息，如果截止時(shí)間內(nèi)葉片失效，則失效標(biāo)志位記為1；否則失效標(biāo)志位記為0。失效標(biāo)志位和截止時(shí)間內(nèi)的葉片壽命信息一同構(gòu)成失效樣本，這類數(shù)據(jù)也稱為截尾數(shù)據(jù)。截尾數(shù)據(jù)中的失效時(shí)間是一個(gè)含有連續(xù)部分和離散部分的混合隨機(jī)變量，綜合式(31)和式(32)可以推導(dǎo)出其對(duì)應(yīng)的似然函數(shù)為

(34)

式中：=(,)′和=(,)′是似然函數(shù)中的待估計(jì)參數(shù)。結(jié)合式(28)、式(31)、式(32)和式(34)，通過迭代求解可以估計(jì)出參數(shù)、、和的值。

3.3 剩余壽命預(yù)測(cè)方法

葉片的剩余壽命是基于當(dāng)前的使用壽命定義的，本身是一個(gè)隨機(jī)變量。若定義葉片剩余壽命分布(Distribution of Remaining Life, DRL)為

(;)=(<≤+|>,())>0

(35)

則第臺(tái)發(fā)動(dòng)機(jī)基于內(nèi)的歷史協(xié)變量過程()=()的條件下，經(jīng)過飛行循環(huán)的剩余壽命分布可表示為

(;)=((<≤+|>,

(),(,+)))

(36)

式中：(,)=(():<<)是～飛行循環(huán)內(nèi)葉片的協(xié)變量歷史；是關(guān)于(,+)|()=()的條件期望。而基于當(dāng)前的歷史協(xié)變信息預(yù)測(cè)未來的協(xié)變過程十分復(fù)雜，且直接得到剩余壽命分布的數(shù)學(xué)表達(dá)式較為困難。對(duì)于一臺(tái)服役狀態(tài)下的某機(jī)隊(duì)的民用航空發(fā)動(dòng)機(jī)而言，一旦投入到航線開始運(yùn)行，通常需要重復(fù)執(zhí)行相同的航班任務(wù)，其歷史飛行航線結(jié)構(gòu)就相對(duì)較為固定，因此可以從歷史航班和對(duì)應(yīng)的協(xié)變量信息中抽取樣本，作為未來飛行航班信息的預(yù)測(cè)，從而使用蒙特卡洛數(shù)值仿真的方法，對(duì)剩余壽命分布進(jìn)行評(píng)估。蒙特卡洛仿真流程如下：

2) 得到更新的協(xié)變量歷史數(shù)據(jù)：

7)重復(fù)過程6)次，得到：

4 仿真驗(yàn)證

要驗(yàn)證提出方法的有效性，不僅需要葉片的壽命統(tǒng)計(jì)數(shù)據(jù)，還需要對(duì)應(yīng)的發(fā)動(dòng)機(jī)工況監(jiān)測(cè)數(shù)據(jù)，目前收集的外場(chǎng)樣本數(shù)據(jù)有限，在某機(jī)隊(duì)航班歷史數(shù)據(jù)基礎(chǔ)上，采用第2節(jié)所述渦輪葉片蠕變壽命仿真方法，通過仿真獲取渦輪葉片的壽命數(shù)據(jù)樣本及其歷史航班的動(dòng)態(tài)協(xié)變量數(shù)據(jù)，驗(yàn)證提出的模型與方法。

4.1 渦輪葉片蠕變壽命周期仿真

收集的某機(jī)隊(duì)的歷史飛行數(shù)據(jù)涵蓋了1～12月的航班數(shù)據(jù)，設(shè)計(jì)轉(zhuǎn)速為14 600 r/min，其中一個(gè)航班的和歷程，如圖10所示。

圖11為渦輪葉片蠕變壽命周期仿真流程，其中FN=1,2,…,100為發(fā)動(dòng)機(jī)編號(hào)。在滿足飛機(jī)起飛安全性的前提下，以相對(duì)于發(fā)動(dòng)機(jī)正常起飛推力較小的推力起飛的方式叫做減推力起飛，也稱靈活推力起飛?，F(xiàn)代高涵道比民航發(fā)動(dòng)機(jī)一般具有較高的推力儲(chǔ)備，因此采用減推力起飛不僅完全可行，這種起飛方式不僅能夠延長(zhǎng)發(fā)動(dòng)機(jī)的使用壽命，還能夠節(jié)省燃油，降低運(yùn)行成本。因

圖10 排氣溫度和高壓轉(zhuǎn)子轉(zhuǎn)速歷程Fig.10 Course of exhaust temperature and high pressure rotor speed history

圖11 蠕變壽命周期仿真Fig.11 Creep life cycle simulation

此考慮到飛機(jī)的減推力起飛，在進(jìn)行仿真時(shí)將發(fā)動(dòng)機(jī)的額定推力等級(jí)考慮在內(nèi)，對(duì)歷史航班按照4個(gè)推力等級(jí)劃分，其燃油流量分別為滿推力的100%、97%、94%和91%，分別對(duì)應(yīng)T0,T1,T2,T3起飛階段減推力的4個(gè)等級(jí)；同時(shí)考慮到壽命評(píng)估的不確定因素，可將不確定性來源總體上劃分為三類：一是葉片設(shè)計(jì)階段的材料參數(shù)分布、計(jì)算機(jī)輔助設(shè)計(jì)誤差等引起的不確定性；二是葉片制造階段的材料晶體缺陷、驗(yàn)收測(cè)量精度等引起的不確定性；三是葉片使用階段的自然環(huán)境載荷、發(fā)動(dòng)機(jī)性能退化等引起的不確定性。這些不確定性相互作用，最終導(dǎo)致了服役葉片的實(shí)際可用時(shí)間與設(shè)計(jì)時(shí)間不符。考慮到以上不確定性來源的影響，使用一個(gè)隨機(jī)變量～(0,)代表這些因素的綜合作用，其中為隨機(jī)變量方差，由于隨機(jī)變量的均值影響的是全壽命樣本整體飛行循環(huán)的偏移，這種數(shù)據(jù)偏移的影響理論上可以通過分布的標(biāo)準(zhǔn)化過程進(jìn)行消除，因此為了方便后續(xù)處理和分析，隨機(jī)變量的均值取為0。壽命周期仿真時(shí)，每次從相應(yīng)的樣本空間隨機(jī)抽取航班，記錄下飛參數(shù)據(jù)和計(jì)算得到的蠕變損傷，則不確定性因素影響下的航班損傷為′=(1+)，將′作為此次抽取的實(shí)際損傷進(jìn)行累積，當(dāng)損傷累積到1時(shí)停止抽樣，每進(jìn)行1次這樣的抽樣循環(huán)便得到了一臺(tái)發(fā)動(dòng)機(jī)的全壽命樣本，按照上述方法，每個(gè)推力等級(jí)下各仿真25臺(tái)全壽命樣本,4個(gè)推力等級(jí)共包含100個(gè)全壽命樣本，對(duì)應(yīng)的壽命分布如圖12所示。

其中1臺(tái)發(fā)動(dòng)機(jī)全壽命周期的蠕變損傷量分布，如圖13所示，可以看到發(fā)動(dòng)機(jī)的單個(gè)航班蠕變損傷量大約分布于10～10之間，同一臺(tái)發(fā)動(dòng)機(jī)不同航班的蠕變損傷呈現(xiàn)較強(qiáng)的分散性，損傷量差距高達(dá)10倍，這主要是不同航班的運(yùn)行載荷不同導(dǎo)致的結(jié)果。

圖12 HPT葉片壽命分布直方圖Fig.12 Distribution histogram of HPT blade lifetime

取其中一個(gè)高損傷與一個(gè)低損傷航班，并分別統(tǒng)計(jì)起飛爬升階段的數(shù)據(jù)，從圖14和圖15中可以看出，蠕變損傷較大的航班對(duì)應(yīng)的相對(duì)排氣溫度和相對(duì)高壓轉(zhuǎn)子轉(zhuǎn)速水平普遍高于蠕變損傷較小的航班，說明服役條件下的載荷水平差異是引起不同航班蠕變損傷差異性的重要因素。

分別統(tǒng)計(jì)1臺(tái)發(fā)動(dòng)機(jī)主要飛行階段的累積損傷以及相應(yīng)的持續(xù)時(shí)間占比，統(tǒng)計(jì)結(jié)果如圖16所示。5階段為飛機(jī)爬升階段，葉片處于高負(fù)載狀態(tài)，此階段持續(xù)時(shí)間所占比重較小，約為總時(shí)間的12.5%左右，但累積蠕變損傷卻占總損傷的80%以上；6階段為飛機(jī)巡航階段，葉片處于低負(fù)載狀態(tài)，此階段的時(shí)間占比達(dá)到了50%以上，但累積損傷占總損傷不到10%，說明葉片在高溫和高應(yīng)力水平下的單位時(shí)間蠕變損傷效應(yīng)已經(jīng)超過了飛行歷程內(nèi)蠕變損傷的時(shí)間累積效應(yīng)。

圖13 航班蠕變損傷統(tǒng)計(jì)情況Fig.13 Statistics of flight creep damage

圖14 起飛爬升階段相對(duì)排氣溫度歷程Fig.14 Relative exhaust temperature course during takeoff and climb phase

圖15 起飛爬升階段相對(duì)高壓轉(zhuǎn)子轉(zhuǎn)速歷程Fig.15 Relative high pressure rotor speed course during takeoff and climb phase

圖16 飛行階段損傷統(tǒng)計(jì)Fig.16 Damage statistics during flight phase

4.2 基于蠕變累積損傷指數(shù)的可靠性壽命評(píng)估

仿真得到的100臺(tái)發(fā)動(dòng)機(jī)樣本中，隨機(jī)抽取50臺(tái)作為失效樣本，這部分發(fā)動(dòng)機(jī)保留壽命周期的飛行航班數(shù)據(jù)和對(duì)應(yīng)的動(dòng)態(tài)協(xié)變量信息，失效標(biāo)志位記為=1(=1,2,…,50)；剩余50臺(tái)作為未失效樣本，對(duì)飛行循環(huán)進(jìn)行隨機(jī)截尾處理，記錄截尾時(shí)刻之前的航班及對(duì)應(yīng)的協(xié)變量信息，失效標(biāo)志位記為=0(=51,52,…,100)，以此模擬機(jī)隊(duì)發(fā)動(dòng)機(jī)葉片實(shí)際使用壽命的統(tǒng)計(jì)情況，其中前20臺(tái)發(fā)動(dòng)機(jī)截尾數(shù)據(jù)分布如圖17所示。

圖17 發(fā)動(dòng)機(jī)截尾數(shù)據(jù)Fig.17 Censored engine data

使用R語言編寫似然函數(shù)，并調(diào)用optim優(yōu)化函數(shù)進(jìn)行估計(jì)參數(shù)的迭代求解，其中，優(yōu)化方法設(shè)置為BFGS，BFGS方法屬于Quasi-Newton方法，其在最優(yōu)值附近的計(jì)算收斂速度快，迭代次數(shù)少。經(jīng)過546次迭代計(jì)算，得到式(28)中的參數(shù)估計(jì)值：

將、反代可以求出每個(gè)航班的累積損傷指數(shù)。將每臺(tái)發(fā)動(dòng)機(jī)所有航班累積損傷指數(shù)累加，即可得到每臺(tái)發(fā)動(dòng)機(jī)截尾時(shí)刻對(duì)應(yīng)的總損傷指數(shù)，計(jì)算結(jié)果的統(tǒng)計(jì)情況如圖18和圖19所示。

由圖18和圖19可見：① 失效樣本的累積損傷指數(shù)普遍大于未失效樣本的累積損傷指數(shù)，且兩種損傷指數(shù)分界較為明顯，可以很好地區(qū)分失效樣本和未失效樣本，說明本文定義的蠕變累積損傷指數(shù)可以很好的表征HPT葉片的實(shí)際損傷狀態(tài)。② 蠕變累積損傷指數(shù)與飛行循環(huán)數(shù)沒有明顯的強(qiáng)相關(guān)性，僅僅借助飛行循環(huán)數(shù)或飛行時(shí)

圖18 蠕變累積損傷指數(shù)散點(diǎn)圖Fig.18 Scatter plot of creep cumulative damage index

圖19 蠕變累積損傷指數(shù)直方圖Fig.19 Histogram of creep cumulative damage index

間來表征HPT葉片的累積損傷并不合適，低飛行循環(huán)數(shù)累積損傷指數(shù)可能較大，葉片已經(jīng)發(fā)生蠕變失效，而高飛行循環(huán)數(shù)累積損傷指數(shù)可能較小，葉片仍可以正常服役。這是因?yàn)镠PT葉片蠕變累積損傷的不僅僅取決于飛行時(shí)間，還取決于經(jīng)歷的溫度、應(yīng)力等環(huán)境載荷參數(shù)。③ 失效樣本的累積損傷指數(shù)呈一定的分散性，失效閾值并不是一個(gè)確定值，這是因?yàn)榧词雇恍吞?hào)的HPT葉片，不同個(gè)體在“設(shè)計(jì)—制造—使用”中不可避免地引進(jìn)不確定性，這些不確定性的相互作用正是導(dǎo)致服役條件下失效葉片累積損傷指數(shù)發(fā)散的重要原因。

式(30)中的估計(jì)參數(shù)值和95%置信區(qū)間邊界值見表4，其中為極值分布的均值，ln()為極值分布標(biāo)準(zhǔn)差以e為底數(shù)的對(duì)數(shù)值。從而得到基于累積損傷指數(shù)的失效概率圖(圖20)，可以看到，累積損傷指數(shù)基本符合威布爾分布。

對(duì)50個(gè)未失效樣本的未來協(xié)變過程進(jìn)行蒙特卡洛仿真，仿真設(shè)置為步長(zhǎng)=80 000，單輪仿真次數(shù)=1 000，仿真輪次=100。仿真計(jì)算得到的其中四臺(tái)發(fā)動(dòng)機(jī)葉片剩余壽命可靠性分布如圖21所示。

表4 蠕變壽命計(jì)算結(jié)果Table 4 Creep life calculation results

圖20 隨累積損傷指數(shù)變化的失效概率圖Fig.20 Diagram of failure probability changing with cumulative damage index

圖21 剩余壽命可靠性分布Fig.21 Distribution of remaining life reliability

上述分布表示了從統(tǒng)計(jì)截止時(shí)間起，未來葉片累積失效率隨飛行循環(huán)數(shù)的變化過程。由于協(xié)變量歷史信息的不同，圖中四臺(tái)發(fā)動(dòng)機(jī)剩余壽命可靠性分布具有明顯差異，這也體現(xiàn)了實(shí)際服役條件下葉片的個(gè)體差異性。針對(duì)剩余壽命分布，工程上常用式(37)所示的平均剩余壽命(Mean Residual Life, MRL)作為剩余壽命預(yù)測(cè)值。

(37)

式中：()為可靠度函數(shù)；為當(dāng)前截止時(shí)間。

在得到了葉片的壽命分布函數(shù)的估計(jì)參數(shù)后，可以根據(jù)分布類型推導(dǎo)出對(duì)應(yīng)的可靠度分布函數(shù)，之后代入式(37)中即可計(jì)算出與分布對(duì)應(yīng)的平均剩余壽命。一般情況下通過直接統(tǒng)計(jì)截尾數(shù)據(jù)的方式也可以得到使用壽命的分布函數(shù)。把這種傳統(tǒng)的處理方法記為基準(zhǔn)模型，分別使用本文定義的累積損傷指數(shù)模型與基準(zhǔn)模型按照式(37) 進(jìn)行平均剩余壽命計(jì)算，部分結(jié)果如表5所示，其中使用壽命指的是截止時(shí)間的實(shí)際使用壽命，剩余壽命是在已知全壽命的基礎(chǔ)上去除使用壽命后剩余的飛行循環(huán)數(shù)，這兩部分?jǐn)?shù)據(jù)對(duì)應(yīng)著全壽命仿真中模擬的未失效樣本。其中30臺(tái)發(fā)動(dòng)機(jī)的預(yù)測(cè)偏差如圖22所示。

結(jié)合表5與圖22中的預(yù)測(cè)結(jié)果可以看出，基于失效時(shí)間分布的剩余壽命預(yù)測(cè)偏差大約在-30%～1 200%的范圍內(nèi)，說明預(yù)測(cè)誤差具有較大的分散性，基本不具備參考價(jià)值；使用累積損傷指數(shù)進(jìn)行的預(yù)測(cè)偏差大約在-20%～20%之間分

表5 平均剩余壽命預(yù)測(cè)結(jié)果對(duì)比Table 5 Comparison of MRL prediction results

圖22 預(yù)測(cè)結(jié)果偏差Fig.22 Forecast result deviation

布，說明預(yù)測(cè)誤差分散性小，結(jié)果總體在可接受范圍內(nèi)，可作為維修計(jì)劃定制的參考剩余壽命。

5 結(jié)論及展望

1) 針對(duì)民航發(fā)動(dòng)機(jī)渦輪葉片的剩余壽命預(yù)測(cè)，提出了一種數(shù)據(jù)驅(qū)動(dòng)與失效物理融合的可靠性剩余壽命預(yù)測(cè)方法，基于服役條件下的歷史工況監(jiān)測(cè)數(shù)據(jù)，結(jié)合發(fā)動(dòng)機(jī)性能仿真模型、葉片應(yīng)力和溫度計(jì)算模型對(duì)葉片表面參考點(diǎn)的溫度和應(yīng)力歷程進(jìn)行計(jì)算；在此基礎(chǔ)上使用拉森米勒參數(shù)法對(duì)葉片蠕變損傷進(jìn)行評(píng)估；之后針對(duì)蠕變損傷定義了一種累積損傷指數(shù)模型，融合歷史協(xié)變量信息排氣溫度和高壓轉(zhuǎn)子轉(zhuǎn)速對(duì)葉片剩余壽命進(jìn)行了可靠性評(píng)估，并通過數(shù)據(jù)仿真分析驗(yàn)證了本文建立模型的可行性，為渦輪葉片可靠性評(píng)估提供了一種新的思路。

2) 針對(duì)服役條件下的HPT葉片蠕變壽命預(yù)測(cè)，使用本文的模型方法得到的預(yù)測(cè)結(jié)果，預(yù)測(cè)偏差大約在±20%范圍內(nèi)分布，而傳統(tǒng)的基于失效時(shí)間分布的可靠性分析方法，得到的預(yù)測(cè)偏差約在-30%～1 200%之間，說明提出的模型具有更好的適用性，在壽命預(yù)測(cè)精度上也更具有優(yōu)勢(shì)。

3) 反映發(fā)動(dòng)機(jī)關(guān)鍵件服役可靠性降級(jí)等行為特性的信息存在于多模態(tài)數(shù)據(jù)中，包括產(chǎn)品設(shè)計(jì)數(shù)據(jù)與使用維護(hù)數(shù)據(jù)，如產(chǎn)品幾何參數(shù)、材料數(shù)據(jù)、狀態(tài)監(jiān)測(cè)數(shù)據(jù)、離線檢測(cè)/檢查報(bào)告、故障記錄和大修報(bào)告等，貫穿產(chǎn)品壽命周期的不同階段，如何充分利用這些多模態(tài)數(shù)據(jù)，實(shí)現(xiàn)多領(lǐng)域物理模型與數(shù)據(jù)融合，是提高服役發(fā)動(dòng)機(jī)關(guān)鍵件壽命預(yù)測(cè)分析準(zhǔn)確性的重要途徑。以民航發(fā)動(dòng)機(jī)HPT葉片蠕變失效為例，在蠕變失效機(jī)理分析基礎(chǔ)上，同時(shí)考慮服役條件下的使用信息、狀態(tài)參數(shù)及截尾失效時(shí)間等數(shù)據(jù)，提出的蠕變累積損傷指數(shù)模型實(shí)現(xiàn)了多模態(tài)信息的融合，對(duì)渦輪葉片服役可靠性進(jìn)行評(píng)估，并預(yù)測(cè)其剩余壽命。

4) 本文中涉及的葉片損傷正向計(jì)算中，每個(gè)計(jì)算環(huán)節(jié)都會(huì)不可避免地引進(jìn)計(jì)算的不確定性，這些不確定性既有認(rèn)知的不確定性，也有材料、設(shè)計(jì)等固有的不確定性，這是正向計(jì)算所不可避免的。因此“數(shù)據(jù)驅(qū)動(dòng)與失效物理融合”是對(duì)可靠性壽命評(píng)估方面創(chuàng)新性的思考，它主要體現(xiàn)在累積損傷指數(shù)模型的建立過程以及外場(chǎng)使用數(shù)據(jù)、截尾壽命數(shù)據(jù)的融合過程中。損傷指數(shù)定義的方式基于一系列的正向的計(jì)算基礎(chǔ)上的推導(dǎo)，而發(fā)動(dòng)機(jī)性能模型、溫度計(jì)算模型、應(yīng)力計(jì)算模型和失效模式對(duì)應(yīng)的損傷計(jì)算物理模型都與失效物理相對(duì)應(yīng)；之后通過分析蠕變失效機(jī)理以及物理模型，推導(dǎo)了一種新的蠕變損傷累積指數(shù)模型，在此初始模型的基礎(chǔ)上進(jìn)一步建立似然函數(shù)，融合外場(chǎng)使用數(shù)據(jù)、截尾壽命數(shù)據(jù)來估計(jì)模型參數(shù)，從而實(shí)現(xiàn)了在失效物理機(jī)理基礎(chǔ)上對(duì)外場(chǎng)數(shù)據(jù)的融合。模型直接用到的是飛行監(jiān)測(cè)數(shù)據(jù)與外場(chǎng)失效時(shí)間數(shù)據(jù)，而由于累積損傷模型基于失效模式定義，因此預(yù)測(cè)結(jié)果在物理層面的解釋性更好。

現(xiàn)階段外場(chǎng)數(shù)據(jù)獲取有限，因此通過仿真數(shù)據(jù)驗(yàn)證了提出方法的可行性和有效性。下一步的工作將基于目前的研究繼續(xù)深入，從以下幾個(gè)方面開展：一是從熱機(jī)械疲勞、高溫氧化等其他失效機(jī)理入手建立對(duì)應(yīng)的累積損傷指數(shù)模型，驗(yàn)證本文方法的通用性和魯棒性；二是研究有效的動(dòng)態(tài)協(xié)變量預(yù)測(cè)方法，來代替目前使用的蒙特卡洛仿真過程，以取得更好的參數(shù)估計(jì)結(jié)果和預(yù)測(cè)精度；三是跟蹤預(yù)測(cè)部件的實(shí)時(shí)狀態(tài)，對(duì)本文的模型方法進(jìn)行實(shí)際數(shù)據(jù)的驗(yàn)證。