航空弧齒錐齒輪輪齒應(yīng)力譜精確獲取及分布估計*

王延忠，吳朝陽，卞　勇，劉　旸

(1.北京航空航天大學(xué) 機械工程及自動化學(xué)院, 北京 100191；2.中國兵器工業(yè)新技術(shù)推廣研究所，北京 100089)

?

航空弧齒錐齒輪輪齒應(yīng)力譜精確獲取及分布估計*

王延忠1，吳朝陽1，卞勇2，劉旸1

(1.北京航空航天大學(xué) 機械工程及自動化學(xué)院, 北京 100191；2.中國兵器工業(yè)新技術(shù)推廣研究所，北京 100089)

航空弧齒錐齒輪在服役過程中，工況惡劣，振動沖擊多，工作過程中的輪齒隨機應(yīng)力譜分布規(guī)律對弧齒錐齒輪的設(shè)計和壽命預(yù)測具有重要作用。根據(jù)弧齒錐齒輪動力學(xué)模型，計算各工況下的動載系數(shù)，建立綜合考慮支持剛度變形和動載荷的弧齒錐齒輪加載接觸分析有限元模型，通過計算和處理得到弧齒錐齒輪工作過程中的齒面接觸應(yīng)力譜和齒根彎曲應(yīng)力譜樣本，進而利用混合正態(tài)分布參數(shù)估計法對輪齒隨機應(yīng)力譜分布進行估計。

弧齒錐齒輪;動載荷;支承剛度變形;加載接觸分析;輪齒隨機應(yīng)力譜分布

航空弧齒錐齒輪是直升機傳動系統(tǒng)中最重要、最復(fù)雜和最薄弱的動力元件。其常處在高速、重載條件下工作，工況惡劣，振動沖擊多，由齒輪失效導(dǎo)致的后果往往是災(zāi)難性的。精確獲取齒面應(yīng)力譜的分布規(guī)律對航空弧齒錐齒輪傳動系統(tǒng)的設(shè)計和壽命預(yù)測具有重大意義。

為了研究弧齒錐齒輪的輪齒隨機應(yīng)力譜分布，就有必要對齒輪應(yīng)力進行分析。目前，齒輪應(yīng)力計算分析通常是利用LTCA有限元法[1-2]；但這是建立在靜力分析基礎(chǔ)上，并未考慮動力因素的影響。由于嚙合剛度的時變性和制造誤差而產(chǎn)生的振動和沖擊，使得齒輪傳遞的載荷遠大于靜態(tài)載荷，國內(nèi)外研究學(xué)者通常利用集中參數(shù)法建立多自由度的齒輪振動模型,進而通過數(shù)值方法進行求解[3-4]。

本文根據(jù)航空弧齒錐齒輪傳動系統(tǒng)的功率譜，通過對弧齒錐齒輪動力學(xué)模型進行數(shù)值求解獲取動載荷系數(shù)；建立弧齒錐齒輪加載接觸分析有限元模型，將動載系數(shù)以及支承剛度和軸系變形量作為有限元模型的邊界條件輸入，精確獲取各個工況下的輪齒應(yīng)力，結(jié)合各工況時間比對輪齒應(yīng)力進行擴展，得到齒面整個工作過程的隨機應(yīng)力數(shù)據(jù)樣本；最后利用混合正態(tài)分布參數(shù)法進行估計，得到弧齒錐齒輪整個工作過程中的輪齒隨機應(yīng)力譜分布規(guī)律。

1　弧齒錐齒輪系統(tǒng)參數(shù)和功率譜

弧齒錐齒輪的參數(shù)見表1。弧齒錐齒輪工作功率譜見表2。其中，時間百分比η為該工況運行時間占所有工況運行時間和的比例。

表1　齒輪相關(guān)參數(shù)

表2　弧齒錐齒輪傳動系統(tǒng)輸入功率譜

2　各工況下的弧齒錐齒輪動載系數(shù)獲取

2.1動力學(xué)數(shù)學(xué)模型的建立

在齒輪動力學(xué)模型中(見圖1)，兩齒輪用集中質(zhì)量和集中轉(zhuǎn)動慣量來模擬，軸模擬成無質(zhì)量的剛體，用彈簧和阻尼器來模擬，輪齒嚙合剛度用一對彈簧和阻尼器模擬。

圖1　弧齒錐齒輪動力學(xué)模型

弧齒錐齒輪傳動系統(tǒng)可以等效為八自由度的動力學(xué)模型，系統(tǒng)的自由度表示為列向量：

q={X1,Y1,Z1,θ1x,X2,Y2,Z2,θ2y}T

式中，Xi、Yi和Zi(i=1,2) 分別為齒輪沿x軸、y軸移動和z軸橫向振動位移；θ1x、θ2y分別為主動齒輪繞x軸、被動齒輪繞y軸的扭轉(zhuǎn)振動位移。

齒面嚙合點間由于振動和誤差產(chǎn)生的沿嚙合點法線方向的相對位移為：

λn=(x1-x2)a1-(y1-y2)a2-

式中，a1=cosδ1sinαn；a2=cosδ1cosαnsinβ1，a3=cosαncosβ1；en(t)為齒輪副的法向靜態(tài)誤差；δ1為主動輪的節(jié)錐角；αn為法面壓力角；β1為主動輪的中點螺旋角；r1、r2為兩齒輪嚙合點半徑。

齒形誤差通常用正弦簡諧函數(shù)表示：

en(t)=e0+ercos(Ωht+φ)

式中，e0、er為輪齒誤差的常值和幅值，通常取e0=0；Ωh為齒輪的嚙合頻率；φ為相位角，通常取φ=0。

齒輪副在嚙合時的法向動載荷及沿坐標方向的分力如下：

嚙合剛度取為：

式中，ch為嚙合阻尼；km為嚙合剛度的平均值；Akl為剛度的l階諧波幅值；Φkl為剛度的l階初相位?；↓X錐齒輪的嚙合剛度可以根據(jù)LTCA原理獲取并進行公式擬合[5]。根據(jù)經(jīng)驗公式，齒輪嚙合阻尼為：

式中，J1，J2分別為主、從動輪的轉(zhuǎn)動慣量；R1，R2分別為主、從動輪的基圓半徑；ξg為輪齒嚙合阻尼比，本文取ξg=0.1。

弧齒錐齒輪傳動系統(tǒng)的振動微分方程組為：

將λn作為新的自由度，消去θ1x和θ2y，系統(tǒng)的自由度數(shù)由8個降為7個，合并后的方程為：

做如下變換：

引入狀態(tài)變量：

將上式中的7個二階微分方程變?yōu)?4個一階微分方程，可以通過龍格庫塔數(shù)值積分法進行求解。

2.2動載系數(shù)的求解

圖2　不同工況下齒輪的動載系數(shù)

從圖2中工況1～工況3可以看出，保持齒輪的轉(zhuǎn)速不變，隨著轉(zhuǎn)矩的增加，弧齒錐齒輪的動載系數(shù)減小，這是因為隨著載荷的增大，弧齒錐齒輪的傳動比增大，傳動誤差變小，傳動更加平穩(wěn)；由工況4和工況5可以看出，弧齒錐齒輪的動載系數(shù)隨著載荷的增大而增大，這是因為當載荷超過一定范圍后，齒輪的動態(tài)沖擊激勵占主導(dǎo)作用，隨著載荷的增大，齒輪的動態(tài)沖擊激勵變大，因而，弧齒錐齒輪的動載系數(shù)增大。

3　考慮動載荷和支承剛度的輪齒應(yīng)力獲取

3.1有限元模型的建立

由于動力學(xué)因素的影響，弧齒錐齒輪的實際工作載荷會比額定載荷增大；同時弧齒錐齒輪傳動系統(tǒng)在載荷作用下，軸承支承和軸系會發(fā)生變形，致使齒輪的嚙合位置發(fā)生變化，這些因素都在不同程度上影響了弧齒錐齒輪的傳動性能。為了提高弧齒錐齒輪的輪齒應(yīng)力的計算結(jié)果準確性，應(yīng)考慮弧齒錐齒輪傳動系統(tǒng)的動力學(xué)因素和支承剛度變形的影響。

結(jié)合弧齒錐齒輪的加載接觸有限元分析原理，建立該對考慮動力學(xué)因素和支承剛度變形的弧齒錐齒輪加載接觸分析有限元模型，大小輪設(shè)置同樣的材料屬性：彈性模量為210 GPa，泊松比為0.3?？臻g問題中對計算對象形狀適應(yīng)較好的是四面體單元，四面體單元對邊界擬合的能力強；因此，本文選擇四面體單元C3D10進行網(wǎng)格劃分。劃分網(wǎng)格后的有限元模型如圖3所示。

圖3　四面體單元C3D10劃分網(wǎng)格后的模型

軸承剛度采用彈簧進行模擬(見圖4)，將軸承支承面與參考點RP-pp-1在x,y和z方向的平移自由度進行耦合，然后將參考點采用“Spring>Creat>Connect points to ground”命令與空間固定點連接起來，設(shè)置彈簧的x、y和z方向的剛度值。

圖4　軸承剛度模擬

本文中的弧齒錐齒輪在傳動過程中發(fā)生二、三齒嚙合交替變化，一個輪齒嚙合周期歷程中的動載系數(shù)變化如圖5所示。根據(jù)輪齒的嚙合位置，將ABAQUS軟件中的Amplitude幅值曲線設(shè)置成相應(yīng)的動載系數(shù)。分別設(shè)定5個工況的載荷大小，提交Job，進行計算分析。

圖5　一個嚙合周期歷程中的動載系數(shù)變化

3.2有限元計算結(jié)果

圖6～圖10分別為5個工況有限元計算結(jié)果的齒面中點接觸應(yīng)力和齒根中點彎曲應(yīng)力在嚙合過程中的曲線變化圖，其中，曲線1、曲線2分別為不考慮和考慮動載系數(shù)的應(yīng)力計算結(jié)果。

圖6　工況1的齒面接觸應(yīng)力和齒根彎曲應(yīng)力

圖7　工況2的齒面接觸應(yīng)力和齒根彎曲應(yīng)力

圖8　工況3的齒面接觸應(yīng)力和齒根彎曲應(yīng)力

圖9　工況4的齒面接觸應(yīng)力和齒根彎曲應(yīng)力

圖10　工況5的齒面接觸應(yīng)力和齒根彎曲應(yīng)力

由圖6～圖10可知，考慮動載系數(shù)的弧齒錐齒輪有限元模型會使輪齒載荷增大，齒面接觸應(yīng)力和齒根彎曲應(yīng)力都會相應(yīng)地提升，受載時間也會變長；由于考慮支承剛度和軸系變形，輪齒在受載情況下會明顯發(fā)生提前嚙合，輪齒嚙合位置發(fā)生改變，進而會影響輪齒的應(yīng)力值，因此，考慮動載系數(shù)和支承剛度變形的應(yīng)力計算結(jié)果會更精確。

將各工況下提取的齒面接觸應(yīng)力值和齒根彎曲應(yīng)力值擴展該工況的時間百分比η×104次，然后將各工況擴展后的齒面接觸應(yīng)力和齒根彎曲應(yīng)力值分別進行疊加，得到弧齒錐齒輪整個工作過程中的齒面接觸應(yīng)力和齒根彎曲應(yīng)力數(shù)據(jù)樣本。

4　輪齒應(yīng)力譜統(tǒng)計及分布估計

為了對弧齒錐齒輪傳動系統(tǒng)的設(shè)計與壽命預(yù)測，就需要得到該弧齒錐齒輪傳動系統(tǒng)的輪齒應(yīng)力譜的總體分布規(guī)律；因此，有必要對應(yīng)力數(shù)據(jù)進行處理和統(tǒng)計。傳統(tǒng)的統(tǒng)計方法一般假設(shè)應(yīng)力譜數(shù)據(jù)來自于某種分布的樣本，然后進行研究推斷；但是由于數(shù)據(jù)形態(tài)面臨越來越復(fù)雜的趨勢，很多情況下不同的局部表現(xiàn)出不同的特性，利用單一參數(shù)分布族函數(shù)將帶來較大的偏差，現(xiàn)在比較成熟的有混合分布模型分析法[6-7]?；旌戏植寄Ｐ褪怯蒏. Person提出的，隨著計算機事業(yè)的蓬勃發(fā)展，計算能力不斷提高，混合分布已經(jīng)成為數(shù)據(jù)處理和統(tǒng)計的一個重要工具，并得到了廣泛的應(yīng)用。

本文先對輪齒應(yīng)力數(shù)據(jù)進行載荷等級分組，然后統(tǒng)計各載荷等級的幀數(shù)以及概率密度，最后利用混合正態(tài)分布法對輪齒應(yīng)力譜進行分布估計?；旌险龖B(tài)分布模型概率密度函數(shù)為：

式中，w-li為各正態(tài)分布的權(quán)重系數(shù)；μi為各正態(tài)分布的期望；σi為各正態(tài)分布的標準差；l為需要估計的正態(tài)分布的個數(shù)。具體數(shù)據(jù)處理流程如圖11所示。

編制相應(yīng)的計算程序，對弧齒錐齒輪整個工作過程中的齒面接觸應(yīng)力譜和齒根彎曲應(yīng)力譜樣本進行處理，得到輪齒隨機應(yīng)力譜的總體分布結(jié)果(見圖12和圖13)。

圖11　流程圖

圖12　齒面接觸應(yīng)力譜分布

由圖12可知，齒面接觸應(yīng)力譜的分布為：w-l1=0.335，μ1=10，σ1=5.4；w-l2=0.079，μ2=129.3，σ2=4.5；w-l3=0.532，μ3=148.6，σ3=107；w-l4=0.054，μ4=298，σ4=3.5。

由圖13可知，齒根彎曲應(yīng)力譜的分布為：w-l1=0.247，μ1=50.6，σ1=20.5；w-l2=0.183，μ2=215.5，σ2=22.7；w-l3=0.165，μ3=500.1，σ3=23.5；w-l4=0.405，μ4=266.9，σ4=95.1。

5　結(jié)語

通過上述分析，得到如下結(jié)論。

1)通過建立弧齒錐齒輪動力學(xué)模型，計算得到了某弧齒錐齒輪系統(tǒng)的各工況下的動載系數(shù)。

2)建立了一種綜合考慮動力學(xué)因素和支承剛度變形影響的弧齒錐齒輪的加載接觸分析有限元模型，計算結(jié)果表明，綜合考慮動力學(xué)因素和支承剛度變形，輪齒在受載情況下會明顯發(fā)生提前嚙合，齒根彎曲應(yīng)力和齒面接觸應(yīng)力會增大，受載時間也會變長。

3)基于某弧齒錐齒輪傳動系統(tǒng)的功率譜，通過計算分析得到弧齒錐齒輪整個工作過程中的齒面接觸應(yīng)力和齒根彎曲應(yīng)力數(shù)據(jù)樣本，并基于混合正態(tài)分布參數(shù)估計法對輪齒應(yīng)力數(shù)據(jù)進行處理，得到了輪齒隨機應(yīng)力譜的分布，為弧齒錐齒輪的設(shè)計和壽命預(yù)測提供了重要依據(jù)。

[1] 王延忠,劉溢溥,張祖智,等. 重載弧齒錐齒輪承載特性分析[J]. 新技術(shù)新工藝，2014(1)：60-63.

[2] 唐進元,陳嫦,蒲太平，等.考慮誤差、軸變形和軸承剛度的弧齒錐齒輪LTCA[J].機械傳動，2010，34(5)：5-8.

[3] 王立華, 黃亞宇, 李潤方, 等. 弧齒錐齒輪傳動系統(tǒng)的非線性振動特性研究[J]. 中國機械工程, 2007，18(3)：260-264.[4] 王三民, 沈允文, 董海軍. 含間隙和時變嚙合剛度的弧齒錐齒輪傳動系統(tǒng)非線性振動特性研究[J]. 機械工程學(xué)報, 2003，39(2)：28-32.

[5] 唐進元, 蒲太平.基于有限元法的螺旋錐齒輪嚙合剛度計算[J].機械工程學(xué)報, 2011，47(11)：23-29.

[6] Pearson K. Contributions to the theory of the mathematieal evolution [J].Philosophieal Transaetions of the Royal Soeiety of London，1894，185：71-110.

[7] Nagode M F M.On a new method for prediction of the scatter of loading spectra[J].International Journal of Fatigue，1998，20(4)： 271-277.

*國家自然科學(xué)基金資助項目(51275020)

國家高科技支撐項目(2014BAF08B01)

航空A計劃項目(ATR-125-02-101)

責任編輯馬彤

Precise Acquisition and Distribution Estimation of Tooth Random Stress Spectrum of Aviation Spiral Bevel Gear

WANG Yanzhong1, WU Zhaoyang1, BIAN Yong2, LIU Yang1

(1.School of Mechanical Engineering and Automation Beihang University， Beijing 100191， China；2.Chinese Advanced Technology Generalization Institute of CNGC， Beijing 100089， China)

Aviation spiral bevel gear is working under severe conditions which are full of vibration and impact. Tooth random stress distribution is important for spiral bevel gear design and life prediction. The loaded contact analysis with finite element model of spiral bevel gear is established, considering deformation of support stiffness and dynamic load, which is calculated based on the dynamic model. Samples of tooth contact stress and tooth root bending stress spectrum are obtained. Then use the mixed normal distribution parameter method to estimate the distribution of tooth stress spectrum.

spiral bevel gear，dynamic load，deformation of support stiffness，loaded contact analysis，distribution of tooth random stress spectrum

TH 132.41

A

王延忠 (1963-),男，教授，主要從事先進精密傳動加工技術(shù)等方面的研究。

2016-04-02