劉越,周廣明,張祖智,杜萬(wàn)里,馬貴葉

(中國(guó)北方車輛研究所車輛傳動(dòng)重點(diǎn)實(shí)驗(yàn)室,北京 100072）

多軸承支撐軸系的精確變形計(jì)算研究

劉越,周廣明,張祖智,杜萬(wàn)里,馬貴葉

(中國(guó)北方車輛研究所車輛傳動(dòng)重點(diǎn)實(shí)驗(yàn)室,北京 100072）

針對(duì)傳動(dòng)系統(tǒng)中典型的多軸承支撐軸系結(jié)構(gòu),提出了軸和軸承剛度耦合建模方法以及軸系結(jié)構(gòu)精確變形迭代求解方法。基于軸承幾何學(xué)、接觸力學(xué)理論,采用Newton-Raphson算法實(shí)現(xiàn)軸承接觸應(yīng)力和變形計(jì)算、軸承內(nèi)部載荷分配計(jì)算、軸承滾動(dòng)體力學(xué)平衡方程組的求解,建立軸承剛度矩陣。實(shí)現(xiàn)軸承與Timoshenko梁的剛度耦合,形成軸系剛度矩陣,并通過(guò)松弛迭代法求解稀疏剛度矩陣。改變軸承剛度矩陣完成迭代過(guò)程,獲得系統(tǒng)的精確變形。通過(guò)算例,驗(yàn)證了該計(jì)算方法收斂快、對(duì)初值要求低,計(jì)算精度較高。

固體力學(xué);矩陣位移法;Timoshenko梁;軸承剛度;SSOR迭代法;Newton-Raphson算法

0 引言

獲取軸和軸承精確的變形量,是正確預(yù)測(cè)軸承和齒輪疲勞壽命的基礎(chǔ),對(duì)輪齒齒面精細(xì)化設(shè)計(jì)具有重要意義[1-2]。

傳統(tǒng)計(jì)算變形量的方法將軸承作為剛性支撐點(diǎn),根據(jù)受力平衡條件、虛位移原理等建立方程組,求解軸承的支撐反力和變形。這種方法雖然簡(jiǎn)便,卻無(wú)法得到軸承的準(zhǔn)確受力情況,也就無(wú)法計(jì)算軸和軸承的精確變形。事實(shí)上,軸和軸承剛度大小和變形量是相互影響的,無(wú)法獨(dú)立計(jì)算得到軸和軸承的精確變形[3]。因此,只能將軸和軸承作為整體進(jìn)行系統(tǒng)分析。

目前,國(guó)內(nèi)有關(guān)軸承剛度的研究較多[4-9],文獻(xiàn)[4-5]采用擬靜力學(xué)方法建立了角接觸球軸承剛度矩陣,文獻(xiàn)[6]通過(guò)實(shí)驗(yàn)驗(yàn)證了擬靜力學(xué)方法求解軸承剛度矩陣的準(zhǔn)確性。但對(duì)傳動(dòng)系統(tǒng)中的多軸承支撐軸系變形的研究較少,仍處于探索階段。文獻(xiàn)[7]對(duì)多支點(diǎn)軸系中滾動(dòng)軸承進(jìn)行了擬靜力學(xué)分析,將常見(jiàn)的階梯軸等效轉(zhuǎn)化為理想的等截面軸,需求解多個(gè)復(fù)雜的非線性方程組,方程數(shù)目龐大,迭代復(fù)雜。文獻(xiàn)[8]需要借助有限元軟件完成軸和軸承的剛度迭代,建立軸系的有限元模型使得計(jì)算效率和通用化程度受到較大的限制。相比而言,國(guó)外這方面的研究比較成熟[9-13],已相繼開(kāi)發(fā)出了相關(guān)的商業(yè)軟件和專用軟件,如Masta、Romax、Kissoft,這些軟件在工程設(shè)計(jì)中得到了良好的運(yùn)用。

1 基本方程

圖1為傳動(dòng)系統(tǒng)中典型的多軸承支撐軸系結(jié)構(gòu)。有3個(gè)滾球軸承作為支撐,軸端受主動(dòng)力矩,齒輪處作用在軸上的載荷有徑向力、切向力和被動(dòng)力矩。軸系采用直角坐標(biāo)系,坐標(biāo)原點(diǎn)O為軸上端面圓心,所在端面為YOZ面,旋轉(zhuǎn)軸為X軸。

圖1 軸系統(tǒng)結(jié)構(gòu)圖Fig.1 Shafting model

1.1 階梯軸的有限元法計(jì)算

將階梯軸簡(jiǎn)化為梁?jiǎn)卧饶芄?jié)省計(jì)算時(shí)間,又能滿足精度要求。由于傳動(dòng)系統(tǒng)中多為短粗軸,為充分考慮剪切力對(duì)軸系變形的影響,選用Timoshenko梁作為分析單元。

考慮剪切變形影響后,Timoshenko梁的應(yīng)變能方程為

將單元節(jié)點(diǎn)變量變成兩部分,即

對(duì)于彎曲變形引起的位移采用Hermite形函數(shù)插值表示,剪切引起的位移采用線性插值法,即有

將(3)式代入(1)式中,彎曲部分的應(yīng)變能Ub:

剪切部分的應(yīng)變能Us:

將Hermite形函數(shù)代入(4)式和(5)式中,得到Timoshenko梁的彎曲剛度矩陣和剪切剛度矩陣.

對(duì)單元?jiǎng)偠染仃囘M(jìn)行疊加,形成軸剛度矩陣,為(6×nShaft)×(6×nShaft)(nShaft為軸上節(jié)點(diǎn)個(gè)數(shù))階的對(duì)稱性帶式稀疏矩陣。圖2為疊加過(guò)程示意圖,其中單元1和單元2分別為梁?jiǎn)卧獎(jiǎng)偠染仃?并且兩個(gè)單元通過(guò)節(jié)點(diǎn)2連接。

1.2 滾動(dòng)軸承的剛度矩陣獲取

由于滾動(dòng)軸承的受載滾動(dòng)體個(gè)數(shù)與所受載荷的大小和方向有關(guān),因此軸承剛度具有非線性特點(diǎn)。以球軸承為例,在無(wú)載荷狀態(tài)下,初始接觸角定義為

式中:A=B×D,B=fi+fo-1,fi、fo分別為內(nèi)外滾道的溝曲率半徑系數(shù),D為滾球直徑;Pd為徑向游隙。

圖2 剛度矩陣耦合Fig.2 Stiffness matrix coupling

如圖3所示,在軸向力Fa和徑向力Fr的作用下,軸承內(nèi)圈產(chǎn)生軸向位移δa和δr,δa和δr為沿接觸線法向位移δn的分量。δn由(7)式確定:

受載后的接觸角αj為

圖3 徑向和軸向力作用下的角接觸球軸承Fig.3 Angular contact ball bearing under the action of the radial and axial forces

根據(jù)Hertz理論,球軸承受載的接觸區(qū)域?yàn)樾D(zhuǎn)橢圓面,計(jì)算橢圓長(zhǎng)半軸和短半軸之比χ需要求解一個(gè)含有第一類和第二類完全橢圓積分的超越方程,如(9)式、(10)式所示[2]:

式中:兩類橢圓積分K(χ)、E(χ)和橢圓接觸區(qū)長(zhǎng)短半軸之比χ滿足方程式(11)式,求解方程可得到χ:

式中:F(ρ)為滾球和滾道的曲率差。

根據(jù)文獻(xiàn)[2]中的理論,滾球軸承中滾動(dòng)體受力Q與變形δ滿足:

根據(jù)滾動(dòng)體的受載—變形關(guān)系,求解軸承所受外力與滾球作用力的平衡方程組(14)式和(15)式:

圖4 軸承變形的幾何參數(shù)Fig.4 Geometric parameters of bearing deformation

在原軸承位移量δ的基礎(chǔ)上取較小位移變量Δ,計(jì)算軸承反力,得到軸承剛度矩陣。

1.3 軸和軸承的剛度耦合計(jì)算

在多軸承支撐的軸系分析中,軸承內(nèi)圈的位移和軸配合節(jié)點(diǎn)處的位移滿足變形協(xié)調(diào)關(guān)系。把軸承單元作為變剛度的彈簧單元,其中一個(gè)節(jié)點(diǎn)(軸承內(nèi)圈)與軸配合節(jié)點(diǎn)耦合;另外一個(gè)節(jié)點(diǎn)(軸承外圈)為輸入的邊界條件。如軸承外圈與箱體配合,若箱體變形可忽略不計(jì),外圈節(jié)點(diǎn)的位移列陣置零;反之,可將箱體的變形量作為外圈的位移。將圖1所示的軸系簡(jiǎn)化為計(jì)算模型,如圖5所示。

圖5 模型簡(jiǎn)化圖Fig.5 Simplified shafting model

按照?qǐng)D2所示的剛度耦合方法,將軸承與軸的剛度矩陣耦合,形成包含軸承——軸的系統(tǒng)剛度矩陣。

2 算法

求解軸承剛度時(shí),需要對(duì)(14)式和(15)式所形成的5個(gè)非線性方程組求解。以軸承內(nèi)圈中心五自由度作為未知量,采用Newton-Raphson算法,利用梯度估計(jì)函數(shù)與自變量的截距,進(jìn)行多維求根。

初值的選取對(duì)非線性方程組的求解至關(guān)重要,為有較好的計(jì)算結(jié)果,預(yù)設(shè)軸承剛度K0,通過(guò)X0= K-10F得到計(jì)算初值。

對(duì)(14)式和(15)式表示的5個(gè)非線性方程組進(jìn)行泰勒展開(kāi):

令Jij=?Fi/?xj,得到偏導(dǎo)數(shù)在i點(diǎn)的雅克比矩陣Ji.

用矩陣形式表示(16)式為

解線性方程組(17)式得到修正量ΔX,然后進(jìn)行第二次計(jì)算,使X2=X+ΔX,至此,迭代過(guò)程產(chǎn)生。通過(guò)設(shè)置允許誤差ε,判斷是否結(jié)束計(jì)算。

根據(jù)疊加生成系統(tǒng)總體剛度矩陣KS,建立系統(tǒng)力學(xué)平衡方程式(18)式:

耦合生成的剛度矩陣KS為稀疏矩陣,采用直接法(如高斯消元法)求解方程組會(huì)破壞稀疏矩陣特性,產(chǎn)生較大的舍入誤差,影響計(jì)算結(jié)果,迭代過(guò)程難以收斂。采用SSOR迭代法取得了良好的計(jì)算效果。設(shè)剛度矩陣KS=I-L-U,SSOR的迭代方法如下:

式中:I、U、L分別是剛度矩陣KS的對(duì)角陣、下三角和上三角矩陣;S為SSOR的迭代矩陣,

計(jì)算方法采用變軸承剛度的計(jì)算思路。算法如圖6所示。

通過(guò)對(duì)軸承剛度的迭代,最終達(dá)到軸系力學(xué)平衡。迭代過(guò)程中通過(guò)軸和軸承的變形協(xié)調(diào)位移大小調(diào)節(jié)軸承的剛度值,不斷逼近目標(biāo)結(jié)果。

3 算例

將圖1所示的結(jié)構(gòu)作為算例,采用傳統(tǒng)方法和開(kāi)發(fā)程序進(jìn)行計(jì)算結(jié)果對(duì)比。軸的材料為合金鋼,材料參數(shù)如下:彈性模量E=206 000 MPa,切變模量G= 79380 MPa,泊松比ν=0.3.載荷條件如下:軸端處輸入扭矩1000 000 N·mm,齒輪所受載荷Fy=-2 426.5 N, Fz=-6 666.7 N,Mx=-1 000 000 N·mm,軸承外圈固定。

圖6 計(jì)算框圖Fig.6 Computation block diagram

軸和軸承的結(jié)構(gòu)尺寸如表1和表2所示。

表1 軸的尺寸參數(shù)Tab.1 Dimension parameters of shafts

采用傳統(tǒng)力學(xué)方法,將結(jié)構(gòu)簡(jiǎn)化為圖7所示的超靜定梁,其中,軸承視為剛性支撐點(diǎn)。運(yùn)用矩陣位移法計(jì)算該結(jié)構(gòu)的支反力,得到支座1的反力Fy1= -250.4 N,Fz1=-688.0 N,支座2的反力Fy2= 1688.9 N,Fz2=4640.1 N,支座3的反力Fy3=988.0 N, Fz3=2 714.5 N.

根據(jù)本文提出的計(jì)算思路和算法,對(duì)算例所示結(jié)構(gòu)作為軸系結(jié)構(gòu)進(jìn)行計(jì)算。將Y方向軸承主剛度KrY的迭代過(guò)程繪制為圖表,如圖8所示。此時(shí),各軸承的徑向支撐反力為:軸承1的支撐反力Fy1= 605.5 N,Fz1=1 299.2 N;軸承2的支撐反力Fy2= 470.8 N,Fz2=1 808.7 N;軸承3的支撐反力Fy3= 1 351.7 N,Fz3=3 558.8 N.

圖7 材料力學(xué)方法簡(jiǎn)化模型Fig.7 Simplified model for mechanics of materials

圖8 軸承剛度迭代過(guò)程Fig.8 Bearing stiffness iteration

對(duì)比上述兩種方法得到的計(jì)算結(jié)果可知,軸承的彈性變形不容忽略。軸承的實(shí)際剛度不僅會(huì)影響外載荷在各軸承上的載荷分配,甚至可以改變軸承的受力方向。

在國(guó)外先進(jìn)的商業(yè)軟件中建立圖1所示的模型,軸的材料選取合金鋼,材料參數(shù)設(shè)置為算例中的參數(shù),軸承選取表2所示的軸承,進(jìn)行變形計(jì)算。與依照本文提出的計(jì)算方法開(kāi)發(fā)的程序進(jìn)行結(jié)果對(duì)比,圖9和圖10為軸系線位移和角位移結(jié)果對(duì)比情況。

圖9 軸的線位移結(jié)果對(duì)比Fig.9 Linear displacement of shaft

從圖9和圖10可以看出,本文計(jì)算方法與商業(yè)軟件的計(jì)算結(jié)果吻合度很高,線位移中3軸承節(jié)點(diǎn)處的位移誤差均在1 μm以內(nèi),其余節(jié)點(diǎn)處的位移也均在1.5 μm以內(nèi)。所有節(jié)點(diǎn)處的角位移誤差均在0.002 mrad以內(nèi)。因此,完全滿足工程應(yīng)用。

4 結(jié)論

通過(guò)上述分析,可以得到如下結(jié)論:

1)在多軸承支撐軸系的變形分析中,軸承的彈性變形對(duì)軸承的受力分配影響很大。采用軸和軸承的剛度耦合聯(lián)合計(jì)算,比單獨(dú)計(jì)算軸和軸承變形更符合實(shí)際情況。

圖10 軸的角位移結(jié)果對(duì)比Fig.10 Angular displacement of shaft

2)在迭代求解軸承變形的過(guò)程中,提出了以軸承內(nèi)圈中心位移為未知量的五維求根算法。相對(duì)于以滾子變形作為未知量的常規(guī)算法,避免了大量非線性方程組的求解,提高了計(jì)算穩(wěn)定性,具有收斂快、初值要求低、準(zhǔn)確度高的特點(diǎn)。

References)

[1] 羅繼偉,羅天宇.滾動(dòng)軸承分析計(jì)算與應(yīng)用[M].北京:機(jī)械工業(yè)出版社,2009.

LUO Ji-wei,LUO Tian-yu.Rolling bearing analysis and calculation[M].Beijing:China Machine Press,2009.(in Chinese)

[2] Harris T A,Kotzalas M N.滾動(dòng)軸承分析:軸承技術(shù)的基本概念[M].5版.羅繼偉,馬偉,楊咸啟,等,譯.北京:機(jī)械工業(yè)出版社,2010.

Harris T A,Kotzalas M N.Rolling bearing analysis:Essential Concepts of Bearing technology[M].5th ed.LUO Ji-wei,Ma Wei,YANG Xian-qi,et al,translated.Beijing:China Machine Press,2010.(in Chinese)

[3] Redmond I.Study of a misaligned flexibly coupled shaft system having nonlinear bearings and cyclic coupling stiffness—theoretical model and analysis[J].Journal of Sound and Vibration,2010, 329(6):700-720.

[4] 謝濤,劉品寬,陳在理.轉(zhuǎn)臺(tái)軸系軸承剛度矩陣的理論推導(dǎo)與數(shù)值計(jì)算[J].哈爾濱工業(yè)大學(xué)學(xué)報(bào),2003,35(3):329-333.

XIE Tao,LIU Pin-kuan,CHEN Zai-li.Theoretical analysis and numerical estimation of stiffness matrix of angular contact ball bearing in a spindle-bearing system[J].Journal of Harbin Institute of Technology,2003,35(3):329-333.(in Chinese)

[5] 黃浩,張鵬順,溫建民.航空發(fā)動(dòng)機(jī)角接觸球軸承剛度的一種實(shí)用分析方法[J].南京航空航天大學(xué)學(xué)報(bào),2000,32(4): 422-427.

HUANG Hao,ZHANG Peng-shun,WEN Jian-min.An analytical method of stiffness for angular contact ball bearing of aeroengine [J].Journal of Nanjing University of Aeronautics&Astronautics, 2000,32(4):422-427.(in Chinese) [6] 方兵,張雷,曲興田,等.角接觸球軸承剛度理論計(jì)算與實(shí)驗(yàn)[J].吉林大學(xué)學(xué)報(bào):工學(xué)版,2012,42(4):840-844.

FANG Bing,ZHANG Lei,QU Xing-tian,et al.Theoretical and experimental research of stiffness of angular contact ball bearing [J].Journal of Jilin University:Engineering Edition,2012,42 (4):840-844.(in Chinese)

[8] 羅祝三,吳林豐.多支點(diǎn)軸系中滾動(dòng)軸承的擬靜力學(xué)研究[J].南京航空航天大學(xué)學(xué)報(bào),1995,27(3):291-298.

LUO Zhu-san,WU Lin-feng.Quasi-static analysis of multi-bearing shaft system[J].Journal of Nanjing University of Aeronautics& Astronautics,1995,27(3):291-298.(in Chinese)

[9] 劉顯軍,洪軍,朱永生,等.多支撐軸系軸承受力與剛度的有限元迭代計(jì)算方法[J].西安交通大學(xué)學(xué)報(bào),2010,44(1):41 -45.

LIU Xian-jun,HONG Jun,ZHU Yong-sheng,et al.Iterative method to solve bearing's force and stiffness for a multi-support spindle system based on finite element analysis[J].Journal of Xi' an Jiaotong University,2010,44(1):41-45.(in Chinese)

[10] 劉衛(wèi)群,羅繼偉,吳長(zhǎng)春,等.滾動(dòng)軸承剛度分析程序[J].計(jì)算力學(xué)學(xué)報(bào),2001,18(3):375-378.

LIU Wei-qun,LUO Ji-wei,WU Chang-chun,et al.Program for the stiffness analysis of rolling bearing[J].Chinese Journal of Computational Mechanics,2001,18(3):375-378.(in Chinese)

[11] Palmgren A.Ball and roller bearing engineering[M].3th ed. Burbank,Philadephia:SKF Industries Inc,1959.

[12] Ding J,Krodkiewski J M.Inclusion of static indetermination in the mathematical model for non-linear dynamic analyses of multibearing rotor system[J].Journal of Sound and Vibration,1993, 164(2):267-280.

[13] Zverv I,Pyoun Y S,Lee K B,et al.An elastic deformation model of high speed spindles built into ball bearings[J].Journal of Materials Processing Technology,2005,170(3):570-578.

Accurate Calculation for Deformation of Multi-bearing Shafting System

LIU Yue,ZHOU Guang-ming,ZHANG Zu-zhi,DU Wan-li,MA Gui-ye
(Science and Technology on Vehicle Transmission Laboratory,China North Vehicle Research Institute,Beijing 100072,China)

A new modeling method of shaft-bearing stiffness coupling and an iterative method for calculating the deformation of multi-bearing shafting system are presented for the typical complex multi-bearing shafting system used in driving system.The bearing stiffness matrix is created by Newton-Raphson method through the following steps:calculation of contact stress and deformation,calculation of internal load distribution,and solution of rolling elements mechanics equilibrium equations,which is based on the theory of bearing geometry and contact mechanics.The shafting system stiffness matrix is constructed by coupling the Timoshenko beam and bearings.The deformation of system can be solved by iterative matrix displacement method of variable bearing stiffness,which is applied to SSOR iterative matrix.It is proven that the method has better convergence,lower demand to initial data,and higher accuracy,which is more convenient to use for engineering design.

solid mechanics;matrix displacement method;Timoshenko beam;bearing stiffness;SSOR iterative matrix;Newton-Raphson method

TP391

:A

1000-1093(2014)03-0305-07

10.3969/j.issn.1000-1093.2014.03.003

2013-09-08

國(guó)防科技工業(yè)技術(shù)基礎(chǔ)科研項(xiàng)目(2009b2542134);車輛傳動(dòng)重點(diǎn)實(shí)驗(yàn)室基金項(xiàng)目(9140C340103120C34124)

劉越(1986—),男,工程師。E-mail:liuyue1620@sina.com