计算机模拟在工程中的应用姓名：曾亚维班级：研机械工程1402班学号：S14030662015年6月30日滑动摩擦的有限元仿真研究引言计算机数值仿真是数值计算理论与计算机相结合而发展起来的一种新技术，在众多工程领域发挥了重要的作用。

随着计算机数值仿真技术的发展和计算机辅助工程(CAE)分析软件的出现，在摩擦学领域，摩擦磨损仿真技术的研究也正蓬勃兴起。

数值仿真的实质是将连续的过程离散化，从而使复杂的动态问题转化成较简单的准静态问题。

该技术应用于复杂的摩擦学磨损问题是十分有效的，并具有广泛的应用前景。

目前人们己经将其用在摩擦副的接触分析、寿命预测、优化设计、状态监测及建立相应的耐磨损强度校核及设计准则等方面，发挥了重要的作用。

采用基于计算机的数值仿真方法进行各种方案的计算，则可以使大难度的工作轻易完成，并可以大大缩短新技术的开发和生产的时间周期。

而且对于不同条件(如材料、载荷、运转速度、温度、工作条件等)下的同类问题，只须改变数学模型中相应的参数，就能够获得对应的仿真结果，给工程设计分析带来质的飞跃。

1、国内外摩擦磨损仿真研究国际上许多学者都对摩擦学仿真进行了大量的研究，为摩擦学在新技术时期的发展做出了突出贡献。

如1977年首次召开专门涉及磨损领域的国际材料磨损会议(Wear of Materials，简称WOM)，专门讨论和交流磨损研究的进展，并且每隔两年召开一次，现在已经召开了14届国际材料磨损会议。

在第14届会议上，约翰·霍普金斯大学的Mark Robbins教授作了“摩擦学的计算机仿真”特邀报告。

Robbin教授指出，利用摩擦学仿真，可以用理想化的试验来对几何参数、化学参数和滑动条件等进行全面控制，使得对摩擦、润滑和磨损等每个变量的影响的控制成为可能；同时，他就违反人们直觉的几种摩擦学现象在纳米尺度上的摩擦学仿真进行了有益而成功的探索，运用表面吸附膜和2个表面间存在的“第三体”解释了Amonton摩擦定理及其它的摩擦学现象等。

计算机仿真技术运用到摩擦学领域，国内的部分专家学者在80年代就开始对此做了深入的研究。

但由于对摩擦副的影响因素较多和磨损机理的复杂性，国内因此在摩擦学的仿真研究方面，多集中在对具体的工程机械摩擦构件上。

如大连大学的江亲瑜教授、大连海事大学的严立教授、西安交通大学方亮教授、洛阳工学院的林纲等人开展了积极的探索，做了大量有意义的研究工作。

江亲瑜教授等人在初期进行有关数值仿真技术在磨损研究中的应用方面的研究，提出了一个研究磨损问题的数值仿真模型，并探讨了该模型在多种磨损问题中的应用前景。

严立教授等人提出一种利用智能仿真技术求解磨损问题的方法。

其思路是在磨损试验、机理研究和系统分析的基础上，采用模糊推理和智能仿真方法对磨损状况进行演绎。

他并不勉强去建立磨损与各影响因素之间的显式数值关系，而是利用有关基础理论、试验资料和经验知识，建立了包含有灰匣和黑匣以及推理机构的框形综合仿真模型；以必要的试验数据及应用工况为初始条件，进行仿真试验并对结果的可信度进行检验，以此评价系统的耐磨特性，预测其发生故障的可能性及使用寿命。

2、目前存在的问题摩擦学系统的复杂性导致人们很难对摩擦副的磨损性能给出定量的预测，尽管研究者们在摩擦磨损数值仿真方面所做的工作，取得了不少成绩，但同时也有一定的问题存在。

许多研究者通过实验分别从不同角度研究了不同磨损机理形式在平滑表面上的磨损模型，但都不能很好地得到验证。

如Podar等学者建立锥形接触时考虑表面形貌的磨粒磨损模型时，用有限元和解析法两种方法对磨损进行了计算，但两种方法计算结果相差较大。

Meakawa，K等人对接触表面的摩擦机理、微观粗糙度和微观形貌进行的模拟和分析也局限于微观机理方面，没有从宏观上把握磨损特性。

在实验室通过磨损试验对摩擦学的性能进行仿真研究的研究者们，只考虑即定参数、即定的磨损阶段的摩擦学行为；在建立仿真模型过程中对一些影响因素的处理有待商榷；对同时改变多参数如材料、载荷、速度等没给予深入研究。

一般说来，在某种特定的条件下通过试验手段所得到的磨损特性，在条件有了一些变化时，就发生很大的改变。

如果想了解该状态下的磨损情况，人们就不得不做试验或者是积累数据以确定磨损曲线等。

然而，达到这个目的需要大量的试验数据的支持。

3、ANSYS对销盘接触问题的仿真实现1970年，JohnS wanson博士发现计算机模拟工程商品化的前景，于是创建了ANSYS公司。

30年来，ANSYS公司致力于设计分析软件的开发，不断吸取新的计算方法和计算技术，使得这一分析软件功能更强大，使用更加便利。

ANSYS软件是融结构、热、流体、电磁、声学于一体的大型通用有限元软件，可广泛的用于核工业、铁道、石油化工、航空航天、机械制造、能源、汽车交通、国防军工、电子、土木工程、生物医学、水利、日用家电等一般工业及科学研究。

该软件提供了不断改进的功能清单，具体包括：结构高度非线性分析、电磁分析、计算流体力学分析、设计优化、接触分析、自适应网格划分及利用ANSYS参数设计语言扩展宏命令功能。

3.1 ANSYS接触分析接触问题是一种高度非线性行为,需要较大的计算资源，为了进行切实有效的计算，理解问题的特性和建立合理的模型是很重要的。

接触问题存在两个较大的难点：其一，在求解问题之前，不知道接触区域，表面之间是接触或分开是未知的，突然变化的，这随载荷、材料、边界条件和其它因素而定；其二，大多的接触问题需要计算摩擦，有几种摩擦和模型供挑选，它们都是非线性的，摩擦使问题的收敛性变得困难。

3.2 ANSYS销盘接触分析步骤ANSYS支持三种接触方式：点一点，点一面，面一面，每种接触方式使用的接触单元适用于某类问题。

为了给接触问题建模，首先必须认识到模型中的哪些部分可能会相互接触，如果相互作用的其中之一是一点，模型的对立应组元是一个结点。

如果相互作用的其中之一是一个面，模型的对应组元是单元。

有限元模型通过指定的接触单元来识别可能的接触匹对，接忽王单元是覆盖在分析模型接触面之上的一层单元。

就销盘接触模式的ANSYS 分析，摩擦表面微凸体放大到半球体接触分析，本文采用点一面接触单元；而对销盘实体接触分析则采用面一面接触单元。

ANSYS 分析不同的接触单元时，由于各接触单元形状以及节点不同,因此有不同的分析过程。

但其基本分析步骤是一致的，就销盘接触所用到的接触单元模式讲，为了保证整个仿真过程的实现和便于采用ANSYS 参数化设计语言(APDL ，ANSYS Parametric Design Language)，需制定求解步骤，并通过实际求解过程来完善，以满足在输入参数改变的情况下快速获取所需要的接触应力和变形值，其设计步骤流程图如图1所示。

图1 ANSYS 对销盘仿真步骤3.3 ANSYS 计算仿真实现采用有限元法对于求解球体接触问题，分析变载荷作用下的变形、应力分布输入初始化参数并建立实体模型-接触体几何形状、相对位置；-材料弹性模量E ，泊松比μ；生成有限元模型（/PREP7）定义并生成接触单元定义边界条件、加载做结构静态分析（/SOLU ）剥离接触区，获取法向相接触压力分布（/POST1） 拾取接触实体的变形值（/POST1）仿真结果结束仿真？ 退出改变实体参数 是否等有比较好的适应性，亦能更真实地反映接触应力的分布情况和输出所需要的参数值等。

对上述销盘接触实体问题利用ANSYS软件分别进行求解。

为了实现接触应力仿真的准确性，对单峰接触且粗糙峰顶的形状都是球体的微凸体，在ANSYS中建模较为困难，所以放大到实体两球体之间的接触，以求看到真实的、可比性的实体接触情况。

3.3.1 两球体模型仿真实现销盘试验所用材料为45号钢，取两球体为同种材料，45号碳钢其弹性模量为200GPa，泊松比为0.3，假设上球体半径R1=10mm，下球体半径R1=15mm,集中载荷作用于上球体。

图2为两球体接触网格有限元模型，用4节点的Quad node 42实体单元，网格采用映射划分，接触单元采用点一面接触方式。

由于对称关系，可取1/4模型进行仿真计算。

同时该问题也是一个轴对称问题，也可按轴对称模型求解。

有限元分析ANSYS结果的精度取决模型网格尺寸。

在划分网格时，接触部位的网格尺寸要小于接触椭圆短半轴的长度，当接触部位的网格尺寸小于接触椭圆短半轴的50%时，有限元分析结果己足够精确。

在本模型中椭圆短半轴和长半轴相等，即接触区为圆形。

模型总节点数1922，单元数1840，其中接触单元为40个，图2中右上角为接触区网格放大图。

由式(10),在P为500N时，代入参数计算得a=0.27mm，而接触区网格为0.1mm，因此已经满足分析的精度要求。

对整个接触仿真模型，为加快求解接触区外的其余部分可用较大的单元表示。

图2 两球体接触网格模型求解后进人ANSYS后处理器，可获取应力、接触区压力和变形等图形仿真结果，图3、图4分别是在P=4kN时仿真出的Y向应力分布云图和接触压力分布云图。

Y向应力分布云图清楚地反映出半球平面上的应力分布情况及大小；而接触区压力分布云图亦能帮助我们理解接触区压力分布情况，其形状为1/4椭圆，当取负值时，即为接触压应力。

由于变形图在图形窗口中需放大才能清楚显示,在此未提供。

图3 Y 方向应力分布云图 图4 接触面上压力分布云图3.3.2 变载荷下球体接触应力、变形计算与仿真结果比较分析两球体模型接触应力和变形计算方法有两种：即以上推导的赫兹接触应力公式、变形计算公式计算与有限元接触仿真分析计算。

由于接触体材料相同，且μ=0.3，简化处理得：3212212R R E R R P 23.1）（+=∆ （1） 322212120R R R R PE 388.0P ）（+= （2） 为了让理论计算结果与仿真结果做确切的对比分析，在此通过改变载荷的大小分别对接触应力、变形做了11组数据。

从公式(1)、(2)可以看出，力与接触应力、变形成指数关系，通过编程计算理论值与仿真值并在Matlab 中实现曲线的三次拟合，得图5和图6。

两图实线为计算拟合线，从中看出，仿真值与计算值存在一定的误差，但其趋势是一致的。

表一为P=4kN 和5.5kN 的计算与仿真数值对比情况，选择户4kN ，可在图4中看到接触压力分布云图，而P=5.5kN ，是在所做的11组数据中误差最大。

从产生的误差分析，存在两种情况，一是从材料上分析，当所加载荷超过材料弹性极限时，材料将会产生塑性变形，继而产生残余应力，由于残余应力的存在，致使接触面的分布压力将发生变化；二是从球体物理模型上分析，误差主要包括单元本身、单元离散、几何、边界范围和计算模型处理方法等的误差。

对第一种情况处理，需对赫兹接触压力分布予以修正，以满足工程分析的需要；对第二种情况，可以使用更多的或更精确的单元、更合适的接触刚度和更精确的几何实体模型解决。