李云涛, 全齐全, 唐德威, 侯绪研

(哈尔滨工业大学 机器人技术与系统国家重点实验室，黑龙江 哈尔滨 150001)



轴孔协调接触建模与试验研究

李云涛, 全齐全, 唐德威, 侯绪研

(哈尔滨工业大学 机器人技术与系统国家重点实验室，黑龙江 哈尔滨 150001)

传统的轴孔接触解析模型通常研究轴与带孔无限大平面的接触问题，未考虑孔壁厚度对轴孔接触特性的影响。本文考虑孔壁厚度对轴孔接触特性的影响，基于Winkler弹性基础模型和几何约束条件，沿接触点的法线方向建立了轴孔协调接触解析模型。分析轴孔协调接触模型发现，接触力和接触刚度受间隙和孔壁厚度影响较大。基于图像法开展的轴孔接触试验结果表明，轴孔的接触刚度随接触深度的增大呈现非线性增大的特征，而且与以往的轴孔静态接触模型相比，本文所建立的模型与试验数据的吻合程度较高，因而本文提出的接触模型能够准确描述轴孔的协调接触特性。

协调接触模型；回转铰链；几何约束；Winkler；弹性层；图像法

由于回转铰链内存在间隙，轴孔接触面间的接触碰撞过程对机构动态特性产生一定影响[1-3]。轴孔静态接触模型可为回转铰链接触碰撞特性和机构动态特性的研究提供理论基础。本文针对小间隙的轴孔协调接触问题开展模型建立与试验研究工作。

解析法是最早被用于建立接触模型的方法。Hertz基于弹性半空间理论建立了物体间的接触模型[4-5]。Goldsmith建立了软材料和高初始接触速度的接触体间的接触模型[6]。Hertz和Goldsmith建立的模型仅适合解决含有大间隙且承受小载荷的轴孔非协调接触问题[7-9]。Steuermann将物体接触表面的几何形状描述为N次多项式的形式，但建立的协调接触模型的计算精度取决于N的取值[10]。Liu沿用Hertz模型的接触应力分布形式，基于几何约束条件细化了接触点距离函数[11]。Persson和Ciavarella忽略孔壁厚度对轴孔接触特性的影响，建立了带孔无限大平面的轴孔接触力模型[12-14]。

数值法已广泛用于研究物体结构对接触特性的影响[15-19]，但在求解接触问题时数值法的计算精度依赖于接触参数的设置以及网格的划分精度。

半解析法具有解析法计算准确和数值法处理边界灵活的优点。由于采用解析法难以完全确定轴孔的接触深度和弹性变形关系式，通常需要结合数值法解决此问题[20-23]。Fang采用有限元法获取了球铰接触应力分布指数与接触半径比的函数式，推导出球铰接触力模型[20]。Liu将孔半径作为弹性层厚度，基于有限元法和Winkler理论建立了间隙球铰和回转铰链接触力模型[21-23]。采用半解析法建模时，解决了物体接触的边界约束问题并节省了计算时间，但中间变量的函数只是近似的解析形式，因而半解析法的计算精度比纯解析法低。

本文建立了基于几何约束条件和Winkler弹性基础理论的轴孔协调接触模型，并采用图像法开展了轴孔接触模型验证试验。试验结果表明所建立的模型能够准确描述轴孔的接触特性。

1 轴孔接触模型建立

1.1 间隙函数

在建立轴孔接触模型之前，提出以下几点假设：1)忽略轴孔间摩擦对接触模型的影响；2)在接触区域内，接触变形前轴孔轮廓线上对应弧长相等的点相互接触；3)在轴的中心处加载并且孔的外轮廓线固定，轴孔变形后的轮廓线连续；4)接触区域的接触半角不大于90°。

在接触区域内，接触变形前轴孔轮廓线上对应弧长相等的两点A、B(图1所示)接触。A、B两点的接触过程可视为两个阶段：1)在载荷作用下，点A沿-Y方向发生刚性位移移动到点A′；2)由于轴孔接触发生弹性变形，A′和B两点分别经过法向un1和un2、切向位移uθ1和uθ2后接触。

图1 轴孔二维接触模型Fig.1 Two-dimensional contact model of shaft and hole

在轴孔刚开始接触时，初始间隙C0(θ)(图1中P1P2距离)为

C0(θ)=R2-r1(θ)

(1)

式中：θ为OP1与x轴夹角，r1(θ)为OP1长度，R2为孔半径。

在△OO1P1中，运用余弦定理建立ΔR与r1(θ)关系式：

式中：R1为轴半径，ΔR为轴孔初始间隙，ΔR=R2-R1。r1(θ)进一步表述为

由于ΔR≪R1，r1(θ)函数中ΔR的二阶及二阶以上的高次项可忽略。式(1)中C0(θ)可描述为

C0(θ)=ΔR(1+sin θ)

在外载荷P作用下，轴沿-y方向移动位移δ。在△OO1′P1′中，运用余弦定理可求出轴孔接触后的间隙函数C1(θ)。由于δ≪R1，所以δ的二阶及二阶以上的高次项可以忽略，则轴孔接触后的间隙函数为

C1(θ)=ΔR(1+sin θ)+δsin θ

(2)

1.2 法向弹性变形函数

接触区域内的法向弹性变形un(φ)为

un(φ)=un1(φ)+un2(φ)， φ∈[-ε,ε]

式中：φ为轴孔接触轮廓线上任意点对应的中心角，un1(φ)和un2(φ)分别为轴和孔的法向弹性变形，ε为接触半角。

考虑轴孔的法向弹性变形，轴孔的间隙函数为

C(θ)=C1(θ)+un(φ)

接触区域内间隙C(θ)为零，可根据上式求得法向弹性变形函数：

un(φ)=-C1(θ)

将上式代入式(2)，并用φ替换θ，则

un(φ)=δcos φ-ΔR(1-cos φ)

(3)

1.3 接触力与接触刚度模型

Winkler弹性基础模型可用于求解在载荷作用下弹性物体发生的弹性变形量[4]，其表达形式如下

式中：p为接触力，E为材料弹性模量，h为弹性层厚度，uz为沿接触深度方向弹性变形。

轴孔均具有一定的柔性，可等效为一系列弹簧的形式。接触角度为φ处的法向接触应力p(φ)与法向弹性变形关系为

(4)

式中：t1和t2(φ)分别为轴和孔的弹性层厚度，E1和E2分别为轴和孔的弹性模量。

当轴孔的材料相同时，E1=E2=E。由于轴孔间隙远小于弹性层厚度，因而在处理弹性层问题时忽略间隙值。图2所示的弹性层厚度满足如下几何关系:

式中：t2min为孔的最小壁厚。

图2 轴孔的弹性层示意图Fig.2 Elastic layer sketch of shaft and hole

联合式(3)和式(4)，获取接触角度为φ处的法向接触应力：

[(δ+ΔR)cos3φ-ΔRcos2φ]

(5)

在接触区域内，y向接触应力沿着接触宽度方向(x方向)分布为

式中：x为接触点的横坐标值。

将上式代入式(5)，可得轴孔y向接触应力:

(6)

在接触区域上任取一弹性微元(图2所示)，则弹性微元承受的外载荷为

dp=p(φ)rcos φdφ

式中：p(φ)为接触角度为φ时轴孔的接触应力，r为接触点与固定坐标系原点距离。由于间隙值ΔR和渗透深度δ且均远小于轴半径R1，因而r≈R1。

对法向接触应力积分获取单位接触长度外载荷

将式(5)代入上式，获取单位接触长度外载荷:

轴孔在接触边界处的弹性变形un(φ)为零，因而由式(3)获得轴孔接触边界条件:

(7)

将式(7)代入上式，获得单位长度的外载荷P与接触深度δ关系式(即轴孔接触力模型)：

(8)

轴孔的非线性接触刚度为

(9)

2 接触模型参数分析

针对所建立的轴孔接触模型，分析不同接触深度下轴孔的接触状态、y向的接触应力分布情况以及轴孔间隙与弹性层厚度分别对接触力和接触刚度的影响。

2.1 接触状态参数

随着接触深度的增大，轴孔的接触宽度和接触半角会相应发生变化，因而接触宽度和接触半角在一定程度上能够反映出轴孔的接触状态。

接触宽度与接触半角存在以下关系：

(10)

式中：a为接触宽度。

为了研究轴孔实时的接触状态，根据式(7)和式(10)绘制出接触宽度和接触半角随着接触深度的变化曲线(如图3所示)。孔的半径、孔的最小壁厚和轴孔间隙分别取12.5 mm、10 mm和0.05 mm。

图3 接触宽度和接触半角曲线Fig.3 Curve of contact semi-angle and contact width

随着接触深度增大，接触半角和接触宽度呈现非线性增大趋势，且接触半角和接触宽度的增长速度逐渐减缓。接触半角随着接触深度增大，逐渐趋于稳定值(90°)，此结果与假设条件4一致。

2.2 y向的接触应力分布

轴孔的y向接触应力分布曲线可按式(6)绘制(图4所示)。在分析接触模型的参数时，轴孔材料为铝合金7075。孔的半径、孔的最小壁厚和轴孔间隙分别取12.5 mm、10 mm和0.05 mm。

图4 y向接触应力分布Fig.4 Contact pressure along y direction

新模型中的y向接触应力分布形式与Hertz不同，未呈现出椭圆形分布。随着接触深度的不断增加，轴孔的接触宽度和峰值接触应力不断增大。轴孔的接触应力呈对称分布且峰值接触压力出现在接触区域中心。

2.3 间隙的影响

为了更直观的研究间隙对接触模型影响程度，按式(8)和式(9)分别绘制不同间隙下接触力和接触刚度随接触深度变化曲线(图5和图6所示)。孔半径和孔的最小壁厚分别取12.5 mm和10 mm。

图5 不同间隙的接触力和接触深度关系图Fig.5 Relation between contact force and contact depth with different clearances

图6 不同间隙的接触刚度和接触深度关系图Fig.6 Relation between contact stiffness and contact depth with different clearances

在相同接触深度下，接触力随着间隙的增大呈现非线性减小趋势。轴孔间隙越小，接触刚度趋于稳定的速度越快，因为间隙越小轴孔的接触宽度趋于轴半径的速度越快。在接触刚度增大还未稳定之前，相同接触深度下不同间隙的轴孔接触刚度值相差较大，但是在稳定后接触刚度值相差较小。

2.4 弹性层的影响

取不同的弹性层厚度值时，根据式(8)和式(9)分别绘制接触力和接触刚度随接触深度变化曲线(如图7和图8所示)。孔半径和轴孔间隙分别取12.5 mm和0.05 mm。

当轴孔接触深度一定时，接触力随弹性层厚度的增加而急剧减小。当弹性层的厚度较小时，接触刚度随着弹性层厚度的增加而迅速上升，呈现出较强的非线性特征。随着最小弹性层逐渐增大，接触刚度的非线性增长趋势逐渐减弱。这是因为当轴孔的材料一定时，Winkler弹性基础模型中接触刚度与弹性层的厚度成反比。当弹性层厚度达到一定数值，轴孔接触刚度的稳定值之间的差异逐渐减小。

图7 不同弹性层厚度的接触力和接触深度关系图Fig.7 Relation between contact force and contact depth with different thicknesses of elastic layer

图8 接触刚度随接触深度和弹性层厚度变化图Fig.8 Relation between contact stiffness and contact depth with different thicknesses of elastic layer

3 轴孔接触试验

目前轴孔接触试验研究主要是采用光弹性技术等技术测试试验件的接触应力分布情况[25-27]。Liu采用静载压力机获取回转铰链的接触刚度曲线，但未考虑轴在载荷作用下发生了挠曲变形[11]。为了验证所建立的轴孔协调接触模型的准确性，本文采用图像法开展了轴孔压缩试验。

3.1 试验原理

使用INSTRON高低温静载试验机在试件的圆环上表面加载(如图9所示)。压力机的位移测量精度(0.01 mm)不能满足试验中测量接触深度的精度要求，采用图像法可以提高位移测量精度。

图9 试验原理图及测试设备Fig.9 Schematic diagram and test setup of experiment

在试验中将厚度为5 mm的1级标准量块作为测量标准块，其实际尺寸与图像中的像素点个数的比值作为测量比例尺。图像中A、B两点的像素变化量可直接获取，然后再乘以比例尺可得到轴和孔的位移变化量ΔSAB。同理可获取C、D两点对应的位移变化量ΔSCD，最后取ΔSAB和ΔSCD的平均值作为轴孔的接触深度。

3.2 试验条件

静载压力机的加载精度为5‰Fmax。试验机的加载速率控制形式分为两种类型：位移速率加载和载荷速率加载。设置位移阈值或载荷阈值可停止加载。试验中采用载荷速率的方式加载，设置载荷阈值停止加载。加载量程、加载速率和加载间隔分别为0～15 000 N、60 N/s和1 000 N。当施加15 000 N的载荷时，由式(6)可知轴孔在y向的接触应力小于铝合金材料的许用接触应力，因而轴孔的材料在试验过程中未发生塑性变形。

每次加载结束后触发快门线按钮，使用1 800万像素相机获取A、B、C、D四处标记点及量块的图像。试验件的尺寸如下：R1=12.46 mm，R2=12.5 mm，L1=18 mm，L2=15 mm，t2min=10 mm，t0=5 mm。

3.3 试验结果及分析

将原始图像边缘处理，得到如图10所示的结果。量块上下边缘在图像中占据的像素点数为1 882，因而单个像素代表的尺寸为0.002 66 mm。

开展5组轴孔压缩试验，获取总载荷与轴孔双边接触深度数据，求出5组试验数据的平均值，绘制如图11所示的曲线。

图10 轴孔接触试验结果Fig.10 Test results when shaft and hole contact

图11 轴孔双边接触试验曲线Fig.11 Experimental curves of bilateral contact between shaft and hole

试验结果反映出轴孔的接触刚度随着接触深度的增大逐渐增大，且呈现出非线性的特征。将Hertz[5]、Persson[12]、Liu CS[23]、Liu RQ[11]以及新建立的接触力模型绘制成单边轴孔的接触刚度曲线(如图12所示)。

Hert接触模型的表达式为

(11)

式中：E*为轴和孔的复合弹性模量。

Persson接触模型的表达式为

(12)

Liu RQ 接触模型表达式为

(13)

Liu CS接触模型表达式为

(14)

图12 接触刚度曲线Fig.12 Curves of contact stiffness

Hertz接触模型中接触刚度随着接触深度的增大呈非线性增大趋势。在Hertz模型中假设轴孔上相互接触的点是沿着接触深度方向的点[5]，存在以下几何关系：

δ-(uz1+uz2)=z1-z2

(15)

式中：δ为轴孔接触深度，uz1和uz2分别为轴孔沿着接触深度方向的弹性变形，z1和z2分别为轴孔接触点距水平切平面距离：

Hertz接触模型假设轴孔间的接触区域很小，上式进一步简化为

式中：z1s和z2s分别为简化后的轴和孔接触点距切平面距离。

简化前后轴孔接触点高度差为

式中：Δz(x)和Δzs(x)分别为简化前和简化后轴孔接触点高度差。

简化前后轴孔接触点高度差的误差为

Δz(x)-Δzs(x)=

(16)

因为接触点沿x轴方向满足x∈[0,R1]，所以Δz(x)-Δzs(x)≥0，即简化后的轴孔接触点的高度差不大于简化前的高度。

由式(15)和式(16)可知，接触深度经过简化相应的减小，因而当接触力取相同值时接触刚度增大，而且随着接触宽度(或接触角度)的增大简化前后接触刚度的差别会更加明显。边界接触点的横坐标等于接触宽度值。当接触力为1 000 N时，简化前后轴孔接触点高度差的误差为0.021 5 mm，因而Hertz模型不适合描述轴孔协调接触问题。

Persson、Liu RQ、Liu CS的接触力模型曲线与试验曲线的逐渐接近。与其他接触模型相比，所建立的轴孔协调接触模型曲线与试验数据最贴近，因而能够准确的描述轴孔协调接触特性。

4 结论

本文提出了一种含间隙的轴孔协调接触解析模型，可为回转铰链的接触碰撞特性的研究提供参考。获得以下结论：

1)轴孔接触试验结果表明，本文提出的模型相对于以往的模型能较准确地描述轴孔间的静态接触特性。

2)在相同接触深度下，接触力随间隙的增大呈现非线性减小趋势，而且间隙越小，接触刚度趋于稳定的速度越快。

3)当弹性层厚度较小时，接触刚度对弹性层厚度的变化呈现出较强的非线性特征，而且随弹性层厚度逐渐增大，接触刚度的非线性增长趋势逐渐减弱。

[1]XU Lixin, YANG Yuhu, LI Yonggang, et al. Modeling and analysis of planar multibody systems containing deep groove ball bearing with clearance[J]. Mechanism and machine theory, 2012, 56: 69-88.

[2]FLORES P, MACHADO M, SILVA M T, et al. On the continuous contact force models for soft materials in multibody dynamics[J]. Multibody system dynamics, 2011, 25(3): 357-375.

[3]BAI Zhengfeng, ZHAO Yang. Dynamic behaviour analysis of planar mechanical systems with clearance in revolute joints using a new hybrid contact force model[J]. International journal of mechanical sciences, 2012, 54(1): 190-205.

[4]JOHNSON K L. Contact mechanics[M]. London: Cambridge University Press, 1985: 84-90.

[5]JOHNSON K L. One hundred years of hertz contact[J]. Proceedings of the institution of mechanical engineers, 1982, 196(1): 363-378.

[6]GOLDSMITH W. Impact: the theory and physical behavior of colliding solids[M]. London: Dover Publications, 2001.

[7]MACHADO M, MOREIRA P, FLORES P, et al. Compliant contact force models in multibody dynamics: evolution of the Hertz contact theory[J]. Mechanism and machine theory, 2012, 53: 99-121.

[8]王庚祥, 刘宏昭. 多体系统动力学中关节效应模型的研究进展[J]. 力学学报, 2015, 47(1): 31-50. WANG Gengxiang, LIU Hongzhao. Research progress of joint effects model in multibody system dynamics[J]. Chinese journal of theoretical and applied mechanics, 2015, 47(1): 31-50.

[9]PEREIRA C M, RAMALHO A L, AMBROSIO J A. A critical overview of internal and external cylinder contact force models[J]. Nonlinear dynamics, 2011, 63(4): 681-697.

[10]NOBLE B, HUSSAIN M A. Exact solution of certain dual series for indentation and inclusion problems[J]. International journal of engineering science, 1969, 7(11): 1149-1161.

[11]LIU Rongqiang, ZHANG Jing, GUO Hongwei, et al. Nonlinear characteristics of revolute joints with clearance[J]. Journal of aerospace technology and management, 2013, 5(3): 341-348.

[12]PERSSON A. On the stress distribution of cylindrical elastic bodies in contact[D]. Gothenburg, Sweden: Chalmers University, 1964.

[13]CIAVARELLA M, DECUZZI P. The state of stress induced by the plane frictionless cylindrical contact. Ⅰ. The case of elastic similarity[J]. International journal of solids and structures, 2001, 38(26/27): 4507-4523.

[14]CIAVARELLA M, DECUZZI P. The state of stress induced by the plane frictionless cylindrical contact. Ⅱ. The general case (elastic dissimilarity)[J]. International journal of solids and structures, 2001, 38(26/27): 4525-4533.

[15]颜东煌, 刘雪锋, 田仲初, 等. 销轴连接结构的接触应力分析[J]. 工程力学, 2008, 25(1): 229-234, 240. YAN Donghuang, LIU Xuefeng, TIAN Zhongchu, et al. Contact stress analysis of pinned connection[J]. Engineering mechanics, 2008, 25(1): 229-234, 240.

[16]KNIGHT M G, DE LACERDA L A, WROBEL L C, et al. Parametric study of the contact stresses around spherical and cylindrical inclusions[C]//Proceedings of the 10th International Workshop on Computational Mechanics of Materials. Galway, Ireland: National University of Ireland, 2002: 115-121.

[17]SUN Zhan'gang, HAO Caizhe. Conformal contact problems of ball-socket and ball[C]//Proceedings of International Conference on Solid State Devices and Materials Science. Macao, China: Elsevier, 2012: 209-214.

[18]LIU Shuangbiao, CHEN W W. Two-dimensional numerical analyses of double conforming contacts with effect of curvature[J]. International journal of solids and structures, 2012, 49(11/12): 1365-1374.

[19]LIU Shuangbiao. Numerical simulation of double conformal contacts involving both interference and clearance[J]. Tribology transactions, 2013, 56(5): 867-878.

[20]FANG Xin, ZHANG Chunhua, CHEN Xun, et al. A new universal approximate model for conformal contact and non-conformal contact of spherical surfaces[J]. Acta mechanica, 2015, 226(6): 1657-1672.

[21]LIU Caishan, ZHANG Ke, YANG Lei. The compliance contact model of cylindrical joints with clearances[J]. Acta mechanica sinica, 2005, 21(5): 451-458.

[22]LIU Caishan, ZHANG Ke, YANG Lei. Normal force-displacement relationship of spherical joints with clearances[J]. Journal of computational and nonlinear dynamics, 2005, 1(2): 160-167.

[23]LIU Caishan, ZHANG Ke, YANG Rei. The FEM analysis and approximate model for cylindrical joints with clearances[J]. Mechanism and machine theory, 2007, 42(2): 183-197.

[24]BURGUETE R L, Patterson E A. A photoelastic study of contact between a cylinder and a half-space[J]. Experimental mechanics, 1997, 37(3): 314-323.

[25]COLOMBO C, GUAGLIANO M. Photoelastic analysis of cylindrical elements with internal cracks under Hertz contact loading[J]. Fatigue & fracture of engineering materials & structures, 2010, 33(12): 885-896.

[26]章伊华, 林丹益, YANG Guoyu. 基于接触理论的螺栓联接接触面力学特性研究[J]. 兵工学报, 2015, 36(5): 946-952. ZHANG Yihua, LIN Danyi, YANG Guoyu. Research on mechanical properties of the contact surfaces in bolt joints[J]. Acta armamentarii, 2015, 36(5): 946-952.

[27]MARSHALL M B, LEWIS R, DWYER-JOYCE R S. Characterisation of contact pressure distribution in bolted joints[J]. Strain, 2006, 42(1): 31-43.

Modeling and experimental research on a coordinated contact between a shaft and hole

LI Yuntao, QUAN Qiquan, TANG Dewei, HOU Xuyan

(State Key Laboratory of Robotics and System, Harbin Institute of Technology, Harbin 150001, China)

In this study, the traditional contact models mainly focus on the contact between a shaft and an infinite plane with a hole, usually ignoring the effect of the hole's wall thickness on the contact characteristics. A coordinated contact model is established on the basis of the geometric constraint condition and the Winkler elastic foundation model, given that the contact characteristics of the shaft and the hole are influenced by the hole's wall thickness. The model is built along a normal direction of contact points. A model analysis shows that the contact force and stiffness are greatly affected by clearance and wall thickness of the hole. The contact experiments based on the image method show that the contact stiffness of the shaft-hole nonlinearly increases following the contact depth increase. In addition, the contact model has good consistency with the test data compared with the other models. The model can more accurately describe the contact characteristics of the shaft and the hole than the previous models, thereby providing a theoretical basis for the analysis of the dynamic characteristics of mechanism with clearance.

coordinated contact model; rotation hinge; geometric constraint; Winkler; elastic layer; image method

2015-09-17.

日期：2016-09-28.

国家自然科学基金项目(61403106).

李云涛(1987-), 男, 博士研究生； 全齐全(1983-), 男，副教授，硕士生导师.

全齐全，E-mail：quanqiquan@hit.edu.cn.

10.11990/jheu.201509043

TH131

A

1006-7043(2016) 11-1546-07

李云涛, 全齐全, 唐德威, 等. 轴孔协调接触建模与试验研究[J]. 哈尔滨工程大学学报, 2016, 37(11): 1546-1552. LI Yuntao, QUAN Qiquan, TANG Dewei, et al. Modeling and experimental research on a coordinated contact between a shaft and hole[J]. Journal of Harbin Engineering University, 2016, 37(11): 1546-1552.

网络出版地址：http://www.cnki.net/kcms/detail/23.1390.u.20160928.1419.050.html