ansys磨损分析.doc

用有限元的方法模拟滑动摩擦磨损用有限元的方法模拟滑动摩擦磨损摘要摘要磨损往往是影响产品寿命的一个主要因素因此磨损预测就成为工程的一个重要部分 这篇论文介绍了用有限元软件 ANSYS 来模拟磨损的方法用线性磨损定律和欧拉解析积 分提出了一个模型化的模拟程序然而， 还要考虑保证模型的正确性和数学方法的收敛性 分别用实验和有限元的方法分析了球形 pin-on –disk 系统在没有润滑条件下的接触问题，使 用了 Lim 和 Ashby 磨损图来区分磨损机理在给定几何尺寸和载荷的条件下，可以用有限 元的方法模拟磨损，得到磨损率对滑动距离的对应关系有限元软件 ANSYS 非常适合解 决接触问题和磨损模拟实际磨损率的分布范围在±40-60%的界限内会导致磨损模拟结果 相当大的偏离因此这些结果必须在一个相对的值上进行估测，从而比较不同的设计 关键词：关键词：磨损模拟；FEA；磨损试验；接触温度1．． 绪论绪论 摩擦副之间最可靠的摩擦学行为的知识可以通过做磨损实验来获得然而，当特别是 设计改变时需要在日常的内部程序基础上进行迅速的估测已经进行了大量的研究工作来 帮助设计者实现这一步 已经证实一个给定系统滑动磨损的主要参数是接触载荷和相对滑动速度。

速度由机构 运动来决定系统载荷怎么影响接触应力是很复杂的一个问题第一个分析两个弹性实体 接触应力的人是赫兹他认为接触体是弹性的，接触部分为椭圆形，而且没有摩擦的这 些假设被用在接触应力的计算中磨损发生在机械构件相互接触时一个重要的实际问题 是在给定的时间里有多少的材料损失由于功能和加工误差等表面的形状是不同的而且 会因为磨损和弹性变形而改变因此压力的分配就依赖于这些条件有限元的方法是一个 通用的工具来解决应力应变的问题这篇论文使用有限元软件 ANSYS5.0A 分析了接触压 力和磨损模拟 2. 磨损模型磨损模型 磨损过程可以认为是动态的，由许多参数决定，这个过程的预测可以看作是一个初始 值的问题从而磨损率就可以由一个总的方程来描述 dh/ds=f( 载荷，速度，温度，材料参数，润滑，….) h 为磨损深度（m） ，s 为滑动距离（m） 文献中可以查到很多磨损模型这些数学表达多 种多样，从简单的经验式到复杂的依赖于物理概念和定义的方程[1]常常包括特定的参数 和变量，只有在特定的情况下是正确的，在手册中是查不到的因此这些模型很少被用于 实际中磨损的估测Lim 和 Ashby 给出了钢板一个更广泛的在大范围载荷和速度下的磨损分类方法。

他们的 工作建立在简化磨损方程并在很多 pin-on-disk 实验数据的基础上进行调整通过这个工作 得到了一个磨损图，图 1 所示给出了磨损机理的轮廓和无量纲化磨损率 Õ， 作为无量纲 化压力 p 和无量纲化速度 v 的作用结果，定义如下（1） V 为磨损体积（m3） ，A 为接触面积（m） ，r0为接触半径（m） ，FN为载荷，H 为接触对中 较软材料的硬度（Pa） ，v 为相对滑动速度（m/s） 表一列出了 Lim 和 Ashby 使用的方程和 参数图 1 所示磨损图中的温度的分析是假定一个简单的温度沿一维方向流动的基础上进行的 更进一步的说就是迅速传播的温度对磨损起了一个很重要作用，热量分配系数 α12=0.5， 如果接触的温度达到 7000C，就会发生严重的氧化磨损在这个温度一下，磨损与载荷成 线性关系，不受速度影响最常用的模型是线性磨损方程 Q=kp，磨损体积率与载荷成比例关系这个模型被认为 是 Archard 的磨损定律，尽管它的基本形式首先由 Holm 发表这个模型是建立在实验观 察基础上的，用公式表示（2） 引入磨损率 K 使理论和实验的结果相吻合Holm 把它作为一个常量，来表述耐磨原子的 数量。

在 Archard’s 工作中，提出了微凸体相互作用导致磨损粒子产生的可能性[4]然而， 那并不是唯一可能的解释Lim 和 Ashby 通过计算认为分层和塑性磨损机理是主要的对 于刚，他们建议用以下的值：然而，对于一个接触实际的 K 值需要通过实验的方法得到对于工程应用来说，相对于磨 损体积，对磨损深度更感兴趣这里 Archard 想通过接触面积 A[4]把方程（2）都分开，给 出H 为磨损深度（m） ，k 为空间磨损率，p 为接触压力磨损过程可以看作是个动态的过程， 它的预测是个初始值的问题从而磨损模型可以通过一个不同的方程来表述，性的情 况下，方程（2）可以用下面公式表示（3） 3. 有限元磨损模拟程序有限元磨损模拟程序 3.1 有限元理论有限元理论用有限元方法计算磨损的主要工作是计算接触应力工程结构被离散成许多单元， 单元与单元之间用节点连接起来在有限元中，可以把单元中一些物理量(位移， 温度等)通过多项式分段拟合来近似描述，通过结点位移来表示[5]可以同时使用不同的单元类型，复杂的载荷和边界条件在结构分析中，把自由的程度定义为节点位 移把各个单元按原来的结构重新连接起来，形成整体的有限元方程[D]{μ}={F} D 为结构刚度阵或总体刚度阵，{μ}为结构的节点位移列阵或变形矢量，{F}为结构节点的 载荷列阵。

节点的应力通过变形得出商用有限元软件 ANSYS 可以处理几种材料和非线 性问题， 例如塑性，粘弹性， 摩擦等[6]有限元磨损计算包括解决全面的接触问题，实际接触面积是未知的，所以分析为非线性 在此情况下，使用了点对面单元 有限元软件配备了多种工具，提高了非线性分析程序， 但选择参数时要格外细心 3.2 磨损模拟程序磨损模拟程序 磨损模拟程序的流程图包括一系列的结构的解决步骤和额外的计算，如图 2 所示初始的参数用来定义模型的几何尺寸，载荷，约束和磨损模型参数以及单元和材料数 据并开发了一些特殊的小程序生成 FE 模型和自动的定义载荷和约束每个几何和载荷 的情况都要很好的离散化模型中采用较多的单元会得到更精确的结果，但是会增加计算 时间和占用硬盘空间当得到了 FEA 迭代应力后，就可以确定接触区域然后可以确定每个接触单元的状态 用接触单元的结点坐标确定接触区域根据接触区域节点应力来确定压力的分配 用欧拉方法求磨损相对于时间的积分系统参数在磨损模拟的每一个载荷步中被假定 为常数，根据下面所描述的磨损模型，在每个单元节点上对磨损深度产生影响4）Δhj,n是在节点 n,j 上的磨损增量根据已知的应力分布，可以估计出节点的磨损增量 Δhj,n（m） 。

如果在一个磨损步中某个节点处的磨损增量过大，模拟的结果就可能变得不确 定接触区域就可能出现间隙因此就根据以往的经验事先假定和引入了一个允许的最大 磨损增量进行时间较短的模拟实验来调整这个值使它尽可能大最初的节点磨损增量经 过一个连续的时间增量 Δt（s）计算出来每个求解步 j 中的时间因数 Mj可以通过下式估 测Δhin,j,max为节点磨损增量 Δhin,j，n的最大值那个求解步 j 中实际的时间间隙被定为 MjΔt然后模型的几何参数改变，根据方程（4）把节点移到新的位置这种方法，没有 用连续的时间步长，提高了 FEA 计算速度在每个求解步后保存输出数据是很重要的这 就要求如果分析由于某种原因中断时，快速数据回顾和保存前一个步中的数据 3.33.3 FEAFEA 结果确认结果确认 也许确认 FEA 结果最可靠的方法是将它与已知的分析结果相对比ANSYS 软件也配备 了能量误差估计方法，建立在 FEA 结构分析结果是一个从一个单元到另一个单元连续的位 移域[6]为了得到更满意的应力，单元节点的应力被平均节点应力误差矢量相应的被估 测出来，作为单元和整个模型能量误差估计的基础当每一个单元能量误差都相等时，这 个离散化的模型才是最有效的。

3.43.4 Sphere-on-planeSphere-on-plane FEFE 模型模型 Pin-on-disk 结构如图 3 所示，用上面所示的 FEA 方法作了分析在这种情况下塑性 变形和摩擦对压力的分配被看作是可以忽略的一个顶端为球形半径 R=5mm 的 pin 通过一 个轴向均质的 Sphere-on-plane FE 模型来表达这个结构采用了两维实体结构单元 PlANE42接触面采用二维点对面接触单元 CONTACT48Pin 和 disk 都被认为是弹性模量 E=210GPa，泊松比 µ=0.3 的钢做成的两个负载为 FN=21N 和 FN=50N接触区单元的 X 方向 的尺寸根据载荷不同分别为 25μm 和 32.5μmANSYS 接触刚度参数定为 KN=5X107N/m为了检验模型的正确性，接触应力的分配分别通过 FEA 和赫兹公式进行了计算 （图 4 所示）E*=E/2（1-µ2）为常态弹性模量这个模型中忽略了弹性变形和摩擦用 FE 和用 Hertz 分别计算的结果相差不超过 5%4.4.实验程序实验程序 一个半径 R=5mm 的 pin 在 disk 上滑动，压力分别为 FN=21N 和 FN=50N。

disk 和 pin 分 别硬化到刚度分别为 HV=4.6GPa 和 HV=3GPaHertz 计算的最大的接触应力被假定为在弹性 范围内实验设备可以直线测量磨损深度和摩擦力矩（图 5） 实验中的滑动速度为 v=25mm/s 实验结果图 6 所示，做了两个载荷下的实验图 6（C）所示的磨损率由图 6 所示两个 实验的的平均磨损深度决定使用了如下方程体积磨损增量通过公式[7]计算i≥1 是样品点的数量，Δs=0.15m 为滑动距离增量Disk 的硬度比 pin 的硬度大平均磨 损率通过试验数据估测， 滑动距离为 s=3m 和 s=4.5m[图 6（c）]当滑动距离分别为 s=3m 和 s=5m，这些值和偏差分别为在 F=21N 时 K=（1.25±0.44）*10-13Pa-1和 K=（2.26±1.44）*10-13Pa-1，在 F=50N 时 K=（1.33±0.54）*10-13Pa-1 和K=（2.01±1.21）*10-13Pa-1所测得的摩擦系数的平均值为 fFr=0.7±0.2通过 Archard 建议的方法分析了接触产生的温度[8]假定摩擦热流量 qn（W/m2）没 有流入接触体内部。

qn=fFrpvP 为接触区域的平均接触压力，计算时没有考虑系统变形通过一个圆形的热源加热，计算出了两个接触体的平均和最大温度没有考虑辐射和转化使用了下面的公式[9]接触温度由[8]计算无量纲参数 Pe=0.5v=vr0/2a0为数计算所得温度没有超过 6K[图 6（d）].在载 荷 F=21N 和 F=50N 时，实验速度和压力分别为 v=0.2,…,1.6, p=0.007,…,0.27 和v=0.3,…,2.0, p=0.009,…0.37在平均赫兹压力的基础上，计算相对于初始情况更高的无量纲压力认为磨损机理为分层或黏着磨损，图 1， 得到了线形磨损规律，公式 （3） 5.5.磨损模拟结果磨损模拟结果 假定符合线性磨损定律，FEA 磨损模拟结果可以用磨损率相对于磨损距离的关系来表 示给定几何尺寸和载荷时，如果 ks 不变，磨损的深度不会改变 5.15.1 Sphere-on-planeSphere-on-plane 滑动接触滑动接触 用以上的模型进行 FE 磨损模拟，假定符合公式（3）的线性磨损定律模拟磨损时估 测空间磨损率 K=（1.33±0.54）*10-13Pa-1最大的磨损增量定为 Δhlim=0.1µm。

求解步磨 损增量通过公式（3）Δh=kpΔs 计算出来Disk 在实验中被假定为硬度较高，因此只有 pin 发生磨损 磨损率做为常数来处理对比实验结果和 FEA 磨损模拟结果，如图 9 所示粗线表示 在平均 k 值下的磨损曲线，细线表示磨损误差的影响在摩擦过程中。