附件機(jī)匣振動(dòng)力學(xué)行為及壽命分析

殷 緣 ，楊 勱 ，信 琦 ，陳云霞

（1.北京航空航天大學(xué)可靠性與系統(tǒng)工程學(xué)院，北京100191；2.中國船舶工業(yè)系統(tǒng)研究院，北京100094；3.中國航發(fā)沈陽發(fā)動(dòng)機(jī)研究所，沈陽110015）

0 引言

航空發(fā)動(dòng)機(jī)附件機(jī)匣是航空發(fā)動(dòng)機(jī)的重要部件，對其力學(xué)行為和壽命進(jìn)行分析以獲得可靠性和壽命指標(biāo)結(jié)果至關(guān)重要。附件機(jī)匣在工作時(shí)，不僅受到軸承載荷、自身和安裝附件重力的影響，同時(shí)由于齒輪系統(tǒng)高速轉(zhuǎn)動(dòng)產(chǎn)熱以及噴油散熱使得其始終處于較復(fù)雜的溫度場中[1]，且由于附件機(jī)匣內(nèi)部齒輪系統(tǒng)嚙合傳動(dòng)以及在安裝、運(yùn)行過程中的誤差和故障等因素導(dǎo)致其處于復(fù)雜的振動(dòng)環(huán)境中[2]?？紤]到附件機(jī)匣工作環(huán)境的復(fù)雜多樣性，以及同時(shí)承受機(jī)械載荷、振動(dòng)載荷以及熱載荷，李錦花等結(jié)合MASTA軟件和ANSYS軟件計(jì)算同時(shí)受到軸承載荷和傳動(dòng)附件安裝彎矩影響下殼體的變形[3]；郭梅等將輪齒嚙合時(shí)變剛度和靜傳遞誤差產(chǎn)生的嚙合力作為系統(tǒng)振動(dòng)激勵(lì)，采用有限元方法分析附件機(jī)匣系統(tǒng)的振動(dòng)響應(yīng)[4]；吳鴻等采用有限元方法建立附件機(jī)匣殼體穩(wěn)態(tài)熱分析模型，給出殼體的溫度場分布[5]。然而，上述研究僅包含機(jī)械載荷、振動(dòng)或者溫度場中的某一方面，不能全面反映其在復(fù)雜工作環(huán)境和熱固條件下的振動(dòng)載荷情況。熱、固、振載荷會(huì)影響附件機(jī)匣應(yīng)力分布等力學(xué)行為，從而影響其疲勞壽命。因此，在熱、固、振條件共同作用下進(jìn)行附件機(jī)匣的力學(xué)行為以及壽命建模分析是十分重要的。

本文基于有限元軟件ANSYS Workbench，同時(shí)考慮附件機(jī)匣殼體自身重力和固定約束條件、軸承載荷、溫度場以及振動(dòng)載荷作用，在不同工況下進(jìn)行附件機(jī)匣殼體力學(xué)行為分析，計(jì)算得到附件機(jī)匣應(yīng)力響應(yīng)PSD譜，并根據(jù)應(yīng)力響應(yīng)PSD譜采用雨流循環(huán)計(jì)數(shù)方法對附件機(jī)匣壽命指標(biāo)進(jìn)行估計(jì)。

1 力學(xué)行為建模及分析

在熱、固、振多載荷耦合作用下，運(yùn)用ANSYS Workbench對其力學(xué)行為進(jìn)行建模及分析。

1.1 模型定義、邊界以及輸入條件

1.1.1 模型定義

在實(shí)際工作環(huán)境中，各種負(fù)載懸掛在附件機(jī)匣兩側(cè)，對其力學(xué)行為分析產(chǎn)生影響。因此，在CAD模型上添加質(zhì)點(diǎn)塊(point mass)來定義負(fù)載對附件機(jī)匣的慣性影響，如圖1所示。

圖1 模型定義

圖1定義模型中各負(fù)載的質(zhì)量以及在局部坐標(biāo)系中的位置見表1。

1.1.2 邊界條件

進(jìn)行ANSYS有限元仿真時(shí)，邊界條件指運(yùn)動(dòng)邊界上方程組的解應(yīng)該滿足的條件。針對附件機(jī)匣的邊界條件為吊耳處固定約束。

表1 負(fù)載附件的質(zhì)量以及坐標(biāo)

1.1.3 輸入條件

進(jìn)行附件機(jī)匣熱、固、振耦合作用下的力學(xué)行為分析，涉及機(jī)械載荷、振動(dòng)載荷和熱，因此輸入條件有：機(jī)匣重力、溫度場、軸承載荷、PSD譜。

根據(jù)ANSYS FLUENT產(chǎn)熱分析，得到附件機(jī)匣殼體上各部位溫度，并采用EXCEL擬合得到附件機(jī)匣溫度分布函數(shù)。環(huán)境溫度設(shè)置為90℃，在附件機(jī)匣中部建立局部坐標(biāo)系，以該坐標(biāo)系為基準(zhǔn)添加溫度，溫度變化圖線如圖2所示。

圖2 附件機(jī)匣不同部位的溫度

因此，溫度變化函數(shù)為

式中：x為附件機(jī)匣與坐標(biāo)原點(diǎn)的相對位置；y為該部位的溫度值。

圖3 附件機(jī)匣溫度場

加載溫度場到附件機(jī)匣上，如圖3所示。附件機(jī)匣內(nèi)部共有10個(gè)齒輪軸，每個(gè)軸兩端各有1個(gè)軸承，僅承受徑向力。采用LMS Virtual Lab.Motion進(jìn)行整體仿真，選取每個(gè)軸承的最大徑向力并加載到附件機(jī)匣上，具體數(shù)值見表2。

表2 軸承載荷

隨機(jī)振動(dòng)輸入的加度素功率譜頻率范圍為0～2000 Hz，各頻段振動(dòng)量值見表3。

表3 該工況下每頻段振動(dòng)量值

1.2 網(wǎng)格劃分

附件機(jī)匣殼體幾何機(jī)構(gòu)較為復(fù)雜，建議采用快速自動(dòng)生成的四面體網(wǎng)格進(jìn)行劃分。此外，在劃分網(wǎng)格時(shí)選擇合適的網(wǎng)格大小對求解的精度以及求解速度十分重要。

結(jié)合附件機(jī)匣實(shí)際尺寸，使用全局8 mm、全局16 mm、附件機(jī)匣兩側(cè)(平板處)32 mm其余8 mm、兩側(cè)16 mm其余8 mm 4種尺寸劃分網(wǎng)格，結(jié)果見表4。

表4 網(wǎng)格劃分結(jié)果相關(guān)參數(shù)

其中，畸變度(skewness)為網(wǎng)格質(zhì)量的重要參數(shù)，指單元相對其理想形狀的相對扭曲的度量，是1個(gè)值在0(極好的)到1(無法接受的)之間的比例因子。根據(jù)表4可知，不同網(wǎng)格尺寸得到的畸變度最大值和平均值基本接近；節(jié)點(diǎn)數(shù)和單元數(shù)為尺寸8 mm的最少，選擇較少節(jié)點(diǎn)和單元的尺寸可以提高運(yùn)算速度。

對于8 mm所得的單元和節(jié)點(diǎn)數(shù)比16 mm得到的少，其原因在于8 mm時(shí)的平滑參數(shù)(transition)設(shè)置為fast，而16 mm設(shè)置為slow。若8 mm網(wǎng)格尺寸設(shè)置為slow則會(huì)輸出過多單元，導(dǎo)致計(jì)算量太大而不能進(jìn)行模態(tài)分析；若16 mm設(shè)置為fast則網(wǎng)格質(zhì)量太低不足以比較。而網(wǎng)格尺寸為8 mm時(shí)畸變度為0.95的單元所占比例最小，因此8 mm的網(wǎng)格尺寸最理想，其對應(yīng)的機(jī)匣局部網(wǎng)格劃分結(jié)果如圖4所示。

圖4 附件機(jī)匣網(wǎng)格

1.3 預(yù)應(yīng)力模態(tài)分析

模態(tài)分析用來計(jì)算結(jié)構(gòu)的振動(dòng)特性，包括固有頻率和振型。由于ANSYS隨機(jī)振動(dòng)分析采用模態(tài)疊加法，因此隨機(jī)振動(dòng)分析前必須進(jìn)行模態(tài)分析。本文中靜力學(xué)分析設(shè)置的邊界和輸入條件使結(jié)構(gòu)中存在預(yù)應(yīng)力，會(huì)導(dǎo)致結(jié)構(gòu)剛度變化進(jìn)而影響模態(tài)分析及后續(xù)振動(dòng)分析結(jié)果。因此，在進(jìn)行預(yù)應(yīng)力模態(tài)分析時(shí)需考慮機(jī)匣載荷和溫度對結(jié)構(gòu)剛度的影響。

設(shè)置PSD譜最高頻率為2000 Hz，進(jìn)行15階模態(tài)分析，固有頻率結(jié)果皆小于500 Hz。因此，為保證響應(yīng)結(jié)果的精確性分析，在前45階模態(tài)中，第1、6、15、45階固有頻率結(jié)果見表5。

1.4 隨機(jī)振動(dòng)響應(yīng)輸出

根據(jù)仿真計(jì)算的薄弱位置的應(yīng)力響應(yīng)PSD譜進(jìn)行壽命分析。根據(jù)隨機(jī)振動(dòng)理論，響應(yīng)PSD譜可以由單自由度系統(tǒng)的傳遞函數(shù)H（ω）根據(jù)輸入（激勵(lì)）PSD譜通過模態(tài)疊加方法計(jì)算得到。第i個(gè)自由度的響應(yīng)PSD譜由以下3部分組成：

動(dòng)態(tài)部分

擬靜態(tài)部分

協(xié)方差部分

式中：n為振型數(shù)；r1、r2分別為節(jié)點(diǎn)數(shù)和基礎(chǔ)PSD譜數(shù)。

單自由度的傳遞函數(shù)根據(jù)不同的輸入PSD譜和響應(yīng)類型有不同形式。本文中輸入為加速度PSD譜，輸出為應(yīng)力PSD譜，其傳遞函數(shù)形式為

式中：ω為外激勵(lì)頻率；ωj為第j階模態(tài)的自然圓弧頻率。

采用ANSYS Workbench得到最大等效應(yīng)力結(jié)果見表6。

表6 最大等效應(yīng)力計(jì)算結(jié)果

從表6中可知，振動(dòng)譜沿Z軸方向等效應(yīng)力最大，如圖5所示。

圖5 等效應(yīng)力

從圖中可見附件機(jī)匣殼體的應(yīng)力最大點(diǎn)（即薄弱位置），從而能夠輸出薄弱位置的正應(yīng)力σx、σy、σz和剪應(yīng)力 τxy、τyz、τzz的響應(yīng) PSD 譜。

應(yīng)力響應(yīng)PSD譜如圖6（a）所示。已知功率譜密度值-頻率值關(guān)系曲線下面積的開方即為均方根值（RMS）[6]。從圖 6（a）中可見，σx響應(yīng) PSD譜的 RMS皆遠(yuǎn)大于其余應(yīng)力，并且各應(yīng)力頻率主要集中于90～150 Hz，將該頻率區(qū)間圖像放大，如如圖 6（b）所示，發(fā)現(xiàn)σz和σyz的RMS遠(yuǎn)小于其余應(yīng)力。

圖6 應(yīng)力響應(yīng)PSD譜

2 壽命建模及分析

2.1 雨流循環(huán)計(jì)數(shù)方法

根據(jù)第1.4節(jié)輸出的響應(yīng)PSD譜，采用雨流循環(huán)計(jì)數(shù)法（RFC）計(jì)算壽命。RFC是循環(huán)計(jì)數(shù)法的1種。循環(huán)技術(shù)法指在已知峰值概率密度函數(shù)的情況下，對應(yīng)力峰值進(jìn)行循環(huán)計(jì)數(shù)得到幅值信息，然后進(jìn)行壽命計(jì)算。常用的循環(huán)計(jì)數(shù)法還有變程計(jì)數(shù)法（RC）和水平穿越計(jì)數(shù)法（NB）[7]。

雨流計(jì)數(shù)法是公認(rèn)最好的疲勞損傷估計(jì)方法，但其計(jì)算過程比較繁瑣，準(zhǔn)確解析式很難給出[8-9]。然而，在任意平穩(wěn)高斯過程中，線性損傷準(zhǔn)則有效的情況下，雨流計(jì)數(shù)損傷總是處于2種損傷值之間[10-11]

式中為雨流計(jì)數(shù)損傷；DRC為變程計(jì)數(shù)損傷；DNB為水平穿越計(jì)數(shù)損傷。

從式（6）中可知，估計(jì)雨流計(jì)數(shù)損傷可以通過在其上下界中找到1個(gè)合適的中間點(diǎn)。因此，D.Benasciutti和R.Tovo提出雨流計(jì)數(shù)損傷可以由和加權(quán)線性組合得到[12]

變程計(jì)數(shù)法和水平穿越計(jì)數(shù)法[13]計(jì)算疲勞損傷的近似公式為

式中：υ0為譜密度函數(shù)Sx（ω）的平均上跨率

式中：λm為單邊譜密度函數(shù)Sx（ω）的譜距

方差由 λ0決定

C和k為S-N曲線中的參數(shù)

權(quán)重b的選擇取決于響應(yīng)的譜密度函數(shù)。D.Benasciutti和R.Tovo給出b的表達(dá)式為

其中

對于譜密度函數(shù) Sx（ω），α1/α2為其帶寬參數(shù)，取值為[0,1][14]。

2.2 線性疲勞累計(jì)損傷理論

根據(jù)線性疲勞累積損傷理論，計(jì)算附件機(jī)匣的總損傷量從而給出其疲勞壽命[15]。

式中：D為總損傷量；Di為單個(gè)循環(huán)造成的損傷。

線性疲勞累計(jì)損傷理論假定損傷量D=1時(shí)試件將發(fā)生疲勞破壞。因此，疲勞壽命為

2.3 疲勞壽命計(jì)算

根據(jù)得到的響應(yīng)PSD譜（圖6），運(yùn)用雨流循環(huán)計(jì)數(shù)法并采用MATLAB編程計(jì)算附件機(jī)匣的疲勞壽命。由于σx響應(yīng)PSD譜RMS遠(yuǎn)大于其余應(yīng)力，根據(jù)σx響應(yīng)PSD譜計(jì)算附件機(jī)匣的損傷以及疲勞壽命，結(jié)果見表7。

根據(jù)6個(gè)應(yīng)力的響應(yīng)PSD譜計(jì)算其單位時(shí)間內(nèi)的損傷，再疊加計(jì)算所有功率譜造成的單位時(shí)間損傷，從而得到附件機(jī)匣總損傷以及壽命結(jié)果，見表8。

表7 損傷及疲勞壽命計(jì)算結(jié)果

從表7、8中可知，2種計(jì)算方式的結(jié)果非常接近，表明 RMS非常小的 σy、σz，τxy、τyz、τzz對附件機(jī)匣的損傷以及壽命計(jì)算結(jié)果影響很小。

3 結(jié)論

本文綜合考慮附件機(jī)匣復(fù)雜的工作環(huán)境，包括自身重力和固定約束條件、軸承載荷、溫度場以及振動(dòng)載荷，分析計(jì)算附件機(jī)匣殼體在機(jī)械載荷、溫度并施加隨機(jī)振動(dòng)譜得到的應(yīng)力響應(yīng)PSD譜，結(jié)合雨流循環(huán)計(jì)數(shù)方法計(jì)算附件機(jī)匣殼體的疲勞壽命。在壽命計(jì)算時(shí)，對比單獨(dú)采用RMS最大的應(yīng)力PSD譜和采用所有應(yīng)力PSD譜計(jì)算得到的壽命結(jié)果分別為3.25828×106min和3.25823×106min，相差僅為50 min。因此，可以選用RMS最大的PSD譜計(jì)算附件機(jī)匣的疲勞壽命。

[1]呂亞國，劉振俠，路彬，等.航空發(fā)動(dòng)機(jī)附件機(jī)匣熱分析研究[J].潤滑與密封,2011,36（10）：62-80.LU Yaguo,LIU Zhenxia,LU Bin,et al.Thermal analysis of aeroengine accessory gearbox[J].Lubrication Engineering,2011，36（10）：62-80.(in Chinese)

[2]王桂華，劉海年.航空發(fā)動(dòng)機(jī)成附件振動(dòng)環(huán)境試驗(yàn)剖面確定方法研究[J].推進(jìn)技術(shù)，2013，34（8）：1101-1107.WANG Guihua,LIU Hainian.Study on formulating method for vibration environment test profiles of aeroengine accessories[J].Journal of Propulsion Technology，2013，34（8）：1101-1107.(in Chinese)

[3]李錦花，史妍妍，張茂強(qiáng)，等.航空發(fā)動(dòng)機(jī)附件機(jī)匣殼體變形分析[J].航空發(fā)動(dòng)機(jī)，2013，39（3）：59-72.LI Jinhua,SHI Yanyan,ZHANG Maoqiang,et al.Analysis of accessory gearbox housing distortion for aeroengine[J].Aeroengine，2013，39（3）：59-72.(in Chinese)

[4]郭梅，陳聰慧，王建軍，等.發(fā)動(dòng)機(jī)附件機(jī)匣結(jié)構(gòu)系統(tǒng)振動(dòng)特性[J].航空動(dòng)力學(xué)報(bào)，2013，28（7）：1607-1612.GUO Mei,CHEN Chonghui,WANG Jianjun.Vibration characteristics of accessory gearbox structure system of engine [J].Journal of Aerospace Power，2013，28（7）：1607-1612.(in Chinese)

[5]吳鴻，李國權(quán)，齊樂華.基于ANSYS的附件機(jī)匣抽油口潤滑油溫度計(jì)算[J].航空發(fā)動(dòng)機(jī)，2006，32（3）：31-35.WU Hong,LI Guoquan,QI Lehua.Lubrication temperature calculation of the oil outlet for accessory gearbox based on ANSYS[J].Aeroengine，2006，32（3）：31-35.(in Chinese)

[6]鄭志國，王宇峰.隨機(jī)振動(dòng)中的參數(shù)介紹及計(jì)算方法[J].電子產(chǎn)品可靠性與環(huán)境試驗(yàn)，2009，26（7）：45-48.ZHENG Zhiguo,WANG Yufeng.Introduction of parameters in random vibration and their calculation[J].Electronic Product Reliability and Environmental Testing，2009，26（7）：45-48.(in Chinese)

[7]王明珠.結(jié)構(gòu)疲勞振動(dòng)壽命分析方法研究[D].南京:南京航空航天大學(xué),2009.WANG Mingzhu.Research on life analysis method for structure vibration fatigue[D].Nanjing：Nanjing University of Aeronautics and Astronautics,2009.

[8]Petrucci G，Zuccarello B.On the estimation of the fatigue cycle distribution from spectral density data [J].Proceedings of the Institution of Mechanical Engineers,Part C:Journal of Mechanical Engineering Science，1999，213（8）：819-831.

[9]Petrucci G,Di Paola M,Zuccarello B.On the characterization of dynamic properties of random processes by spectral parameters[J].Journal of Applied Mechanics，2000，67（3）：519-526.

[10]Rychlik I.Note on cycle counts in irregular loads[J].Fatigue and Fracture and Engineering Materials and Structures，1993，16（4）：377-390.

[11]Frendahl M,Rychlik I.Rainflow analysis:Markov method[J].International Journal of Fatigue，1993，15（4）：265-272.

[12]Benasciutti D,Tovo R.Comparison of spectral methods for fatigue analysis of broad-band Gaussian random processes[J].Probabilistic Engineering Mechanics，2006，21（4）：287-299.

[13]Rychlik I.On the narrow-band approximation for expected fatigue damage[J].Probabilistic Engineering Mechanics，1993，8（1）：1-4.

[14]Lutes L D,Sarkani S.Stochastic analysis of structural and mechanical vibrations[M].Oxford：Butterworth-Heinemann,2003：261-272.

[15]姚起杭，姚軍.工程結(jié)構(gòu)的振動(dòng)疲勞問題[J].應(yīng)用力學(xué)學(xué)報(bào)，2006，23（1）：12-15.YAO Qihang,YAO Jun.Vibration fatigue in engineering structures[J].Chinese Journal of Applied Mechanics，2006，23（1）：12-15.(in Chinese)