速度对聚四氟乙烯摩擦系数影响的分子动力学模拟*

潘登 刘长鑫 张泽洋 高玉金 郝秀红†

1) (燕山大学机械工程学院,秦皇岛 066004)

2) (燕山大学自润滑关节轴承共性技术航空科技重点实验室,秦皇岛 066004)

3) (旭硝子汽车玻璃（中国）有限公司,秦皇岛 066004)

1 引 言

摩擦系数是衡量材料减摩性能的主要参数.速度作为影响摩擦系数的重要参数一直受到学者们的广泛关注.研究表明,接触界面的势垒可以影响摩擦力进而影响摩擦系数[1-3].摩擦速度增加引起的温升有助于界面上的原子翻越势垒从而降低摩擦系数[4].若在适当的摩擦速度下生成了起润作用的物质[5]也会降低摩擦系数.此外,摩擦速度增加引起的温升还可改变高分子材料的黏-弹性,从而改变高分子材料与其对磨副的相互作用力,最终影响摩擦系数[6-8].然而,聚四氟乙烯 (polytetrafluoroethylene,PTFE)作为一种典型的自润滑高分子材料,在低速摩擦时接触温度没有明显变化[9],但其摩擦系数随着速度的增加会先上升后下降[10].针对该现象,本文从微观角度分析接触区内PTFE分子形状的变化,从而揭示了速度对其摩擦系数的影响.

PTFE因具有摩擦系数低、耐高温等特点而得到广泛应用[11].学者们通过实验方法研究了速度对PTFE摩擦系数的影响.实验研究表明: 低速滑动时,PTFE可在对偶副的表面形成分子定向排列的转移膜,该转移膜可起固体润滑剂作用[12]; 而高速滑动时PTFE转移膜破裂,磨削呈块状和片状,因此表现出较高的摩擦系数[13].黄传辉[14]还发现速度较高时材料回弹不充分,从而导致摩擦系数升高.尽管实验研究在宏观上对PTFE的摩擦机理进行了一定的解释,但受实验仪器精度及功能的限制,仅用实验方法难以深入地揭示PTFE材料在实际摩擦过程中的微观作用机理.因此需要从分子水平上对PTFE的摩擦行为进行研究.

由于PTFE是一种半晶体高分子材料,既包含规则分布的PTFE分子,也包含不规则分布的PTFE分子.上述PTFE摩擦模型已对规则分布的PTFE的摩擦行为进行了深入研究,但不规则分布的PTFE的摩擦行为尚未提及.此外,在上述分子动力学模拟中,需在PTFE摩擦模型中将PTFE分子中的全部原子表示出来,从而增加了模拟所需时间,降低了计算效率.为提高分子动力学模拟的计算效率,可采用粗粒化方法分析PTFE的摩擦过程.粗粒化是指将全原子体系中的多个原子或原子团按照一定的映射关系转化成粗粒化粒子.由于粗粒化过程忽略了粒子内部各原子间的相互作用,从而提高了模拟效率,因此常用于分析材料的物理及化学性能[22-24].本文采用粗粒化分子动力学方法研究不规则分布的PTFEPTFE摩擦过程,分析不同速度下接触区内下层PTFE分子键长、键角分布、分子形状变化及接触过程中摩擦力和正压力的变化,得出了速度对PTFE摩擦系数的影响规律.由于分子动力学模拟时并未考虑空气中的分子与PTFE材料的相互作用,因此本文参考真空条件下PTFE的摩擦实验工况[10]设置模拟参数,如环境温度、外载荷及速度等.本文的研究结果可为改善PTFE及其复合材料的减摩性能提供理论参考.

2 计算方法

2.1 粒子间作用势

PTFE是由多个重复单元组成的链状高分子材料,其分子结构如图1(a)所示.当重复单元数n取10时[25],其分子模型如图1(b)所示.为增加计算效率,本文将一个重复单元简化成一个粒子,简化原理如图1(c)所示.由于分子链两端的单体较中间单体多一个氟原子,因此两端的单体被简化成2个端部粒子.中间的8个单体被简化成8个中间粒子.简化后的分子链如图1(d)所示.

图1 (a) 化学式; (b) 全原子模型; (c) 映射原理; (d) 粗粒化模型Fig.1.(a) Chemical structure; (b) all-atom model; (c) mapping scheme; (d) coarse-grained model.

粒子间相互作用势包括非成键相互作用势和成键相互作用势(图2).非成键相互作用势包括分子间粒子的相互作用及同一分子链上不相邻2个粒子间的相互作用.非成键相互作用势由Lennard-Jones (L-J)表示,其表达式为

其中,ε为势阱深度,σ为势能为零时2个粒子间的距离,r为粒子之间的距离,rc(～20 Å)为截断半径.

为了确定(1)式中U和r的关系,本文采用逆玻尔兹曼法(iterative Boltzmann inversion,IBI)[26]确定ε和σ的值.IBI法利用径向分布函数g(r)与目标径向分布函数gtarget(r)之差来优化并最终确定ε和σ的值.迭代修正过程可由下式表示:

图2 粒子间相互作用势Fig.2.Schematic of interaction among beads.

经过迭代后最终得到ε和σ的值分别为0.3 kcal/mol和5 Å.

粒子间成键相互作用势包括键伸缩势Ubond和键角弯曲势Uangle(图2).基于图1所示的全原子模型与粗粒化模型间的映射关系,粒子间的成键相互作用势可能包含多个峰值.为了提高计算精度,本文采用multicentered Gaussian-based势[27]来描述粒子间的成键相互作用势.键伸缩势及其相应的力表达式如下:

HER-2蛋白的表达与患者的性别、患病年龄和病变部位无关，与胃癌的分化程度和Lauren分型相关。见表1。

其中,l为成键的两个粒子间的距离;nb为高斯函数数量;lbci,Abi和wbi分别为第i个高斯函数的中心、面积及宽度.各参数如表1所列.键角弯曲势Uangle及其相应的力表达式如下:

表1 键伸缩势参数[28]Table 1.Parameters of the bond strength potential[28].

其中,θ为键角;na为高斯函数数量;θaci,Aai和wai分别为第i个高斯函数的中心、面积和宽度.各参数如表2所示.

表2 键角弯曲势参数[28]Table 2.Parameters of the bond angle potential[28].

2.2 PTFE摩擦模型

当PTFE与其对偶副摩擦时,PTFE会转移到其对偶副表面,使PTFE与其对偶副的摩擦变成PTFE-PTFE摩擦[29].因此,本文采用LAMMPS软件[30],基于2.1节中所示的粗粒化过程及参数,构建了用于分析PTFE-PTFE摩擦的分子动力学模型.模型构建过程如图3所示.图3(a)为上、下层PTFE分子的初始分布.上层PTFE包含120条分子链,下层PTFE包含1000条分子链.在摩擦磨损实验中,下试件在xoy平面内的尺寸通常比上试件大.为了使所构建的PTFE摩擦模型更接近真实的摩擦实验,在本文所构建的PTFE摩擦模型中,下层PTFE在xoy平面内的尺寸(135 Å×65 Å)较上层PTFE (55 Å×27.5 Å)大.上、下层PTFE的高度分别为142和135 Å.此外,在图3(a) 中,上、下层PTFE均为规则分布.为了避免规则分布的分子链影响最终的模拟结果,本文通过调整L-J势中截断半径rc的值,使上、下层PTFE分子链呈无规则分布.由于上层PTFE在x和y方向上的尺寸较下层PTFE小,为避免无规则化过程中上层PTFE的尺寸发生明显变化,本文将上层PTFE分子置于固定粒子中.固定粒子与PTFE中间粒子及端部粒子的相互作用势均由L-J势表示.此外,为避免下层PTFE分子在无规则化过程中沿z方向的尺寸发生明显变化,采用“fix wall/reflect”命令限制下层PTFE分子在z方向的运动范围(z＜140 Å).模拟时沿着x,y方向取周期边界,沿z方向取非周期边界.模拟过程基于微正则系综(NVE),时间步长为1 fs.经过60 ps后,上、下层PTFE呈无规则分布(图3(b)).

图3 PTFE-PTFE摩擦的粗粒化模型 (a)初始分布; (b)随机分布; (c)平衡态分布; (d)最终分布Fig.3.The coarse-grained PTFE-PTFE friction model: (a) Initial distribution; (b) random distribution; (c) equilibrium distribution;(d) final distribution.

为了获得PTFE的平衡态分布,将rc的值设为20 Å,基于(1)-(6)式求解粒子的分布.考虑到PTFE分子中各粒子间的摩擦作用,本文基于Langevin热浴,采用“fix langevin”命令将上、下层PTFE的温度控制在300 K.经过150 ps后,系统达到平衡状态(图3(c)).最终删去固定粒子得到双层PTFE的摩擦模型(图3(d)).上层PTFE的尺寸为34 Å×20 Å×142 Å,比固定粒子形成的框架尺寸小.这是因为PTFE粒子间的成键及非键相互作用减小了粒子间的初始分布距离(图3(a)).上、下层PTFE间的距离约为35 Å,下层PTFE的底部固定不动(z=0-15 Å).模拟PTFE-PTFE摩擦时,对上层PTFE分子沿-z方向施加外载荷,并使其以固定的速度沿x方向移动.通过调节上层PTFE的线速度从微观角度分析速度对PTFE摩擦学性能的影响.

3 结果与讨论

PTFE的摩擦学性能通常用摩擦系数和磨损量表示.鉴于分子动力学模型与实际摩擦试件在尺度上有较大差异,本文重点研究速度对PTFE摩擦系数的影响.

3.1 计算摩擦系数

在模拟PTFE-PTFE摩擦时,上层PTFE沿-z方向受外载荷后将与下层PTFE接触,并以一定的速度沿x方向移动.本文将接触过程中下层PTFE分子对上层PTFE分子沿x方向及z方向的作用力分别定义为摩擦力和正压力.在接触的任一时刻t,摩擦力和正压力的数学表达式分别如(7)式和(8)式所示.摩擦系数为摩擦力与正压力的比值,如(9)式所示.

其中,为摩擦力,N为下层PTFE所含粒子数,为t时刻粒子i与上层PTFE沿x方向的相互作用力.

其中µt为t时刻的摩擦系数.

由于分子动力学模拟并未考虑空气中的气体分子与PTFE的相互作用,本文所选取的模拟工况均基于真空条件下钢球与PTFE膜的摩擦实验工况[10].根据赫兹接触理论((10)式),当半径为5 mm的钢球与PTFE膜在真空条件下接触时,其平均接触应力P约为31 MPa.因此本文将上层PTFE所受外载荷设为31 MPa.

其中,F(6 N)为钢球所受外载荷,R0(5 mm)为钢球直径,µ1(0.4)与E1(1.14-1.42 GPa)分别为PTFE的泊松比和弹性模量,µ2(0.3)与E2(210 GPa)分别为钢球的泊松比和弹性模量.

当上层PTFE所受载荷为31 MPa,速度为0.3 m/s时,摩擦力、正压力和摩擦系数随模拟时间的变化如图4所示.随着模拟时间的增加,摩擦力先缓慢增加后趋于平稳(图4(a)).当上层PTFE与下层PTFE刚开始接触时,PTFE分子间的吸引力起主导作用.随着模拟时间的增加,上层PTFE浸入到下层PTFE中,上、下层PTFE的吸引力逐渐减小,排斥力逐渐增加.因此正压力随着模拟时间的增加先短暂下降后急剧上升(图4(b)).上、下层PTFE接触过程中,摩擦系数随着载荷的变化如图4(c)所示.当模拟时间小于90 ps时,摩擦系数略有波动; 当模拟时间超过90 ps时,摩擦系数趋于平稳,约为0.06.在不同的速度下,当模拟时间超过90 ps后,摩擦系数均趋于平稳.因此本文取模拟时间为90-150 ps之间的摩擦系数的均值作为相应速度下的摩擦系数,从而分析速度对PTFE摩擦系数的影响.

为避免上、下层PTFE的初始接触位置影响模拟结果,本文通过改变上层PTFE相对于下层PTFE沿x方向的位置,使上、下层PTFE在3个不同位置接触并产生移动.本文将3次模拟所得摩擦力及正压力的平均值作为该速度下的摩擦力与正压力.

3.2 速度对摩擦系数的影响

图5为摩擦系数随着速度的变化规律.可见,当外载荷保持不变时(31 MPa),随着速度的增加,摩擦系数先升高后降低.这与已发表文献[10]的实验研究结果一致.此外,图5中误差棒的变化范围为0.0024-0.0407.该范围介于平行滑动的最小误差(0.002)与垂直滑动的最大误差(0.089)之间[21].这是因为上、下层PTFE分子接触时,上层PTFE分子与下层PTFE分子之间的夹角介于0°到90°之间,即介于平行滑动与垂直滑动2种极限情况之间.

图4 (a) 摩擦力; (b)正压力; (c)摩擦系数随模拟时间的变化Fig.4.(a) Friction force; (b) normal force; (c) friction coefficient as a function of simulation time.

图5 当外载荷为31 MPa时,摩擦系数随速度的变化Fig.5.Effect of velocity on friction coefficient when the externally applied load is 31 MPa.

为了分析摩擦系数随速度先上升后下降的原因,本文分析了3种速度下(0.1,1.2和2.0 m/s),接触区内(25 Å＜x＜100 Å,z＞ 40 Å)下层PTFE分子形状变化及接触过程中摩擦力和正压力随模拟时间的变化.

3.3 速度对PTFE分子形状的影响

图6 不同速度下接触区内下层PTFE的键长 (a)和键角分布 (b)Fig.6.(a) Bond length and (b) bond angle distributions of the bottom PTFE molecules within the contact area under different velocities.

本文分析了3种速度下接触区内下层PTFE分子键长分布、键角分布及回转半径Rg(radius of gyration)的变化.图6为不同速度下接触区内下层PTFE分子键长和键角的数量分布.图6中数量越多意味着处于该键长和键角的粒子越多.由图6(a)可见,尽管接触区内大多数PTFE分子的键长在平衡距离附近(～2.7 Å),但随着速度的增加,处于平衡距离的粒子数先减少后增加,间距处在2.37 Å附近的粒子数量先上升后下降.接触区内的键角分布(图6(b))与键长分布相似,大多数PTFE粒子处于平衡键角附近(～175°),但随着速度的增加,处于平衡键角的粒子数先减少后增加,处于90°-140°范围内的粒子数先增加后减小.这是因为随着速度的增加,接触区内的PTFE分子被压缩,导致接触区内PTFE分子的键长及键角均变小.速度进一步增大时,键长和键角多处在2.7 Å和175°位置,产生这种现象的原因可能是上、下层PTFE分子的相互作用势导致接触区内下层PTFE分子沿上层PTFE分子移动方向倾斜[31],分子链内粒子间的相互作用势使键长和键角多处在平衡值附近.

回转半径可以反映聚合物分子形状的变化.本文分析了3种速度下接触区内PTFE分子回转半径沿x方向((11)式)及z方向((13)式)分量的变化.

其中,M和N1分别为一条PTFE分子链的质量和所含粒子数,mi为粒子i的质量,xi为粒子i的x方向坐标,xcm为该分子链质心的x方向坐标,其表达式如下:

其中,mj为一条分子链内第j个粒子的质量,xj为粒子j的x方向坐标.

同理可得Rgz的表达式如下:

其中,zi,zj分别为一条分子链内粒子i和j的z方向坐标.

图7为Rgx与Rgz的概率密度分布.概率密度值越大说明接触区内处于该状态的分子数越多.本文将上、下PTFE层未发生接触时回转半径的概率密度分布作为参照.对比Rgx与Rgz的参照曲线可见,Rgx的峰值介于6-8 Å之间,而Rgz的峰值介于0-3 Å之间.这说明上、下PTFE层未发生接触时,下层PTFE分子大多沿x方向分布,沿z方向略微弯曲.由图7(a)可见,当上层PTFE与下层PTFE接触且沿x方向移动时,Rgx的峰值明显下降.随着移动速度的增加,Rgx峰值下降幅度先变大后变小.这说明,随着移动速度的增加,接触区内PTFE分子沿着x方向的变形量先增大后减小.这与图6中接触区内PTFE分子键长和键角的分布结果相符.由图7(b)可见,当上层PTFE所受外载荷不变时,Rgz的峰值虽有下降但随速度的变化并不明显.这说明,速度对接触区内PTFE分子沿z方向变形量的影响很小.

图7 接触区内下层PTFE分子沿(a) x方向及(b) z方向的回转半径分布Fig.7.Distributions of radius of gyration along (a) x and(b) z directions of the bottom PTFE molecules within the contact area.

3.4 速度对摩擦力及正压力的影响

图8为3种速度下摩擦力与正压力随模拟时间的变化.由图8(a)可见,随着模拟时间的增加,摩擦力逐渐增加,且随着速度的增加,摩擦力先增大后减小.由图8(b)可见,尽管正压力随着模拟时间的增加而明显增加,但不同速度下正压力几乎保持不变.因此根据(9)式可得摩擦系数随着速度的增加而先增大后减小(图5).

通过对比图6,7和8可见,摩擦力随速度的变化趋势与接触区内下层PTFE分子键长和键角的分布,以及下层PTFE分子沿x方向的形状变化对应.正压力随速度的变化趋势与接触区内下层PTFE分子沿z方向的形状变化对应.这意味着速度的增加导致接触区内PTFE分子沿x方向变形,从而使粒子间的键长变短,键角变小.而PTFE分子链具有一定刚度,使键长及键角恢复到平衡位置,在此过程中产生一定的回弹力.接触区内下层PTFE分子链的变形量越大,回弹力就越大,导致其与上层PTFE分子之间的相互作用力增加,从而增大接触过程中的摩擦力(图8(a)).随着速度的进一步增大,接触区内PTFE分子的键长和键角分布趋于平衡值附近,且沿x方向的变形量变小.导致这种现象的原因很可能是速度进一步增大时接触区内下层PTFE分子链沿移动方向发生倾斜[31],改变了上、下层PTFE分子的夹角,使上、下层PTFE分子趋于平行滑动,从而降低了摩擦力,最终导致摩擦系数降低.

图8 (a) 摩擦力和(b) 正压力随模拟时间的变化Fig.8.(a) Friction force and (b) normal force as a function of simulation time.

4 结 论

本文采用粗粒化分子动力学方法构建了PTFE-PTFE摩擦模型,通过分析摩擦过程中接触区内下层PTFE分子键长分布、键角分布、分子形状变化及接触过程中摩擦力和正压力的变化,从微观角度探究了速度对PTFE摩擦系数的影响规律.随着速度的增加,接触区内下层PTFE分子的键长及键角变小,分子链弯曲量增大导致回弹力增大,从而增加摩擦力.当速度大于临界速度1.2 m/s时,接触区内下层PTFE分子链键长和键角多位于平衡值处,分子链弯曲程度变小,这很可能是由于分子链发生倾斜,使上、下层PTFE分子链的滑动趋于平行滑动,从而降低了摩擦力.然而,不同速度下接触区内PTFE分子沿z方向的变形量及正压力几乎不变,因此摩擦系数随着速度的增加而先增大后减小.本文加深了速度对PTFE摩擦系数影响的理解,为改善PTFE及其复合材料的减摩性能研究提供了理论参考.