基于Matlab的球轴承接触应力与变形和负荷分布的计算

Ξ

陈锦江,任成祖,徐燕申

(天津大学机械学院,天津　300072)

摘　要:提出了在MATLAB 环境下计算滚动轴承中接触应力与变形和负荷分布的方法,具有编程简洁高效,计算精度高,通用性好等

特点:给出了金属球轴承无量纲接触参数的计算曲线和两个计算实例。

关键词:球轴承;接触应力与变形;负荷分布;超越方程;MATLAB

中图分类号:TH133.33　　　文献标识码:A 　　　文章编号:1007-4414(2004)01-0059-02

　　众所周知,对于滚动轴承的分析和计算是相当复杂的,往往需要借助于相应的专用程序[1]。如接触应力与变形和负荷分布的计算是滚动轴承分析的基础,需要求解一个含有第一类和第二类完全椭圆积分等特殊积分的超越方程[2]。问题是此类方程无法采用解析方法,因而有的用C 或FORTRAN 等语言进行编程计算[3,4]、或采用筒化方程[5]、或寻找其它各种替代算法

[6～8]

、甚至直接查表进行插值等等,但都有不理想的地方。

由于MA TLAB 语言具有强大的计算和绘图等功能,笔者尝试基于MA TLAB 计算出了金属球轴承无量纲接触参数和径向负荷分布积分,并通过实例以图形曲线的方式给出了接触应力与变形以及负荷分布之间的函数关系。

1　球轴承的超越方程

1.1　接触应力和接触变形的超越方程

根据分析滚动轴承的赫兹空间弹性点接触问题时的赫兹理论

[1]

,为计算轴承球与套圈沟道相互接触时的接触应力与

变形,需要求解一个超越方程:

F (p )=

(1+k 2

)L (e )-2k 2

K (e )

(1-K 2)L (e )

(1)

也可化成:21-e 2

K (e )-[2-e 2

-F (p )e 2

]L (e )=0

(2)

式中:e 为待求的接触椭圆的偏心率,0≤e ≤1;k 为系数,k =(1

-e 2)1/2,0≤k ≤1;K (e )、L (e )为第一和第二类完全椭圆积分。

K (e )=

∫

π

2

(1-e 2sin 2φ

)-(1/2)

d φ(3)L (

e )=

∫

π

2

(1-e 2

sin 2

φ

)1/2

d φ

(4)

　　F (ρ

)为主曲率函数,当接触副形状确定后为已知量;1.2　负荷分布的超越方程

对于径向游隙不为零的深沟球轴承,受载最大的滚动体

负荷Q 0是:Q 0=F r /ZJ r (ε

)(5.a )此时的最大接触变形是:δmax =K n Q 0

2/3

(5.b )

负荷分布参数是:ε=[1-u r /(2δmax +u r )]/2(5.c )将上述3个式子联立后可得1个方程:

2K n (F r /ZJ r (ε

))2/3

-u r /(1-2ε

)+u r =0(5)

式中:K n 为系数,与轴承材料、结构参数和第一类完全椭圆积分等有关,为已知量;F r 、Z 、u r 为分别是轴承承受的径向负荷、

滚动体的数量、径向游隙,均已知;J r (ε

)为径向负荷分布积分。J r (ε

)=1

2

π∫

φ

L

-φ

L

1-

1-cos φ2ε

3/2

cos φd φ(6)

　　φL 为表征轴承负荷分布区的大小,

φL =cos -1

[u r /(2δmax +u r )],由式(5.c )可得

φL =cos

-1

(1-2ε)(7)

　　在式(1)或(2)中,未知的变量只有一个e ,隐含在这两类完全椭圆积分中;式(5)的未知变量也只有一个ε,隐含在积分(式(6))和上下限(式(7))中。因此,如何求解超越方程(2)和(5),则是轴承计算的关键,也是难点。

2　使用MAT LAB 求解超越方程

2.1　完全椭圆积分的计算

MA TLAB 关于K (e )和L (e )的定义表达式为:

K (m )=∫π2

(1-m sin 2φ

)-(1/2)

d φ;　　　　L (m )

=

∫

π

2

(1-m sin 2

φ

)-(1/2)

d φ　　　　

　　为获得式(3)和(4)的计算结果,只需编程满足条件m =

e 2并输入自变量e 的初始值后,直接调用一个计算完全椭圆

积分的内部函数指令

[K ,L ]=e llipk e (e )

(8)

　　这里的e 表示需要输入的自变量;K 表示输出第1类完全椭圆积分的计算值,L 是第2类的。作者以1°间隔输入e 时,计算值与数表[9]完全一致;图1表示由绘图指令plot 得到的结果。

2.2　径向负荷分布积分的计算

(1)建立一个文件名是,′funname ′,的函数式M 文件,以

便定义方程(6):

function J r =funname (ε)(9)

J r =(1-(1cos φ

)/2ε)3/2cos φ　　(2)另建一个命令文件,以调用计算数值积分的内部功能函数指令。这里仅给出采用自适应递推Lobatto 方法求数值积分的quadl ,调用格式为:

J repsilon =quadl (′funname ′,φL ,-φL )

(10)

　　J repsilon 表示式(6)的计算结果,得到的J r (ε

)与ε的关系见图2,与文献[2]的一样。

·

95·Ξ项　　目:教育部先进陶瓷与加工技术重点实验室资助项目

收稿日期:2003-05-06

作者简介:陈锦江(1962-),男,河北秦皇岛人,燕山大学在职博士研究生,主要研究机械设计。

2.3　求解方程的两个步骤

(1)建一个文件名是。′funn-nanme ′,的函数式M 文件,以

定义方程(2)(或(5))中出现的各个变量,同时调用函数指令

(8)(或(10))。

(2)求根。另建一个命令文件,以调用一个求单变量的

非线性方程零点的内部功能函数指令

z =fzero (′funn-name ′,[x 0,x 1]

(11)

Z 表示解出超越方程式(2)或(5)的根,也就是e 或ε值;[x 0,x 1]表示需要预先指定e 或ε可能存在的范围(初始区域)。

此外,由于式(2)和(5)无法求导,为防止求解过程中出现无解或丢根等现象,可通过绘图的函数指令等方法检验超越方程的数学特性,还能进一步确定根存在的初始区域

。

图1　K (e )、L (e )与e 的函数关系曲线　图2　J r (ε)与ε

的函数关系曲线　图3　m a 、m b 、

δ3及e 与F (ρ)的函数关系曲线3　球轴承无量纲接触参数与F (ρ)的函数关系

根据参考文献[2]中的式(3-20),可令δ3为;

δ3=[2K (e )/π][πK 2/2L (e )]1/3

(12)以及m a =[2L (e )/(

πK 2)1/3(13)m b =[2L (e )kJ

π]1/3(14)

　　式中的m a 、m b 和δ3是仅与接触椭圆偏心率e 有关的无量纲参数,通过式(1)和(8),可计算出这些参数与主曲率函数

F (ρ

)之间的关系,当结果以具体的数值形式输出时,与文献[2]的表3-1完全一致;

根据e 的定义还可得到它与F (ρ

)的关系e =[1-(m b /m a )2]1/2

(15)

　　这些计算结果见图3。

4　球轴承接触应力与变形的计算实例

4.1　所用的原始参数

取单列角接触球轴承中的几何结构参数为[2]:球中心圆直径D m =60mm ,球径D b =12.7mm ,内圈滚道沟曲率半径r i

=6.54mm ,内圈接触点滚道直经D i =48.6mm ,接触角α=

26°;并由此可以计算出主曲率函数F (ρ

)=0.954。轴承钢的弹性模量E =207(GPa ),泊松比v =0.30。

4.2　Q 与p 0、δ、a 和b 的因数关系

当F (ρ

)=0.954,球和套圈的材料均选择为轴承钢,使两接触物体相互压紧的法向负荷Q 的变化范围取为0～5000

(N ),依据上述所介绍的方法,可计算出此时Q 与最大接触应力p 0、球与套圈之间的相互趋近量δ、接触椭圆的长轴a 和短轴b 的函数关系。

5　球轴承径向负荷分布的计算实例

某型单列深沟球轴承的几何结构参数为:内外圈沟底直径分别是D i =52.291mm 、D 0=77.706mm ,内外圈滚道沟曲率半径r i =r o =6.6mm ,球径D b =12.7mm ,球数Z =9,径向载荷F r =00N 。则相应的计算结果为ε=0.4409。

若令式(5)的右端等于y 可构造一个新方程

y =2K n (F r /ZJ r (ε

))2/3-u r /(1-2ε)+u r (16)

　　通过作图的函数指令,能得到y 与ε之间的函数关系(图

4),从中进一步验证方程求根的正确性。

图4　参数y 与ε的函数关系曲线

6　结论

本文提出的基于MA TlAB 语言的计算方法,能方便地求解轴承中赫兹接触问题的超越方程,结果显示出该方法编程简单、容易实现、通用性好、精度高,可将轴承研究与设计人员从繁重的编程计算中出来,也为解决很多类似的工程实际问题提供了一种非常适用的方法和工具;而且若利用该软件强大的计算功能和众多的工具箱,则将具有十分广泛的应用前景,值得推广使用。

参考文献:

[1]　罗祝三,吴林丰.任意载荷作用下高速滚动轴承的寿命计算[J ].

机械工程学报,1990,26(5):57-65.

[2]　万长森.滚动轴承的分析方法[M ].北京:机械工业出版社,

19871

[3]　何高明.在Matlab 计算机语言程序中调用勒让德第一类完全椭

圆积分[J ].广东机械学院学报,1996,14(4):24-29.

[4]　陈家庆,信朝华.弹性点接触问题的数值求解与应用[J ].轴承,

2001(1):1-4.

[5]　杜迎辉,邱　明,蒋兴奇,等.高速角接触球轴承刚度计算[J ].轴

承,2001(11):5-8.

[6]　丁长安,张　雷,周福章,等.线接触弹性接触变形的解析算法

[J ].摩擦学学报,2001,21(2):135-138.

[7]　裴礼清,张　峰1滚动轴承的波动体与内圈接触部分的相对滑

动分析[J ]1机械设计与研究,2000(2):53-55.

[8]　周　海,陈家庆,陈　琪.牙轮钻头滚动轴承滚子接触问题的数

值模拟[J ].石油矿场机械,2000,29(6):5-9.

[9]　《数学手册》编写组.数学手册[M].北京:人民教育出版社,19791

·

06·