劉越,周廣明,張祖智,杜萬里,馬貴葉

(中國北方車輛研究所車輛傳動重點實驗室,北京 100072）

多軸承支撐軸系的精確變形計算研究

劉越,周廣明,張祖智,杜萬里,馬貴葉

(中國北方車輛研究所車輛傳動重點實驗室,北京 100072）

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

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

0 引言

獲取軸和軸承精確的變形量,是正確預測軸承和齒輪疲勞壽命的基礎,對輪齒齒面精細化設計具有重要意義[1-2]。

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

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

1 基本方程

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

圖1 軸系統(tǒng)結構圖Fig.1 Shafting model

1.1 階梯軸的有限元法計算

將階梯軸簡化為梁單元既能節(jié)省計算時間,又能滿足精度要求。由于傳動系統(tǒng)中多為短粗軸,為充分考慮剪切力對軸系變形的影響,選用Timoshenko梁作為分析單元。

考慮剪切變形影響后,Timoshenko梁的應變能方程為

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

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

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

剪切部分的應變能Us:

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

對單元剛度矩陣進行疊加,形成軸剛度矩陣,為(6×nShaft)×(6×nShaft)(nShaft為軸上節(jié)點個數)階的對稱性帶式稀疏矩陣。圖2為疊加過程示意圖,其中單元1和單元2分別為梁單元剛度矩陣,并且兩個單元通過節(jié)點2連接。

1.2 滾動軸承的剛度矩陣獲取

由于滾動軸承的受載滾動體個數與所受載荷的大小和方向有關,因此軸承剛度具有非線性特點。以球軸承為例,在無載荷狀態(tài)下,初始接觸角定義為

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

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

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

受載后的接觸角αj為

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

根據Hertz理論,球軸承受載的接觸區(qū)域為旋轉橢圓面,計算橢圓長半軸和短半軸之比χ需要求解一個含有第一類和第二類完全橢圓積分的超越方程,如(9)式、(10)式所示[2]:

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

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

根據文獻[2]中的理論,滾球軸承中滾動體受力Q與變形δ滿足:

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

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

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

1.3 軸和軸承的剛度耦合計算

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

圖5 模型簡化圖Fig.5 Simplified shafting model

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

2 算法

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

初值的選取對非線性方程組的求解至關重要,為有較好的計算結果,預設軸承剛度K0,通過X0= K-10F得到計算初值。

對(14)式和(15)式表示的5個非線性方程組進行泰勒展開:

令Jij=?Fi/?xj,得到偏導數在i點的雅克比矩陣Ji.

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

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

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

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

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

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

通過對軸承剛度的迭代,最終達到軸系力學平衡。迭代過程中通過軸和軸承的變形協(xié)調位移大小調節(jié)軸承的剛度值,不斷逼近目標結果。

3 算例

將圖1所示的結構作為算例,采用傳統(tǒng)方法和開發(fā)程序進行計算結果對比。軸的材料為合金鋼,材料參數如下:彈性模量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 計算框圖Fig.6 Computation block diagram

軸和軸承的結構尺寸如表1和表2所示。

表1 軸的尺寸參數Tab.1 Dimension parameters of shafts

采用傳統(tǒng)力學方法,將結構簡化為圖7所示的超靜定梁,其中,軸承視為剛性支撐點。運用矩陣位移法計算該結構的支反力,得到支座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.

根據本文提出的計算思路和算法,對算例所示結構作為軸系結構進行計算。將Y方向軸承主剛度KrY的迭代過程繪制為圖表,如圖8所示。此時,各軸承的徑向支撐反力為:軸承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 材料力學方法簡化模型Fig.7 Simplified model for mechanics of materials

圖8 軸承剛度迭代過程Fig.8 Bearing stiffness iteration

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

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

圖9 軸的線位移結果對比Fig.9 Linear displacement of shaft

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

4 結論

通過上述分析,可以得到如下結論:

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

圖10 軸的角位移結果對比Fig.10 Angular displacement of shaft

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

References)

[1] 羅繼偉,羅天宇.滾動軸承分析計算與應用[M].北京:機械工業(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.滾動軸承分析:軸承技術的基本概念[M].5版.羅繼偉,馬偉,楊咸啟,等,譯.北京:機械工業(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] 謝濤,劉品寬,陳在理.轉臺軸系軸承剛度矩陣的理論推導與數值計算[J].哈爾濱工業(yè)大學學報,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ā)動機角接觸球軸承剛度的一種實用分析方法[J].南京航空航天大學學報,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].吉林大學學報:工學版,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] 羅祝三,吳林豐.多支點軸系中滾動軸承的擬靜力學研究[J].南京航空航天大學學報,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].西安交通大學學報,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)群,羅繼偉,吳長春,等.滾動軸承剛度分析程序[J].計算力學學報,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

國防科技工業(yè)技術基礎科研項目(2009b2542134);車輛傳動重點實驗室基金項目(9140C340103120C34124)

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