图1 物理模型
Fig.1 Physical model
孟凡净1,2，刘 焜3，吴华伟2
(1．河南工学院机电工程系，新乡 453003；2．湖北文理学院纯电动汽车动力系统设计与测试湖北省重点实验室；襄阳 441053；3．合肥工业大学摩擦学研究所，合肥 230009)
摘 要：颗粒流润滑是可用于苛刻环境下的新型润滑方法，颗粒流的剪切膨胀特性是全面建立颗粒流润滑理论亟需解决的关键问题。该文基于离散单元法建立了颗粒流润滑的分析模型，研究了润滑启动阶段剪切膨胀现象的变化规律，剪切膨胀的速度适应机制和压力适应机制，并在总结剪切速度和压力载荷对剪切膨胀影响的基础上分析得出了剪切膨胀现象的力学机制。研究结果表明：润滑启动阶段，颗粒润滑系统的剪胀率随时步的变化服从二次多项式规律；颗粒润滑系统的剪切膨胀幅度随剪切速度，压力载荷的增大均呈现逐步递减的变化趋势；颗粒润滑系统的剪切膨胀特性与颗粒润滑介质间细观力链的波动密切相关，当细观力链的波动幅度变大时，颗粒润滑系统剪切膨胀的变化幅度亦增大，反之则减小。颗粒润滑介质间的接触力分布按幂函数规律变化，且颗粒润滑介质间力链的强度随着剪胀率绝对值的增大而变强。
关键词：颗粒流润滑；剪切膨胀；离散单元法；细观力链；力学机制
颗粒流润滑是根据颗粒物质力学理论提出的一种新型的润滑方式[1―3]，其基本思路是将硬质的、不黏附的固体颗粒作为润滑介质导入到摩擦副间隙, 达到减小作相对运动两表面接触的目的，从而保护表面免于损伤[4]。颗粒物质具有强耗散特性和温度对体系的影响可以忽略的特点，因此该润滑方式环境适应性强，可以适应超高温、高真空等极端工况环境，克服了油与脂润滑温度适应性差的问题，从而使得颗粒流润滑特别适用于一些特殊工况环境，如高温成形[5]、航空航天[6]等。
当处于摩擦副间隙的颗粒润滑介质形成颗粒流润滑时，随着剪应变的不断变化，颗粒润滑介质的空间位置和体积发生变化，为了适应这一变化，必然伴随着颗粒润滑介质沿剪切法方向的上下移动，这种现象称为剪切膨胀[7―8]，它是颗粒润滑介质状态变化过程中所特有的一种现象[9]。对颗粒流润滑的剪切膨胀特性进行研究，是进一步阐明颗粒流润滑动力学机制和承载机理的重要组成部分，也是颗粒流润滑理论进一步发展必须面对的关键问题。
现今，颗粒流润滑的研究处于起步阶段，还没有掌握颗粒流润滑剪切膨胀的力学机理，因此很有必要对颗粒流润滑剪切膨胀现象开展系统研究，通过研究掌握颗粒流润滑剪切膨胀的速度适应机制、压力适应机制并总结得出产生该变化的力学机理，从而为颗粒流润滑系统的参数选择提供指导。
本文利用基于非连续介质力学的离散单元法理论，建立了颗粒流润滑的离散元分析模型，研究了剪切运动初始阶段颗粒流润滑剪切膨胀的变化规律，速度适应机制和压力适应机制，并根据剪切膨胀随剪切速度和压力载荷的变化规律总结出了颗粒流润滑剪切膨胀现象的力学机制。
为了分析颗粒流润滑状态下剪切膨胀的力学机理，建立了如图1所示的颗粒流润滑离散元模型，模型的长度L和宽度B分别为1.7 mm和0.85 mm，接触模型的建立采用Hertz-Mindlin方法。
图1 物理模型
Fig.1 Physical model
模型的下摩擦副为一条无限长的线性墙体，在模型建立的初始阶段，通过对模型反复的初始化测验, 选择合适的墙体参数，在模拟运算过程中，下摩擦副以速度 u自左向右作剪切运动。 模型的上摩擦副是用CLUMPS命令粘接成板块的50个球形颗粒组成，允许上摩擦副沿坐标轴y向的微量移动，对上摩擦副添加x向和转动约束，并添加沿-y方向的压力载荷 P。模型中的颗粒润滑介质是由 1200个直径大小在27 μm ～ 41 μm之间的球形颗粒集合组成的，颗粒润滑介质直径范围的分布状态为均匀分布，且颗粒润滑介质间以及颗粒润滑介质与表面之间没有任何粘接, 处于完全松散堆积状态。
为了使颗粒润滑介质从沿坐标轴x向右侧出来时，在坐标轴x向左侧以同样的速度和位置出现，特在模型L×B的区域内添加了周期性空间。如图1所示，模型中的上摩擦副材料为1Cr18Ni9Ti，下摩擦副材料为 45#钢，颗粒润滑介质的材料为 SiC，详细参数如表1所示。
表1 模拟参数及数值大小
Table 1 Simulation parameters and value
颗粒流润滑中的颗粒润滑介质是离散和不连续的，无法用传统的连续性介质力学理论来分析解决。但是，Cundall等[10]等首先提出的离散单元法应用于颗粒流润滑系统中却能够很好的解决该散体系统的分析问题，国内外研究者利用该方法对颗粒物质系统进行了广泛的研究和探索[11―15]。
离散单元法的求解是采用动态显式算法。动态显式算法在差分求解过程中无需转置刚度矩阵，不用直接求解切线刚度，不用平衡迭代，使动态显式算法计算简单，计算速度快；在时间步长取的足够小的情况下，无收敛问题的存在；且可以很容易的进行并行计算，编制程序也比较简单。动态显式算法在求解需要分成许多时间增量来达到高精度分析的复杂非线性问题时具有很大优势，如散体接触、碰撞、冲击等高度非线性问题。
图2为离散元数值计算的详细流程图。如图2所示，在模型初始化测试通过后，对模型添加周期性空间，同时在上摩擦副上面添加压力载荷 P，颗粒润滑介质在压力 P和重力的共同作用下下压至初始的平衡状态。之后，下摩擦副以剪切速度 u自左向右作剪切运动，在摩擦的作用下把运动和力依次传给颗粒润滑介质和上摩擦副。
图2 计算流程图
Fig.2 Numerical calculation flow chart
图3为用简化方法表示的颗粒物质剪切流动过程中剪切膨胀现象的演示图形，图中的1表示下摩擦副，2表示上摩擦副， A～D表示颗粒润滑介质，黑色箭头表示下摩擦副和颗粒润滑介质的位移方向。
如图3(a)所示，当下摩擦副1自左向右做剪切运动时，在摩擦力的作用下颗粒A和B亦自左向右运动，同时推动颗粒C与D向右上方移动和上摩擦副2向上移动，当上摩擦副2上升到如图3(b)所示的位置时，到达最高点并停止上移；之后，如图3(c)所示，由于下摩擦副1持续运动，颗粒C和D开始向右下方移动，同时上摩擦副2也开始向下移动，如此循环往复出现。以上颗粒沿剪切法方向的上下移动现象称之为剪切膨胀现象。
图3 剪切膨胀现象的过程演示
Fig.3 The shear dilatancy process
对于颗粒流润滑系统来讲，其剪切启动阶段的瞬态行为非常复杂，且至今还未掌握在该阶段时系统的动力学特性特别是剪切膨胀规律，从而阻碍了研究工作的进一步深入开展。基于此原因，特对颗粒润滑系统在剪切运动初始阶段的剪切膨胀规律进行了如下研究。
图4 剪切膨胀现象的变化规律
Fig.4 The variation law of shear dilatancy
图4为剪切运动初始阶段(根据轴承的启动阶段来确定)颗粒流润滑系统剪胀率λK随时步K的变化规律以及它们的拟合规律。模拟分析时下摩擦副的剪切速度 u始终保持为4.5 m/s，剪胀率 λK的计算如式(1)所示。
其中：y0表示K=0时步(静态)时上摩擦副的y向坐标；yK表示K(K＞0)时步时上摩擦副的y向坐标，其坐标值可以通过离散元程序中的数据采集命令直接提取出来。
从图4可以看出，颗粒润滑介质在剪切运动初始阶段的剪胀率λK始终为负值，这说明在剪切运动初始阶段随着颗粒润滑介质剪切膨胀过程中的位移调整，使得上摩擦副的y向位置与起始时步(K=0)相比均出现了不同程度的下降, 但是下降的幅度较小。产生此现象的根本原因是，颗粒流润滑系统在静态时，压力载荷还未施加到颗粒物质系统中，颗粒润滑介质之间接触程度相对较为松散，颗粒润滑介质整体的弹性变形量较小，因此在静态时上摩擦副的y向位置最高；颗粒润滑系统在剪切运动状态时，在压力载荷和剪切运动的共同作用下，颗粒润滑介质之间的接触程度变得较为紧密, 颗粒润滑介质整体的弹性变形量相比静态时也进一步增大。原因是由于：颗粒流润滑系统在静态时，压力载荷还未通过上摩擦副施加到颗粒流润滑系统中。颗粒流润滑在剪切运动阶段时，压力载荷已经完全通过上摩擦副施加到了颗粒流润滑系统中，且在压力载荷的作用下，上摩擦副的y向位置会下降。因此与静态时相比较，上、下摩擦副间的间隙减小，且作用在颗粒润滑介质上的接触载荷明显增大，所以颗粒润滑介质之间的接触程度变得较为紧密。又由于颗粒所受到的接触力明显增大，且颗粒间的挤压程度变得紧密，故颗粒润滑介质整体的弹性变形量相对静态时也进一步增大。因此在剪切运动状态时上摩擦副的y向位置必定低于静态时。图4中的实体圆形为剪胀率散点图，为了从理论上分析剪切运动初始阶段剪胀率 λK随时步 K的变化规律，利用剪胀率散点并通过数值拟合的方法得到了剪胀率 λK随时步K变化的拟合曲线，如图4中的实体曲线所示。图4中颗粒润滑介质在剪切运动初始阶段的剪胀率Kλ随时步K的拟合规律服从形如式(2)所示的二次多项式规律变化。
其中：拟合系数α、β和γ的数值分别为0.00535、5.8917×10-8和 8.4977×10-14。
从图4中的拟合曲线可以看出，在上摩擦副压力P和重力的共同作用下，剪切运动初始阶段颗粒润滑介质整体的y向位置在剪切膨胀调整过程中首先出现了较大幅度的下移，当下降到一定的程度时又会出现小幅度的反弹上升。
颗粒流润滑是一种新型的润滑方式，若投入工业应用，比如应用颗粒流润滑理论设计新型的耐高温轴承，则该轴承在启动阶段时，上摩擦副位移随时步的变化特性则是该新型颗粒润滑轴承设计中需要着重考虑的一个技术难题。如果通过研究掌握了启动阶段上摩擦副位移随时步的变化规律，则亦可以掌握在稳定工作阶段时上、下摩擦副之间的间隙情况，从而可以为颗粒流润滑轴承设计中摩擦副间隙的合理调整提供指导。基于上述原因，对剪胀率 λK随时步 K的变化规律进行研究，不仅具有重要的理论意义，而且具有实际的应用价值。
与油润滑和脂润滑系统相比较，剪切膨胀现象是颗粒流润滑系统的特有现象。如图4所示，虽然从剪胀率与时步的拟合规律来看，剪胀率是随着时步呈收缩态势。但是，从图4中的E-F-G的单个剪切膨胀循环过程不难看出，颗粒流润滑系统在 E-F阶段是膨胀上升的，在F-G阶段则是收缩下降的。因此，颗粒流润滑系统摩擦副间的容积始终是处于膨胀然后收缩再膨胀再收缩的变化状态。图中 E-F段的膨胀阶段是颗粒流润滑系统的特有现象，油或脂润滑系统在受到压力载荷收缩后则不会出现此膨胀现象。
Campbell[16]通过对颗粒物质系统剪切流动行为的研究认为，滑动板的剪切速度是促使颗粒物质的流动行为发生转变的一个非常重要的外部输入条件。对于相同的颗粒物质系统，当剪切速度较小时，颗粒物质的流动行为以准静态流或慢速流为主，颗粒物质之间的接触行为主要为粘滑接触；随着剪切速度的逐渐增大，颗粒物质的流动行为逐渐转变为快速流，颗粒物质之间的接触行为也转变为滑动接触。粘滑接触和滑动接触的区别是：在粘滑接触状态下，颗粒间的接触状态由界面间的粘滞和滑移两种接触状态组成，而在滑动接触状态下，颗粒间的接触状态则主要是界面间的滑移。同样对于本文所涉及的颗粒流润滑系统，颗粒润滑介质的流动状态也通过下摩擦副剪切速度的变化而发生转变，而颗粒润滑介质流动状态的变化又与颗粒介质通过剪切膨胀过程中的位置调整密切相关，但是由于至今对颗粒流润滑研究尚未涉及剪切膨胀现象速度适应机制方面的研究工作，导致不能把剪切速度变化-剪切膨胀-流动状态转变之间形成有效串联，因此有必要通过研究掌握下摩擦副剪切速度对颗粒润滑介质剪切膨胀现象的影响规律，以期全面了解剪切速度影响的颗粒润滑介质流动状态的转变机制。
表2 拟合系数数值
Table 2 Numerical fitting coefficient
图5为下摩擦副的剪切速度从1.5 m/s增大到15 m/s时颗粒润滑介质剪切膨胀过程中的剪胀率λK随时步K的拟合规律图。
图5 剪切速度对剪切膨胀影响的拟合规律
Fig.5 The fitting law of shear dilatancy influencedby shear velocity
如图5所示；不同的下摩擦副剪切速度状态下，剪胀率λK随时步N的拟合规律均服从形如式(2)所示的二次项式规律变化，拟合系数α、β和γ的数值如表2所示。图6为颗粒润滑介质在时步K内的平均剪胀率λ随下摩擦副剪切速度u的变化图形，平均剪胀率λ的计算如式(3)所示。
式中，N为统计分析时所取的总时步数，在这里N=5×105。
从图5和图6可以看出，当下摩擦副剪切速度较小时(u=1.5 m/s，4.5 m/s，7.5 m/s)，颗粒润滑介质剪胀率λK随时步K变化的拟合规律为抛物线型，且在下降阶段的变化幅度较大，平均剪胀率λ的绝对值数值较大，分别为0.0125、0.0129和0.0122，但彼此之间的差值较小；当下摩擦副剪切速度增大到10.5 m/s以上时，颗粒润滑介质剪胀率λK随时步K变化的拟合规律变为近直线型；且在下降阶段的变化幅度逐渐减小，平均剪胀率λ的绝对值数值与下摩擦副剪切速度为 7.5 m/s时相比则出现了突然的下降，降幅分别为0.00299、0.003535和0.003989，但是它们之间的差值亦较小。这说明当下摩擦副剪切速度较小时(u=1.5 m/s、4.5 m/s、7.5 m/s)，对于颗粒流润滑系统来讲其剪切膨胀的幅度较大，但是当下摩擦副剪切速度增加到 7.5 m/s以上时，系统剪切膨胀的幅度则得到了抑制，变得较小。
图6 平均剪胀率随剪切速度的变化规律
Fig.6 The variation law of average dilatancy rate with increased shear velocity
从对图5和图6的分析结果，再结合作者在文献[17](模型与本文相同)中关于颗粒润滑介质动力学状态随下摩擦副剪切速度变化的研究结论可以得出，当下摩擦副剪切速度小于 7.5 m/s时，颗粒润滑介质的动力学状态为慢速流或准静态流，颗粒润滑介质之间的接触行为以粘滑接触为主，此状态下颗粒润滑介质彼此之间的接触时间较长，接触过程中颗粒润滑介质自身的弹性变形也较为充分，这为颗粒润滑介质沿坐标轴-y方向的位置调整提供了空间，同时也直接影响到了颗粒流润滑系统的剪切膨胀程度，使系统剪切膨胀的幅度变得较大。当下摩擦副剪切速度增大到10.5 m/s以上时，颗粒润滑介质的动力学状态为快速流，颗粒润滑介质之间的接触行为也演变为滑动接触，此状态下颗粒润滑介质之间的碰撞加剧，颗粒润滑介质彼此之间接触时间变得非常短，颗粒润滑介质自身的弹性变形与粘滑接触时相比也变得非常有限，致使颗粒润滑介质没有沿坐标轴-y方向的位置调整空间，这也直接抑制了颗粒流润滑系统的剪切膨胀程度，使剪切膨胀的幅度变得较小。
颗粒润滑介质弹性变形的程度与施加在其上压力载荷的大小密切相关，在弹性变形范围内对于单个颗粒润滑介质单体来讲，弹性变形量基本上随着压力载荷的增加而变大。当施加在上摩擦副上的压力载荷增大时，作用在每个颗粒润滑介质上的压力载荷的平均值也随之增大，使得颗粒润滑介质整体的弹性变形量也随之增加，进而会影响颗粒流润滑系统的剪切膨胀现象。通过以下的研究工作，分析剪胀率 λK随作用在上摩擦副上压力载荷的变化规律，掌握剪切运动初始阶段颗粒流润滑系统剪切膨胀现象的压力适应机制，从而为颗粒流润滑系统设计中压力参数的合理选择提供指导。
图7为剪切运动初始阶段施加在上摩擦副上的压力载荷P从3 MPa均匀增加到8 MPa时剪胀率λK随时步K变化的拟合图。图8为平均剪胀率λ随压力载荷P的变化规律图。图9为颗粒流润滑系统在压力和重力共同作用下达到初始平衡时，摩擦副上表面的y向坐标随压力载荷P的变化图。
如图7所示，当上摩擦副压力载荷P从3 MPa均匀增加到8 MPa时，剪胀率λK随时步N的拟合规律亦服从形如式(2)所示的二次项式规律变化，拟合系数α、β和γ的数值如表3所示。
图7 压力载荷对剪切膨胀影响的拟合规律
Fig.7 The fitting law of shear dilatancy influenced by the pressure load
从图7和图8可以明显的看出，当作用在上摩擦副的压力载荷P的范围为3 MPa～6 MPa时，剪胀率 λK随时步 K的拟合规律为开口向上的抛物线，平均剪胀率λ的绝对值数值较大，分别为0.01446、0.0146、0.01259、0.012，但彼此之间的差值较小；当作用在上摩擦副的压力载荷 P的范围为7 MPa～10 MPa时，剪胀率λK随时步K的拟合规律突变为开口向下的抛物线，平均剪胀率λ的绝对值数值也变得较小，分别为0.0057、0.0054和0.0052。这说明当作用在上摩擦副的压力载荷P低于6 MPa时，对于颗粒流润滑系统来讲其剪切膨胀的幅度较大，但是当作用在上摩擦副的压力载荷P增加到7 MPa以上时，系统剪切膨胀的幅度则得到了抑制，变得较小。
表3 拟合系数数值
Table 3 Numerical fitting coefficient
图8 平均剪胀率随压力载荷的变化规律
Fig.8 The variation law of average dilatancy rate with increased pressure load
图9 摩擦副上表面的y向位置随压力载荷的变化图
Fig.9 The position of upper plate in y direction with increased pressure load at initial pressure balance
从以上的分析结果并结合图9可以得出该颗粒流润滑系统剪切膨胀现象的压力适应机制如下所述。当施加的压力载荷P ＜ 6 MPa时，颗粒流润滑系统经过初压达到平衡状态后上-下摩擦副表面之间的空间较大(如图9所示)，同时颗粒润滑介质整体的弹性变形量较小，这为剪切运动过程中颗粒润滑介质整体的弹性变形往更深层次发展提供了充足的空间和条件，同时也使得系统的剪切膨胀幅度变得较大；当施加的压力载荷P增加到7 MPa以上时，颗粒流润滑系统经过初压达到平衡状态后上-下摩擦副表面之间的空间较小，同时颗粒润滑介质整体的弹性变形量较大，使得剪切运动过程中颗粒润滑介质整体的弹性变形变得非常有限，也抑制了系统的剪切膨胀幅度，使其变得较小。
图10为下摩擦副剪切速度为4.5 m/s，作用在上摩擦副的压力载荷P=2.5 MPa时颗粒流润滑系统的宏观流变图形，该图形是在时步K=100000时获得的。图10中的黑色网状结构为颗粒润滑介质之间的接触力形成的细观力链，该细观力链贯穿于颗粒润滑介质之间形成网络，且黑色线条越粗则说明细观力链越强，反之则越弱。
图10 颗粒润滑系统的宏细观图形
Fig.10 The macroscopic and mesoscopic figure of granular matter system
颗粒润滑介质之间的弹性变形量是微观参量，其大小和变化幅度可以通过图10所示细观力链的变化反映出来，也可以通过计算构成细观力链的接触力变化进行定量性描述，除此以外，颗粒流润滑系统的剪切膨胀与颗粒润滑介质整体的弹性变形密切相关，因此完全可以通过定量性研究颗粒流润滑系统力链变化与剪切膨胀之间的关系，掌握颗粒流润滑系统剪切膨胀现象的力学机制。
为了定量性研究颗粒流润滑系统力链变化与剪切膨胀之间的关联机制，引入了参数(D(FK))1/2如式(4)所示[18]。
式中：FK为K时步时所有颗粒润滑介质之间接触力大小之和；＜FK>为时步K=1到K=N内所有颗粒润滑介质之间接触力大小的平均值；N(D(FK))1/2为FK偏离＜FK>的均方差，并称该均方差(D(FK))1/2为“波动力链”，(D(FK))1/2越大，则说明力链的波动幅度越大，反之力链的波动幅度越小。
图11为剪切运动初始阶段下摩擦副剪切速度从1.5 m/s变化到15 m/s时，“波动力链”(D(FK))1/2随下摩擦副剪切速度的变化图。
图11 “波动力链”随剪切速度的变化规律
Fig.11 Variation of ‘fluctuation force chains’ with increased shear velocity
图12为剪切运动初始阶段施加在上摩擦副上的压力载荷P从3 MPa增加到15 MPa时，“波动力链”(D(FK))1/2随压力载荷P的变化图。
从图11可以看出，剪切运动初始阶段“波动力链”(D(FK))1/2随下摩擦副剪切速度的增大而剪小，与平均剪胀率 |λ| 随下摩擦副剪切速度的变化趋势相同(如图6所示)。从图12可以看出，剪切运动初始阶段“波动力链”(D(FK))1/2随压力载荷P的增大而减小，与平均剪胀率 |λ| 随压力载荷的变化趋势亦相同(如图8所示)。从图11和图12可以看出，平均剪胀率 |λ| 和“波动力链”(D(FK))1/2随下摩擦副剪切速度和压力载荷P的变化趋势是一致的。
图12 “波动力链”随压力载荷的变化规律
Fig.12 Variation of ‘fluctuation force chains’ with increased pressure load
图13为下摩擦副剪切速度为4.5 m/s时，颗粒润滑介质间通过相互接触挤压形成的接触力的分布概率。图13中的横坐标f表示无量纲接触力，具体的定义方法如式(5)所示：
式中：FIK表示第I(I=1-1200)个颗粒在K(K=1,2,3,…,N)时步时的接触力大小。图13中的P( f )表示无量纲接触力f的分布概率，把fmin≤f≤fmax区间从大到小等分为 15个等区间，统计得到无量纲接触力的概率分布如式(6)所示
式中：NJ表示无量纲接触力f在第J个区间的接触力数目。如图13所示，无量纲接触力的概率分布规律服从形如式(7)所示的幂函数规律变化，其中拟合系数a、b的数值分别为0.431和-2.037。从无量纲接触力概率分布的拟合规律可以清楚的得出，接触力越大，则无量纲接触力的概率分布P(f)越小。
为了进一步研究颗粒流润滑系统的剪切膨胀现象与力链之间的关系，根据如图13所示的无量纲接触力f的大小和分布特性，把颗粒润滑介质间的力链按照无量纲接触力f的大小分为弱力链和强力链，0≤f≤fmax/2范围内的力链定义为弱力链，fmax/2≤f≤fmax范围内的力链定义为强力链。
图13 接触力的概率分布规律
Fig.13 The probability distribution of contact forces
图14(a)和图14(b)分别为下摩擦副剪切速度为4.5 m/s时，颗粒润滑介质间的弱力链和强力链随剪胀率λK的变化规律。从图14可以看出，弱力链和强力链的分布概率P( f )均随着无量纲接触力f的增大而减小。当剪胀率 |λK| 从0.006增大到0.014时，弱力链所占的份量则从93.56%下降为77.3%，强力链所占的份量则从6.44%上升为22.7%。所以，当剪胀率|λK|增大时，强力链在颗粒流润滑系统中的数量会明显增多，而弱力链的数量则会明显减少。在如图14所示的五种不同的剪胀率状态下，弱力链所占的份量约为总力链的 78%～94%，而强力链仅占总力链份量的6%～22%。结果表明，颗粒润滑介质间力链的强度随着剪胀率|λK|增大而变强。在颗粒流润滑系统中，弱力链分布最广，而强力链的分布则较为稀疏。颗粒润滑介质间形成的少数的强力链阻碍了颗粒润滑介质的宏观流动速度，而弱力链则对强力链起到了支撑作用。
从以上的分析结果可以看出，颗粒流润滑系统的剪切膨胀现象是颗粒润滑介质整体弹性变化的宏观外在反映，而颗粒润滑介质间的细观力链变化又反映了颗粒润滑介质整体弹性变形量的大小和幅度，据此可以总结出颗粒流润滑系统剪切膨胀现象的力学机制如下所述。 颗粒润滑介质在剪切运动过程中，当颗粒润滑介质间的接触力发生变化时，会使细观力链发生变化并通过如图10所示的图形显现出来，同时由于接触力的变化会影响颗粒润滑介质的弹性变形量，而该弹性变形量的变化系统又以宏观剪切膨胀的形式表现出来，从而会使摩擦副上表面上下移动以适应颗粒润滑介质的位置调整。
图14 强、弱力链概率分布规律
Fig.14 The probability distribution of force chains
本文建立了颗粒流润滑的离散元分析模型，研究了润滑启动阶段剪切膨胀现象的变化规律，速度适应机制和压力适应机制，分析得出了剪切膨胀现象的力学机制。研究结果表明：
(1)剪切运动初始阶段剪胀率λK随时步K的变化特性服从形如式(2)所示的二次多项式规律，颗粒润滑介质整体的y向位置在剪切膨胀调整过程中首先出现了较大幅度的下移，在后期又出现小幅度的反弹上升。
(2)颗粒润滑介质的动力学状态为慢速流或准静态流时，系统剪切膨胀的幅度较大。当颗粒润滑介质的动力学状态演变为快速流状态时，系统剪切膨胀的幅度则变得较小。
(3)当压力载荷P较小时，系统初压平衡后颗粒润滑介质的弹性变形量较小，上-下摩擦副之间的空间较大，这为颗粒润滑介质的弹性变形往更深层次发展提供了充足的空间和条件，使得系统的剪切膨胀幅度较大；当压力载荷P较大时，系统初压平衡后颗粒润滑介质的弹性变形量较大， 上-下摩擦副之间的空间较小，从而使得系统的剪切膨胀幅度变得较小。
(4)当颗粒润滑介质间的接触力发生变化时，会使细观力链发生变化，同时由于接触力的变化会影响颗粒润滑介质的弹性变形量，而该弹性变形量的变化系统又以宏观剪切膨胀的形式表现出来，从而会使摩擦副上表面上下移动以适应颗粒润滑介质的位置调整。
参考文献：
[1]Wang W, Gu W, Liu K. DEM simulation on the startup dynamic process of a plain journal bearing lubricated by granular media [J]. Tribology Transactions, 2014,57(2): 98―205.
[2]王光谦, 孙其诚. 颗粒物质及其多尺度结构统计规律[J]. 工程力学, 2009, 26(增刊II): 1―7.Wang Guangqian, Sun Qicheng. Granular matter and the scaling laws [J]. Engineering Mechanics, 2009, 26(Sup.II): 1―7. (in Chinese)
[3]杨智勇, 曹子君, 李典庆, 等. 颗粒接触摩擦系数空间变异性对颗粒流双轴数值试验的影响[J]. 工程力学,2017, 34(5): 235―246.Yang Zhiyong, Cao Zijun, Li Dianqing, et al. Effect of spatially variable friction coefficient of granular materials on its macro-mechanical behaviors using biaxial compression numerical simulation [J].Engineering Mechanics, 2017, 34(5): 235―246. (in Chinese)
[4]王伟, 刘小君, 焦明华, 刘焜. 颗粒运动和形位参数对颗粒流润滑特性的影响规律[J]. 机械工程学报, 2009,45(7): 101―107.Wang Wei, Liu Xiaojun, Jiao Minghua, Liu Kun.Influences of several slider's parameters on particles flow lubrication characteristics [J]. Journal of Mechanical Engineering, 2009, 45(7): 101―107. (in Chinese)
[5]Timm K, Myant C, Spikes H A, Grunze M. Particulate lubricants in cosmetic applications [J]. Tribology International, 2011, 44(12)1695―1703.
[6]Heshmat H, Walton J F. High-temperature powderlubricated dampers for gas-turbine engines [J]. Journal of Propulsion and Power, 1992, 8(2): 449―456.
[7]常在, 杨军, 程晓辉. 砂土强度和剪胀性的颗粒力学分析[J]. 工程力学, 2010, 27(4): 95―104.Chang Zai, Yang Jun, Cheng Xiaohui. Granular mechanical analysis of the strength and dilatancy of sands [J]. Engineering Mechanics, 2010, 27(4): 95―104.(in Chinese)
[8]季顺迎, 孙珊珊, 陈晓东. 颗粒材料剪切流动状态转变的环剪试验研究[J]. 力学学报, 2016, 48(5): 1061―1072.Ji Shunying, Sun Shanshan, Cheng Xiaodong. Shear cell test on transition of shear flow states of granular materials [J]. Chinese Journal of Theoretical and Applied Mechanics, 2016, 48(5): 1061―1072. (in Chinese)
[9]Marzougui D, Chareyre B, Chauchat J. Microscopic origins of shear stress in dense fluid–grain mixtures [J].Granular Matter, 2015, 17(3): 297―309.
[10]Cundall P A, Strack A. A discrete numerical model for granular assemblies [J]. Geotechnique, 1979, 29: 47―65.
[11]井国庆, 王子杰, 施晓毅. 多围压下三轴压缩试验与不可破裂颗粒离散元法分析[J]. 工程力学, 2015,32(10): 82―88.Jing Guoqing, Wang Zijie, Shi Xiaoyi. Ballast triaxial tests and non-breakable particle discrete element method analysis under different confining pressures [J].Engineering Mechanics, 2015, 32(10): 82―88. (in Chinese)
[12]焦杨, 章新喜, 孔凡成, 刘海顺. 湿颗粒聚团碰撞解聚过程的离散元法模拟[J]. 物理学报, 2015, 64(15):154501.Jiao Yang, Zhang Xinxi, Kong Fancheng, Liu Haishun.Discrete element simulation of impact disaggregation for wet granule agglomerate [J]. Acta Physica Sinica, 2015,64(15): 154501. (in Chinese)
[13]葛良, 桂南, 徐文凯, 闫洁. 波形滚筒内颗粒混合和导热分布形态特性的研究[J].中国科学: 技术科学, 2014,44(1): 62―70.Ge Liang, Gui Nan, Xu Wenkai, Yan Jie. Numerical study of the distribution characteristics on particle mixing and thermal conduction in wavy drum [J].Scientia Sinica Technologica, 2014, 44(1): 62―70. (in Chinese)
[14]Yang Y, Chen Y M, Wang J A. Exploring the contact types within mixtures of different shapes at the steady state by DEM [J]. Powder Technology, 2016, 301: 440―448.
[15]Li M, Yu H, Wang J, et al. A multiscale coupling approach between discrete element method and finite difference method for dynamic analysis [J]. International Journal for Numerical Method in Engineering, 2015,102(1): 1―21.
[16]Campbell C S. Granular material flows-an overview [J].Powder Technology, 2006, 162(3): 208―229.
[17]Meng F J, Liu K, Wang W. The force chains and dynamic states of granular flow lubrication[J]. Tribology Transactions, 2015, 58(1): 70―78.
[18]Hsiau S S, Shieh Y M. Fluctuations and self-diffusion of sheared granular material flows [J]. Journal of Rheology,1999, 43(5): 1049―1066.
RESEARCH ON MECHANICAL MECHANISM OF SHEAR DILATANCY FOR GRANULAR FLOW LUBRICATION
MENG Fan-jing1,2, LIU Kun3, WU Hua-wei2
(1. Department of Mechanical and Electrical Engineering, Henan Institute of Technology, Xinxiang 453003, China;2. Hubei Key Laboratory of Power System Design and Test for Electrical Vehicle, Hubei University of Arts and Science Xiangyang 441053, China;3. Institute of Tribology, Hefei University of Technology, Hefei 230009, China)
Abstract:Granular lubrication is a new lubrication method which can be used in severe environments, while the shear dilatancy characteristics of granular flow lubrication are the key obstacles that remains unsolved in fully setting up the granular lubrication theory. A granular flow lubrication analytical model based on the discrete element method is proposed. Using this model, the variation rule, speed adaptation mechanism, and stress adaptation mechanism of shear dilatancy in the lubrication startup phase are studied, and the mechanical mechanism of shear dilatancy was obtained by analyzing the influence of shear velocity and pressure load on the shear dilatancy. The results show that in the lubrication startup phase, the shear dilatancy rate of granular flow lubrication system exhibits a quadratic polynomial law with time steps. The shear dilatancy amplitude of granular flow lubrication system all gradually decreases with increased shear velocity and pressure load, respectively. The shear dilatancy characteristics of granular flow lubrication is closely related to the fluctuation of force chains between granular lubrication medium. When the fluctuation amplitude of force chains increases, the variation amplitude of shear dilatancy for granular flow lubrication system also increases, and vice versa. Distribution of contact forces between granular lubrication medium is found to obey a power law, and the strength of force chains increases with the increased absolute value of shear dilatancy rate.
Key words:granular flow lubrication; shear dilatancy; discrete element method; mesoscopic force chain;mechanical mechanism
中图分类号：TH117
文献标志码：A
doi:10.6052/j.issn.1000-4750.2017.04.0282
文章编号：1000-4750(2018)07-0236-09
收稿日期：2017-04-11；修改日期：2017-07-28
基金项目：国家自然科学基金项目(51605150，11472096)
通讯作者：孟凡净(1981―)，男，山东人，副教授，博士，硕导，主要从事机械摩擦学力学机制研究(E-mail: mengfanjing0901@126.com).
作者简介：
刘 焜(1963―)，男，陕西人，教授，博士，博导，主要从事摩擦学与表面工程研究(E-mail: liukun@mail.hf.ah.cn);
吴华伟(1979―)，男，湖北人，副教授，博士，硕导，主要从事有限元分析与优化设计研究(E-mail: wuhuawei1130@163.com).