潘五九
, 李小彭
, 王雪
, 李木岩
东北大学 机械工程与自动化学院, 辽宁 沈阳 110819
收稿日期：2016-05-09
基金项目：国家自然科学基金资助项目（51275079，51575091）；中央高校基本科研业务专项资金资助项目（N160306003）。
作者简介：潘五九(1986-), 男, 安徽当涂人, 东北大学博士研究生;
李小彭(1976-), 男, 江西宁都人, 东北大学教授, 博士生导师。
摘要：考虑三维结合部形貌的W-M函数，推导了结合部静摩擦因数非线性分形理论模型.数值模拟了考虑三维形貌的结合部静摩擦因数与法向载荷
P、分形维数
D、分形尺度系数
G的关系，以及在二维分形和三维分形模型中的差异.结果表明：结合部静摩擦因数与法向载荷成单调递增关系，与分形尺度系数成单调递减关系.当
D小于2.5时，静摩擦因数随分形维数的增大而增大；当
D大于2.5时，静摩擦因数随分形维数的增大而减小；三维分形下的静摩擦因数小于二维分形下的静摩擦因数.
关键词：三维形貌结合部静摩擦因数分形理论非线性
Nonlinear Fractal Model for Static Friction Coefficient Considering Three-Dimensional Topography of Joint Surfaces
PAN Wu-jiu
, LI Xiao-peng
, WANG Xue
, LI Mu-yan
School of Mechanical Engineering & Automation, Northeastern University, Shenyang 110819, China
Corresponding author: LI Xiao-peng, professor, E-mail:
xpli@me.neu.edu.cnAbstract: Considering three-dimensional topography W-M function of joint surfaces, a nonlinear fractal model for static friction coefficient considering three-dimensional topography of joint surfaces was established. Relations between the static friction coefficient considering three-dimensional topography of joint surfaces and normal load, the fractal dimension, fractal roughness were numerically simulated, as well as the difference between static friction coefficients in the two-dimensional fractal model and three dimensional fractal model. The results showed that the static friction coefficient increases monotonically with the increase of normal load, and decreases monotonically with the increase of topothesy. When
D is less than 2.5, the static friction coefficient increases with the increase of fractal dimension; when the
D is greater than 2.5, the static friction coefficient decreases with the increase of fractal dimension; the three-dimensional fractal static friction coefficient is smaller than that of the two-dimensional.
Key Words: three-dimensional topographyjoint surfacesstatic friction coefficientfractal theorynonlinear
经过装配的部件间的接触面称之为结合部.结合部大量存在于机械系统中, 且系统的60%刚度和90%的阻尼皆产生于结合部
[1-2].正常工作的机器, 其结合部由于外载传递将不可避免地产生摩擦现象, 反过来摩擦又对系统动态稳定性产生影响
[3-5].因为工业界约有30%~50%的能量会通过热传导、振动、噪声等形式消耗在结合部摩擦上
[6]，因此, 预测和研究结合部的静摩擦特性是必要的.摩擦因数的概念首次由达芬奇给出并引入理论研究范畴, 至此很多学者进行了结合部静摩擦因数特性的研究.Chang等
[7]考虑结合部上粘着力的影响, 并基于统计的Greenwood-Williamson(G-W)模型, 给出了结合部静摩擦因数预测模型.田红亮等
[6]改进了分形理论公式推导, 得出了最大静摩擦力及摩擦因数等与微凸体实际接触面积间的函数预测模型; 兰国生等
[8]根据传统Majumdar-Bhushan (M-B)模型, 给出了静摩擦因数与法向载荷之间的函数变化趋势.
经典的G-W模型受限于测量仪器分辨率和采样长度, 所以推导出的相关模型
[7]具有局限性.M-B分形模型弥补了这一局限, 因其对表面形貌进行分形表征具有尺度独立性, 将不再受仪器分辨率和采样长度限制.但是以上基于分形理论建立的静摩擦因数模型
[6-8]均是建立在二维分形曲线(1 <
D < 2) 基础上, 用二维分形模型计算预测三维形貌(2 <
D < 3) 的接触特性, 本身就存在理论上的不合理性
[9-10].因此本文旨在修正的M-B分形接触理论和赫兹接触理论基础上, 考虑三维分形接触特性, 建立了新的结合部静摩擦因数非线性数学模型.
1 结合部建模由Hertz理论可知, 结合部上两微凸体在外载下相互接触可等效为一刚性光滑平面和一球形等效微凸体相接触, 如
图 1所示.
图 1(Fig. 1)
图 1中三角形
obc根据勾股定理易得
 | (1) |
式(1) 变形可得
 | (2) |
由于微凸体变形量
ω远小于
R, 故设
, 可得
 | (3) |
图 1中没有变形的等效微凸体和刚性光滑平面接触的横截面积为
 | (4) |
由经典赫兹弹性接触理论知, 微凸体处于弹性变形阶段时所承受的法向载荷为
[6] | (5) |
当微凸体接触时, 实际变形的微凸体与刚性平面所形成的接触面的半径为
[11] | (6) |
 | (7) |
式中:
E为微凸体的等效弹性模量;
E1和
E2分别为两相互接触的微凸体的弹性模量;
ν1和
ν2则分别为它们的泊松比.
将式(5) 代入式(6) 可得
 | (8) |
将式(3) 代入式(8) 可得实际变形的微凸体与刚性平面所形成的接触面半径
r与等效刚性平面截得等效微凸体后所形成的横截面半径
r′之间的关系:
 | (9) |
由式(4) 和式(9) 可知,
图 1中等效微凸体和刚性光滑平面实际接触面积为
 | (10) |
1.1 各变形阶段微凸体的临界接触面积Yan等
[11]改进单变量的W-M函数, 得到了能模拟三维表面形貌的曲面函数.修改后的双变量W-M函数为
 | (11) |
式中:
L为形貌取样长度;
D为结合部形貌分形维数(2 <
D < 3);
G为分形尺度系数;
γ(
γ > 1, 一般取1.5) 为决定频率密度的参数;
M为组成粗糙表面形貌时峰脊叠加的数目;
n为频率指数,
nmax=Int[lg(
L/Ls)/lg
γ];
Ls为最低截断长度;
x, y为表面微凸体的直角坐标系;
φm, n为乱数相位.
由式(11) 知, 该余弦函数的波峰与波谷幅值差可用来定义微凸体受载后的变形量, 即
 | (12) |
由式(4) 和式(10) 将式(12) 变形为面积
a的函数:
 | (13) |
将式(13) 代入式(10) 中, 得等效微凸体曲率半径:
 | (14) |
弹性变形阶段微凸体的临界变形量为
[12] | (15) |
式中:
为材料特征参数;
σy为较软材料的屈服强度;
kμ为摩擦力的修正参数,
μ为动摩擦系数.当0≤
μ≤0.3时,
kμ=1－0.228
μ; 当0.3 <
μ < 1时,
kμ=0.932 e
－1.58(μ－0.3).
由式(13)~式(15) 可推导出微凸体的临界弹性变形面积:
 | (16) |
1.2 各变形阶段微凸体的法向接触载荷由式(5)、式(13) 与式(14) 可将微凸体弹性变形阶段法向载荷表示成接触面积
a的函数:
 | (17) |
微凸体塑性变形阶段法向载荷为
 | (18) |
式中:
λ=H/σy, 为定义的系数,
H为较软材料的硬度.
1.3 结合部的实际接触面积结合部上所有微凸体的面积分布函数
n(
a)
[11]与其中最大的微凸体接触面积
al间的关系为
 | (19) |
对面积分布函数
n(
a)在微凸体各变形区间进行连续积分, 可得整个结合面上的实际接触面积
Ar.它应包含结合部上微凸体经历从弹性变形到塑性变形阶段的接触面积总和:
 | (20) |
2 法向接触载荷三维分形模型对单个微凸体在其弹性和塑性变形阶段与面积分布函数
n(
a)进行连续积分, 由式(17)~式(19) 可得总法向载荷.
1) 当2 <
D < 3,
D≠2.5时,
 | (21) |
2) 当2 <
D < 3,
D=2.5时,
 | (22) |
对式(21) 无量纲化可得
 | (23) |
式中:
P*=
P/(
EAa);
a1*=
a1/
Aa;
ac*=
ac/
Aa;
;
Aa为结合面名义接触面积.
对式(22) 进行无量纲化可得
 | (24) |
3 切向接触载荷三维分形模型根据汉米尔顿假设
[13], 若微凸体的屈服发生在接触点边缘, 接触边缘的3个主应力分别为
 | (25) |
 | (26) |
 | (27) |
式中:
Q, p分别为单个微凸体受到的切向与法向载荷.
由Tresca屈服条件知:
 | (28) |
将式(25)、式(26) 代入式(28), 得
 | (29) |
由式(29) 和式(19) 可知所有处于弹性变形阶段微凸体的切向载荷.
1) 当2 <
D < 3,
D≠2.5时,
 | (30) |
2) 当2 <
D < 3,
D=2.5时,
 | (31) |
对于已经发生完全塑性变形的微凸体, 没有承受切向载荷的能力, 故在计算结合部切向载荷时不涉及微凸体完全塑性变形阶段.那么结合部总切向接触载荷为
 | (32) |
对式(30) 和式(31) 分别进行无量纲化可得
 | (33) |
式中,
T*=
T/EAa.
 | (34) |
4 静摩擦因数非线性分形模型随着切向力的增大, 结合部将处于临界滑动状态, 此时结合部上微凸体所受的总切向载荷
T就是最大静摩擦力
[6-8].
因此, 由式(21)、式(22) 及式(30)~式(32) 可得考虑真实三维形貌静摩擦因数的非线性分形模型:
 | (35) |
根据式(23)、式(24)、式(33)、式(34) 可得其无量纲化模型:
 | (36) |
5 仿真与讨论分别分析静摩擦因数在不同
D和不同
G下的变化规律; 对比分析了二维分形和三维分形下静摩擦因数模型间的差异, 以及理论计算值与试验值间的误差.
5.1
D对三维静摩擦因数的影响令
φ=1,
ν=0.3,
G*=10
－5,
λ=2.8, 由式(35) 可得
D对静摩擦因数的影响, 如
图 2所示.可知:静摩擦因数
f均随无量纲法向载荷
P*的增大, 单调增大; 当
D小于2.5时, 静摩擦因数随分形维数的增大而增大; 当
D大于2.5时, 静摩擦因数随分形维数的增大而减小.
图 2(Fig. 2)
5.2
G*对三维静摩擦因数的影响令
D分别为2.3, 2.5, 2.7时, 分形尺度系数
G*对
f的影响如
图 3所示.可知:
f随着
G*的减小而增大;
f的增长率随
G*的增大而减小.
图 3(Fig. 3)
5.3 本文模型的验证为了验证文中所提模型的有效性，将其与现有试验数据进行对比.田红亮等
[14]以YK31320型滚齿机床的静压导轨结合面为研究对象, 给出了结合面上静摩擦因数随不同法向载荷的变化值.结合面相关参数如
表 1所示.结合面的分形参数
D, G分别为2.503 3, 2.006 6×10
-9 m.当结合面法向载荷
P=100 kN时, 由
表 1中数据及式(15) 和式(20) 可得此时微凸体的最大接触面积
a1=2.100 676×10
－4 m
2.再由式(29)、式(31) 和式(34) 可计算出此时的
f=0.523 0.如
图 4所示, 当
P=100 kN时的试验数据为
f=0.466 93.从
图 4a可知, 静摩擦因数与试验数据的变化趋势一致, 均随总法向载荷呈增长趋势.在
图 4b和
图 4c中, 分别给出了理论计算值与试验值的绝对和相对误差.绝对误差和相对误差的变化范围分别为-0.139 6~0.056 07和-38.78%~12.01%.另外, 试验值的最大与最小差值:
fmax－
fmin=0.466 93－0.36=0.106 93, 整体波动范围较小, 因此分析两者之间的绝对误差较为合理.基于以上分析, 本文模型对于工程结合面上静摩擦因数预测有一定的合理性.
图 4(Fig. 4)
表 1(Table 1)
 表 1 结合面参数Table 1 Parameters of joint surfaces | 参数 | 表面1(35钢) | 表面2(45钢) | | E/GPa | 197 | 205 | | ν | 0.31 | 0.3 | | ρ/kg·m-3) | 7 826 | 7 833 | | σy/MPa | 346 | 353 | | H/MPa | 478 | 500 | | γ | 1.5 | 1.5 |
| 表 1 结合面参数 Table 1 Parameters of joint surfaces |
6 结论1) 静摩擦因数与分形维数
D、分形尺度系数
G有着较复杂的非线性关系.
2) 静摩擦因数随
D的增大先增大后减小并以
D=2.5为界; 静摩擦因数随
G的减小而增大.
3) 对比本文理论模型计算值和试验值, 绝对误差范围为-0.139 6~0.056 07.
参考文献 | [1] | Ren Y, Beards C F. Identification of effective linear joints using coupling and joint identification techniques[J].ASME Journal of Vibration and Acoustics, 1998, 120(2): 331–338.DOI:10.1115/1.2893835 |
| [2] | Ibrahim R A, Pettit C L. Uncertainties and dynamic problems of bolted joints and other fasteners[J].Journal of Sound and Vibration, 2005, 279(3/5): 857–936. |
| [3] | 李小彭, 梁友鉴, 孙德华, 等. 系统参数对自激振动系统动力学稳定性的影响[J].东北大学学报(自然科学版), 2015, 36(5): 690–694.
( Li Xiao-peng, Liang You-jian, Sun De-hua, et al. Impact of the system parameters on self-excited vibration system dynamic stability[J].Journal of Northeastern University(Natural Science), 2015, 36(5): 690–694.) |
| [4] | Baramsky N, Seibel A, Schlattmann J. Modeling of friction-induced vibrations during tightening of bolted joints[J].Proceedings in Applied Mathematics and Mechanics, 2016, 16(1): 259–260.DOI:10.1002/pamm.v16.1 |
| [5] | Iroz I, Hanss M, Eberhard P. Transient simulation of friction-induced vibrations using an elastic multibody approach[J].Multibody System Dynamics, 2017, 39(1/2): 37–49. |
| [6] | 田红亮, 赵春华, 方子帆, 等. 金属材料表面静摩擦学特性的预测研究——理论模型[J].振动与冲击, 2013, 32(12): 40–44.
( Tian Hong-liang, Zhao Chun-hua, Fang Zi-fan, et al. Predication investigation on static tribological performance of metallic material surfaces-theoretical model[J].Journal of Vibration and Shock, 2013, 32(12): 40–44.DOI:10.3969/j.issn.1000-3835.2013.12.008) |
| [7] | Chang W R, Etsion I, Bogy D B. Static friction coefficient model for metallic rough surfaces[J].ASME Journal of Tribology, 1988, 110(1): 57–63.DOI:10.1115/1.3261575 |
| [8] | 兰国生, 张学良, 丁红钦, 等. 基于分形理论的结合面静摩擦因数改进模型[J].农业机械学报, 2012, 43(1): 213–218.
( Lan Guo-sheng, Zhang Xue-liang, Ding Hong-qin, et al. Modified model of static friction coefficient of joint interfaces based on fractal theory[J].Transactions of the Chinese Society for Agricultural Machinery, 2012, 43(1): 213–218.) |
| [9] | Liou J L, Lin J F. A modified fractal microcontact model developed for asperity heights with variable morphology parameters[J].Wear, 2010, 268(1/2): 133–144. |
| [10] | Luen L J, Fin L J. A new microcontact model developed for variable fractal dimension, topothesy, density of asperity, and probability density function of asperity heights[J].Journal of Applied Mechanics, 2007, 74(4): 603–613.DOI:10.1115/1.2338059 |
| [11] | Yan W, Komvopoulos K. Contact analysis of elastic-plastic fractal surfaces[J].Journal of Applied Physics, 1998, 84(7): 3617–3624.DOI:10.1063/1.368536 |
| [12] | Li X P, Yue B, Zhao G, et al. Fractal prediction model for normal contact damping of joint surfaces considering friction factors and its simulation[J].Advances in Mechanical Engineering, 2014, 9(1): 1–5. |
| [13] | Hamilton G M. Explicit equations for the stresses beneath a sliding spherical contact[J].Journal of Mechanical Engineering Science, 1983, 197(1): 53–59.DOI:10.1243/PIME_PROC_1983_197_076_02 |
| [14] | 田红亮, 刘芙蓉, 赵春华, 等. 金属材料表面静摩擦学特性的预测研究——实验佐证[J].振动与冲击, 2014, 33(1): 209–220.
( Tian Hong-liang, Liu Fu-rong, Zhao Chun-hua, et al. Prediction of static friction performance of metallic material surfaces with experimental proof[J].Journal of Vibration and Shock, 2014, 33(1): 209–220.) |
