摩擦模型（精选8篇）
摩擦模型 篇1
摘要：双机架平整机组湿平整轧制过程中影响摩擦因数的因素较多, 采用常规方法无法回归问题。在大量的现场试验与理论研究的基础上, 通过对辊缝中润滑油膜厚度与摩擦因数之间关系的定量研究, 建立了一套双机架平整机组湿平整轧制过程中的摩擦因数模型, 解出了双机架平整机组机架间的延伸率分配系数, 并将其应用到某1420连退平整机组的生产实践, 取得了良好的使用效果, 具有进一步推广应用的价值。
关键词：湿平整,摩擦因数,润滑油膜
0引言
在平整轧制过程中, 摩擦因数是轧制压力预报的基础, 其计算精度直接影响到轧制压力、板形、延伸率等参数的预报与设定精度, 对成品带材质量具有举足轻重的影响。对于双机架高延伸率湿平整轧制过程而言, 影响摩擦因数的主要有来料因素 (包括钢种、规格等) 、轧制参数 (包括延伸率、轧制速度等) 、轧辊参数 (包括轧辊原始粗糙度、换辊后轧制公里数等) 以及平整液参数 (包括粘度等) 等4类因素。上述4类因素对摩擦因数又往往是交叉影响、综合作用的。这样, 如何考虑上述因素的综合作用, 建立相应的摩擦因数模型就成为现场攻关的重点。以往现场对于摩擦因数模型的开发主要集中在冷连轧机, 而对于平整轧制尤其是双机架湿平整轧制过程中摩擦因数模型的开发则未有文献涉及。为此, 经过大量的现场试验与理论分析, 综合考虑到高延伸率湿平整轧制的工艺特点, 充分利用现场实际数据, 通过对辊缝中润滑油膜厚度与摩擦因数之间关系的定量回归分析, 首次建立了一套双机架湿平整轧制过程中的摩擦因数模型, 并将其应用于国内某1420连退平整机的生产实践, 取得了良好的使用效果。
1摩擦因数模型的建立
1.1 摩擦因数影响模型的构造
在湿平整轧制过程中, 变形区润滑油膜的厚度与摩擦因数之间存在着一一对应的关系, 而来料因素、轧制参数、轧辊参数以及平整液参数等主要是通过影响润滑油膜的厚度来改变摩擦因数的。这就是说, 如果能够通过现场实际数据回归出润滑油膜厚度与摩擦因数之间的定量关系, 实质就已经建立起了双机架湿平整轧制过程中的摩擦因数机理模型。根据相关文献[1,2], 构造出湿平整轧制过程中润滑油膜厚度对摩擦因数的影响模型:
μi=a+beBξξ0i 。 (1)
其中:μi为双机架平整机第i机架的摩擦因数;a为干摩擦影响系数;b为液体摩擦影响系数;Bξ为摩擦因数衰减指数;ξ0i为第i机架的润滑油膜厚度。根据式 (1) 可以知道, 只要求出a、b、Bξ就可以得到相应的湿平整轧制过程中的摩擦因数模型。
1.2 摩擦因数模型参数的求解
根据基本轧制理论[3]可以知道, 对于双机架平整机而言, 在轧制过程中任意第i机架轧制压力可表示为:
Pi=F0 (hi, εi, σ1i, σ0i, k, μi, B, Ri) 。 (2)
其中:hi为双机架平整机第i机架出口带材厚度;εi为双机架平整机第i机架的延伸率;σ1i为双机架平整机第i机架的前张力;σ0i为双机架平整机第i机架的后张力;k为带材的变形抗力;Ri为双机架平整机第i机架工作辊半径;B为带材宽度。
与单机架平整机及冷连轧机不同, 双机架平整机在机架之间没有测速仪与测厚仪, 往往只能实测出总的延伸率, 而延伸率在各个机架之间的实际分配则是未知的。为此, 引入机架间延伸率分配系数ζ这一变量, 同时用ε0来表示总延伸率, 那么有:
ε1=ζε0 。 (3)
undefined。 (4)
h1=H0 (1-ε1) =H0 (1-ζε0) 。 (5)
其中:H0为来料厚度。
undefined。 (6)
将式 (3) ～式 (6) 代入式 (2) , 又可以得到双机架平整机轧制过程中任意第i机架轧制压力为:
Pi=F1 (H0, ε0, k, B, Ri, ζ, σ1i, σ0i, μi) 。 (7)
考虑到湿平整轧制过程中带材的前滑值比较小, 因此可以认为轧辊表面线速度vr与带材的入口速度v0近似相等, 这样湿平整轧制过程中, 润滑油膜厚度[4,5]可以表示为:
undefined。 (8)
其中:αi为双机架平整机第i机架的咬入角;Krs为轧辊传递到带材上的复印系数;η0为润滑剂的动力粘度;θ为润滑剂的粘度压缩系数;Rar0i为双机架平整机第i机架的工作辊原始粗糙度;Li为双机架平整机第i机架工作辊换辊后的轧制公里数;BL为工作辊粗糙度衰减系数;vi为双机架平整机第i机架带材出口速度;krg为第i机架带材的变形抗力。
由于双机架平整机的宽展量比较小, 所以可认为:
v1h1=v2h2 。 (9)
将式 (5) 、式 (6) 代入式 (9) 得:
undefined。 (10)
综合式 (1) 以及式 (8) ～式 (10) , 可以将摩擦因数表示为:
μi=F2 (H0, ε0, k0, B, Ceq, Ri, ζ, v2, Rar0i, Li, η0, θ, a, b, Bξ, σ1i, σ0i, Pi) 。 (11)
其中:k0为来料的变形抗力;Ceq为材料的碳当量。
实际上, 对于特定的双机架平整机组而言, 根据其工艺及设备特点, 在一个特定的平整轧制过程中, H0、ε0、k0、B、Ceq、Ri、v2、Rar0i、Li、σ1i、σ0i、Pi等参数是已知的。如果工艺润滑制度也确定, 那么η0、θ的值也是已知的。这样, 未知的参数主要有a、b、Bξ、ζ 4个参数。而式 (7) 、式 (11) 代表4组方程, 显然, 4组方程解4个未知数, 通过一定的算法是可以求出解的。基本求解思路如下:
(1) 收集一组特定的轧制工艺参数H0, ε0, k0, B, Ceq, Ri, v2, Rar0i, Li, σ1i, σ0i, Pi及润滑工艺参数η0、θ。
(2) 令X={a, b, Bξ, ζ}, 并给定X的初始值X0。
(3) 将X0代入到方程组 (11) , 求出摩擦因数的计算值μ′i。
(4) 将μ′i代入到方程组 (7) , 求出轧制压力的计算值P′i。
(5) 计算目标函数值undefined。
(6) 判断POWELL条件是否成立?如果不成立, 调整X的设定值, 直到POWELL条件成立为止。
2摩擦因数模型的建立
针对国内某1420双机架平整机组, 为了建立湿平整轧制过程中的摩擦因数模型, 根据现场相关实际数据, 采用本文所介绍的相关方法, 回归计算出该双机架平整机组在湿平整模式下的摩擦因数模型:
μ=0.07+0.3e-1.625ξ0 。 (12)
利用式 (12) 作出摩擦因数与油膜厚度之间的关系曲线, 如图1所示。
3结论
本文在大量现场实验与理论分析的基础上, 充分考虑来料因素、轧制参数、轧辊参数以及平整液参数的交叉影响, 建立了一套适合双机架湿平整高延伸率轧制的摩擦因数模型, 给出了相应的模型参数计算方法, 并将其应用到国内某1420双机架平整机组的生产实践, 取得了良好的使用效果, 具有进一步推广应用的价值。
参考文献
[1]白振华, 周庆田, 窦爱民, 等.冷连轧高速轧制过程中摩擦因数机理模型的研究[J].钢铁, 2007, 42 (5) :47-50.
[2]白振华.冷连轧机高速生产过程核心工艺数学模型[M].北京:机械工业出版社, 2009.
[3]白振华.薄带平整轧制时轧制压力模型的研究[J].机械工程学报, 2004, 40 (8) :63-66.
[4]W.L.罗伯茨.冷轧带钢生产[M].王延博, 译.北京:冶金工业出版社.1985.
[5]李小玉, 顾正秋.轧制工艺润滑[M].北京:冶金工业出版社, 1981.
摩擦模型 篇2
关键字出口模式贸易摩擦引力模型中德比较
一、 背景
商务部发布的《中国对外贸易形势报告》（2013以及2014年秋季）显示，中国俨然成为全球贸易保护主义的最大受害者。尽管世界经济增长有所起色，但国际市场竞争更加激烈，全球贸易摩擦依然高发，贸易摩擦形势依然严峻复杂。在国际贸易保护主义回潮的背景下，针对中国产品的贸易摩擦有增无减。2014年前三季度，共有21个国家（地区）对中国出口产品发起救济调查75起，同比增长17%，其中不少摩擦针对中国战略性新兴产业《2015年中国对外贸易发展环境分析》，中国商务部网站，http：//www.mofcom.gov.cn/article/i/dxfw/gzzd/201411/20141100810825.shtml，访问日期：2014年11月26日。。同时，发展中国家对中国产品发起的贸易救济案件数量和金额均超过发达国家。据英国智库经济政策研究中心的“世界贸易预警”项目监测，2008年国际金融危机爆发以来，全球40%的贸易保护主义措施针对中国。《述评：中国的继续强大有利实现“人类命运共同体”》，新华网，http：//news.xinhuanet.com/politics/2012-11/13/c_113680750.htm，访问日期：2014年11月20日。同时，对华贸易保护手段的打击对象正在从具体产品到整个产业，从劳动密集型产业到新能源与高科技产业，且有越来越多的发展中国家迈入对华贸易保护阵营。
中国在日益成长为世界贸易大国的过程中，成为遭受贸易救济措施最严重的国家之一，而同为出口大国的德国，却鲜少陷入与其他各国的反倾销争端。虽然中德两国的出口规模相当，以2012年为例中国出口总额22 213亿美元，德国17 826亿美元，且德国的出口占GDP比重高达52%，中国约27%，但两国却有十分不同的遭遇，相比较德国的出口，中国出口的产品有什么不同？是否这种差异性导致中国的产品遭受全球贸易壁垒？德国对中国发展贸易强国有什么借鉴？
本文选择德国作为中国出口产品的比较对象，试图从中德比较的角度分析中国出口产品与德国的差异，揭示中国产品频遭贸易壁垒的内在原因。本文利用2013年中国和德国出口的HS92六分位贸易数据，以引力模型为基础进行计量实证 。
二、 相关文献综述
对于中国出口频遭贸易壁垒的问题，已有大量文献进行分析，主要集中在对反倾销的研究上。谢建国做了美国对华贸易反倾销实证研究，显示美国国内工业产出的波动与对华贸易逆差显著提高了美国对华的反倾销调查频率。王孝松和谢申祥通过多年度跨国跨行业的数据，发现中国频遭反倾销的重要因素有出口激增、人民币贬值等，其他还有一国反倾销的总体状况、中国总体反倾销能力、各国的关税减让幅度以及各国加入FTA的数量等。
中国出口的大幅提升，引起了国内外学者的广泛关注，有不少针对中国出口的研究。区别于上述文献对外部环境等的分析研究，中国出口自身存在的特殊性成为造成贸易问题的内因。杨汝岱等研究发现，中国出口产品的技术复杂程度不断上升，出口结构与OECD国家逐渐接近，这种有限的赶超对于经济增长具有重要的意义。Amiti等认为，在中国出口贸易中剔除加工贸易后，中国出口结构并不存在优化趋势，中国出口仍然集中在劳动密集型行业，产品技术复杂度也不存在明显上升。由此，施炳展基于引力模型对美国的进口数据进行计量，发现中国出口数量大、价格低的特殊性。
近些年，传统引力模型已被很好地扩展和补充，部分学者扩充了引力模型的解释变量，如Andersen等。随着新新贸易理论的发展，贸易广度与深度的概念被引入传统模型，地理距离等因素不仅关系到贸易的发生，更对贸易的数量和贸易的价格有重要的影响，Helpman、Bernard、Baldwin及Harrigan等对此进行了深入的研究和推进。
本文将利用2013年最新数据，基于扩展的引力模型，分析新的形势下中国的出口与贸易摩擦。依据引力模型的最新发展，研究2013年中国与德国分别出口产品的数量、价格与引力模型变量之间的关系，将引力模型、贸易方式与贸易摩擦相结合。比照德国的出口，通过对比研究中国出口是否仍存在其特点，且这种特点是否构成中国频遭贸易保护的内在原因。
三、 模型、指标与数据
1. 计量模型
在引力模型的基础上扩展出回归计量模型，分别将双边贸易的数量与价格对双边相应变量进行回归，如（1）式：
lnxik=α+β1lngdp+β2lngdpper+β3lndist+γContrlk+uik（1）
式（1）中，i表示进口国，k为HS92版本六分位产品的标识;lnxik为进口国i进口k产品的数量或价格;核心变量包括lngdp、lngdpper、lndist，分别代表进口国的经济发展水平及规模、出口国与进口国的地理距离;Controlk是一系列控制变量，包括进口国与出口国是否接壤、是否共同语言、进口国的地理特征等;uik代表其他影响出口价格或数量的因素。本文将分别以中国和德国作为出口国，运用式（1）进行回归分析，进而比较中德出口的差异性。
2. 指标数据
表1列出式（1）中所有指标的含义及来源，并对其进行统计分析。
其中，进口国的需求水平影响其进口需求，用lngdpi、lngdpperi代表进口国的经济发展规模和水平，一般来说，经济规模大的国家较规模小的国家，进口商品的品种和数量相对较大，且一国经济规模越大，本国生产产品的种类和数量越多，进口商品需要有较强的价格优势才能与本地商品竞争，因此价格对lngdp的预期为负，数量对lngdp的预期为正。而经济发展水平越高，对高品质的进口品需求越大，进口价格也越高，因此，价格对lngdpper的预期为正。另外，传统引力模型常用的指标还包括衡量双边贸易成本的一系列指标，一般来说贸易成本越高，贸易数量越低，贸易价格越高，数量对lndist的预期为负，价格预期为正。
从基本统计分析看，中国出口价格均值096，低于德国出口价格均值122;中国出口数量均值1083，高于德国出口数量均值931;中国出口的价格/数量明显低于德国的水平，初步可以看出中国的出口与德国有明显不同，中国出口的平均价格较低，对此进行进一步的计量分析。
四、 实证结果分析
（一） 初步回归
分别采用2013年中国和德国出口数量与价格的数据进行初步分析，回归结果汇总在表2中，左右两边分别是根据中德两国各自的出口价格和出口数量对解释变量所做的回归。其中回归（1）是针对核心变量的回归，包括lnwp、lngdp、lngdpper、lndist;回归（2）则加入其他控制变量，有contig、landlocked、langoff等。
首先，关注出口价格的回归下lngdp变量的系数。两种回归方法下德国的lngdp系数都显著为负，而中国出口数据得出的该系数为正。德国的出口与理论预期相符合;而中国则不然，经济规模越大的国家，对中国商品的进口价格越高。但从基本统计分析中可以看出，中国出口商品的总体价格水平（即平均价格）低于德国，说明两国出口商品的种类有很大差异。事实上，由于以美国为主的大部分高GDP国家的制造业产业主要是高科技和尖端制造，德国出口与高GDP国家之间的竞争关系明显，在这些竞争产品上须有一定的价格优势才能进口，因而德国的出口价格与进口国的GDP呈逆向关系;而中国出口的商品中有较多的纺织服装和低技术的机电产品，仍以低端制造为主，与高GDP国家本身的产品竞争不大，故中国出口价格与进口国的GDP水平并非负相关。
从出口数量对lngdpper的回归结果看，两种回归方法下德国的lngdpper系数都显著为正，中国出口数据得出的lngdpper系数显著为负。这说明经济发展水平越高的国家对德国商品的进口越多，经济发展水平越低的国家则进口更多的中国商品。
所有代表贸易成本的变量contig、landlocked、langoff在中国和德国的样本回归中均符合预期。再观察lndist，德国样本下价格对lndist的系数均显著为正，数量对lndist的系数显著为负，符合理论预期;但中国样本计量结果显示，价格对lndist的系数为负，距离越远，价格越低，同时，中国出口产品数量对距离的反应弹性较小，表明在控制贸易成本增加相同的情况下，中国出口数量减少量远低于德国，价格甚至更低，反映了中国出口的优势仍是低价。
通过初步回归可以比较得出，德国和中国虽然都出口大量产品，但两国在出口的产品种类及定价方面存在差异：德国以高端制造类产品为主，而中国出口中纺织、机电等低技术产品比重较大;在定价方面，撇去两国等条件可以看出，中国价格偏低。
（二） 稳健性检验
1. 按不同发展水平国家
初步回归分析中，比较了中国与德国出口世界各国的数量和价格，发现中国与德国的出口方式及产品存在差异，中国的产品相对低质。如此，中国与德国产品的差异性在发达国家中应有更明显的表现。因此将原样本分成高收入国家、中等收入国家和落后地区，按照世界银行划分标准。进一步回归。
表3中，回归（1）是针对高收入国家数据分别就价格和数量所做的回归，回归（2）是中等收入国家，回归（3）是落后地区。结果显示，在各类型国家德国lngdp对价格的系数仍为负，中国在所有类型的地区中都是正的。这与中国是制造业大国，且与出口商品种类有关，既有纺织、服装、塑料制品等劳动密集型产品，也有大型机电和高新技术产品，与不同类型国家都能有差异化的产品贸易。同时，中国出口数量在发达地区数据回归下对lngdpper的系数仍为负，且反应程度远远高于对其他类型国家的反应，表明中国在与发达国家的贸易中，价格敏感度较大，仍以低价取胜。德国的出口数量对lngdpper的系数仍为正，从发达国进一步考虑贸易结构，用中德出口产品的技术层次分析稳健性。前述的初步回归看出中国出口产品以低端制造为主，而德国出口多高技术水平的产品，且德国出口品的价格高过中国。因此仅对高技术水平产品的出口数据做回归，结合贸易结构和中德出口比较。
表4是针对高新技术产品所做的回归。按《中国高新技术产品目录》，选用HS92六分位产品中的高新技术产品样本。回归结果中，中德两国lngdp对价格和数量的系数与初步回归相同，即两国出口的高新技术产品符合前述对所有出口品的分析，但lngdpper、lndist和contig对价格的系数变成了负值，lngdpper的影响与预期相反说明高收入国家进口的中国的高新技术产品与德国相比，并无技术或垄断优势。同时，以距离和是否接壤为代表的贸易成本的增加，反而降低了中国产品的进口价格，更加说明中国出口的高新技术产品主要是低价竞争策略。事实上，中国的高新技术产品制造大多缺少核心竞争力，以简单加工制造为主，价值增加值较低。由此，按产品技术类型的检验结果也进一步说明了基本回归的结论。
3. 按不同时间数据
此外，由于金融危机后世界经济形势有了新的变化，国际贸易壁垒也呈现出新的特点，危机后中国和德国出口方式间的差异性可用2009年数据检验结论的稳健性。表5中对2009年数据的回归得到的结果仍然与上述结论相同，从时间维度验证了结论的稳健性。
五、 结论
本文从出口价格和数量两个角度比较了中国与德国的出口差异性，利用2013年两国HS92版六分位出口数据，基于引力模型进行实证分析。通过引力模型的理论分析，可以看出中国出口的模型拟合在很多方面不符合预期，这可能就是中国在贸易出口中遭遇强摩擦的原因，而中国的低价商品以及以加工贸易为主的出口方式背后，隐藏着以下三个原因。
1. 生产率的提高
根据戈莫里和鲍莫尔的观点，一个国家生产率的提高往往会损害其他国家的整体福利，国际贸易可能会导致各贸易国之间的利益冲突，特别是当新兴贸易国家在全球市场上开始占有重要位置时，该国将发展出更多不利于发达国家的产业，因此发达国家为维持其贸易优势展开激烈竞争，导致两国间的贸易摩擦，这就是中国与很多发达国家之间高发贸易摩擦的根本原因。已有不少学者通过各种方法对中国参与全球性生产和贸易后全要素生产率的增长进行估够以较低的生产成本在国际商品市场始终保持低价的优势。
2. 产业结构调整
根据国际贸易理论，国际产业结构之间的结构性互补是国际贸易双方共赢的前提，而产业趋同则容易发生贸易摩擦。随着中国产业结构的逐渐升级，出口产品的技术含量也逐渐提高，中国与发达国家之间产业结构的趋同不断深化，同类产品之间的冲突和竞争也越来越激烈。中国近些年的主要贸易商品中，包括汽车及其零部件、有机化学品、装备制造业等所占份额越来越高，产业结构调整和升级使得中国和美国等发达国家在一定程度上产业结构趋同，造成了激烈的市场竞争，成为国家之间易摩擦产生的直接原因。同时，东南亚一些国家近年来重点发展的纺织业、加工制造业等，也构成与我国初级产品的竞争。因此，中国的贸易摩擦不仅由初级产品向工业制成品和高新技术产品转移，还在纺织等劳动密集产品上有所增加。
3. 产能过剩
长期以来的考核体制导致不少地方政府以追求GDP增长速度为首要目标，通过各种政策促进投资，以扩大投资规模，一些相关行业产能严重过剩。数据显示，2012年底，我国钢铁、水泥、电解铝、平板玻璃、船舶产能利用率分别仅为72%、73.7%、71.9%、73.1%和75%，明显低于国际通常水平，仍有一批在建、拟建项目，产能过剩呈加剧之势，并扩大到风电、光伏、碳纤维等包括新兴产业在内的19个行业。国内的过剩产能转移到国际市场，也成为一些国家对中国发起贸易保护的借口。以钢铁贸易为例，20世纪60年代美日、美欧之间频发的钢铁贸易就是因为美国钢铁产量下滑，而欧盟和日本的钢铁产量大幅上升;20世纪80年代开始，随着中国的钢铁产能急剧攀升，中国也成为发生钢铁贸易摩擦最多的国家。
参考文献：
[1] 谢建国.经济影响、政治分歧与制度摩擦：美国对华贸易反倾销实证研究[J].管理世界，2006（12）.
[2] 施炳展.金融危机后中国频遭贸易壁垒的内因分析：以中美贸易为例[J].财贸研究，2011（8）.
[3] 王孝松，谢申祥.中国究竟为何遭遇反倾销：基于跨国跨行业数据的经验分析[J].管理世界，2009（12）.
[4] 朱诗娥，杨汝岱.中国本土企业出口竞争力研究[J].世界经济研究，2009（1）.
[5] 闫克远.中国对外贸易摩擦问题研究.东北师范大学博士论文，2012.
[6] 施炳展.中国出口结构在优化吗：基于产品内分类的视角[J].财经科学，2010（5）.
[7] 国务院发展研究中心“区域协调发展和优化全国生产力布局”课题组.生产力布局的内涵及我国生产力布局存在的问题[J].发展研究，2014（12）.
[8] 林峰，戴磊，林珊.从国际服务贸易摩擦透视自由化谈判的利益差异：兼论中国服务贸易发展的战略选择.亚太经济，2014（11） .
（责任编辑：张晓薇）
摩擦模型 篇3
干摩擦出现在许多工程问题中,如制动系统[1,2]、加工机床[3]、轴承和齿轮[4,5]等,许多学者对其进行了研究,但在处理时几乎都将物体之间的摩擦因数μ视为常数,如单颖春等[6,7,8]利用库仑摩擦模型提出了平面接触和复杂接触运动条件下求解摩擦力的数值轨迹跟踪法,并将其用于计算凸肩复杂接触面内的非线性摩擦力;Hamraoui等[9]研究了两个滚筒在库仑干摩擦条件下的传热问题。文献[10]提出了Stribeck摩擦模型,认为干摩擦条件下两接触体间的相对切向速度与摩擦因数存在着一定的非线性关系。由此可见,用传统的库仑摩擦模型研究物体的干摩擦接触势必会影响到问题的真实求解,这对高精度的分析是不利的。
目前,在求解摩擦类接触问题时,大多借用商业有限元软件,其优点是省略了复杂的程序编制工作,可以得到比较理想的离散网格,这在处理复杂的接触问题时其优势尤其明显。但在处理摩擦类问题时,只能使用软件自带的几种摩擦模型,对Stribeck等非线性摩擦模型不能或难以处理,限制了其处理特定问题的能力。
本文所研究的干摩擦滑块模型源于工程中一系列接触场合的简化(接触目标体可假定为刚体),结合文献[10]中论述的Stribeck摩擦模型,用有限元方法分析滑块—皮带接触模型的平面应力特征,探讨接触体之间相对速度的变化对接触应力的影响。研究表明:用有限元方法准确求解Stribeck非线性摩擦模型在简化处理条件下的二维接触问题是成功的,接触体之间的相对速度对接触问题的求解有影响。
1 干摩擦滑块接触问题的分析模型
分析模型如图1所示,为使问题简单化,本文在准静态条件下研究干摩擦滑块接触问题,假定皮带刚性(即单边接触),此时滑块底边竖直方向上的位移为零。由于滑块在准静态时与左边固定约束的相对速度为零,所以该模型中的阻尼元件可以略去,将滑块与皮带之间的摩擦因数模型取为Stribeck模型,见图2。图2中,横坐标vr为皮带与滑块间的相对速度,纵坐标μ为摩擦因数,相对速度为0和0.5m/s时对应的摩擦因数分别为0.4和0.25,用点(0,μs)和(vm,μm)表示,此时点(0,μs)和(vm,μm)分别对应库仑摩擦和摩擦因数极小值点。
2 干摩擦滑块接触问题的有限元数值求解
2.1 模型的离散化
由于滑块模型比较简单,选取三节点三角形单元划分网格,共有35个节点、48个单元,节点与单元编号如图3所示。x轴方向、y轴方向分别表示滑块的水平方向与竖直方向,接触面为滑块与皮带的交界面,接触面法向即y轴负方向。滑块的几何参数为0.3m×0.2m,厚度取0.05m,单元几何尺寸为0.05m×0.05m。
2.2 求解步骤
2.2.1 选择位移模式
三节点三角形单元共有6个节点,由此构造的单元位移模式为
式中,u、v为节点位移;(x,y)为节点坐标;α1、α2、…、α6为待定系数。
2.2.2 建立单元刚度矩阵
由虚位移原理,将位移模式、几何方程与物理方程联合可推导单元e的刚度矩阵Ke:
式中,D为弹性矩阵;B为应变矩阵,B=[BiBjBm];t为滑块厚度;A为单元面积。
由于选取的是常应变三角形单元,由式(2)可得
S=DB=[SiSjSm]
对于平面应力问题,有
式中,E、ν分别为弹性模量和泊松比。
2.2.3 等效节点载荷
本文所讨论的模型涉及重力(体力)、分布摩擦力(面力)的移置,设单元单位体积的体力F=(X,Y)T,其单元边界上作用的面力
体力的移置
面力的移置
式中,N为形函数矩阵;s为积分面积。
2.2.4 整体分析
在同一整体坐标系下,由单元刚度矩阵Ke按角标编号“对号入座”叠加,最后集合整体刚度矩阵
2.2.5 求解总刚度平衡方程
求得总刚度矩阵K后,采用因子化法求解
式中,δ为位移矢量。
2.3 自制软件包与ANSYS的对比分析
在上述分析的基础上,对图3所示离散模型在库仑摩擦条件下分别采用自制软件包与ANSYS进行对比求解。求解相关参数如下:库仑摩擦因数μ=0.4,弹性模量E=200GPa,泊松比ν=0.3,滑块质量m=1kg。求解模型的工况如下:滑块上端中部作用1kN的集中力;滑块左端中部的弹簧力与滑动摩擦力大小相等,形成平衡;约束滑块底面节点y方向自由度。用FORTRAN语言编制的软件包与ANSYS分析的结果分别见图4和图5。对比自制软件包所得结果与ANSYS求解的结果可知,笔者自制软件包求解问题所得结果是可信的。
3 Stribeck摩擦滑块接触模型的有限元求解
用自制软件包研究Stribeck摩擦滑块接触模型时,取滑块和皮带之间的切向相对速度vr为一系列离散值(vr=0.2m/s,0.5m/s,0.7m/s),从图2 Stribeck摩擦曲线中将对应摩擦因数代入式(8)进行求解,求解的应力数据经MATLAB绘制成contour图,见图6～图8。
3.1 接触应力的contour图
由图6～图8可知,滑块在顶部中心受到竖直集中力的作用,受力点附近出现比较大的压应力,相应地在滑块左端中间位置,集中力作用使得该处附近的x方向正应力比较大。在Stribeck摩擦条件下求解滑块接触模型时,不同相对速度下的接触应力分布大致一样,但仍有一些细微的差别,局部地区的相对误差甚至不可忽略。库仑摩擦结果(已在图4、图5中给出)与Stribeck非线性摩擦模型计算结果的差值在相对速度为0.5m/s时达最大,31号节点x方向上的正应力相对误差达到100.08%,y方向正应力的最大误差在1号节点上,为82.73%,xy方向剪应力最大相对误差在3号节点上,为42.43%。
3.2 边界上接触节点结果对比
将滑块与皮带接触的边界节点1、6、11、16、21、26、31在不同相对速度下x方向位移(切向位移)、x方向正应力(接触面切向应力)、y方向正应力(法向应力)和xy方向剪切应力的分布进行对比,结果见图9。图9中,vr=0代表库仑摩擦时的情形,其余的3种速度分别代表各自相对速度下的情形。
Stribeck摩擦模型认为摩擦因数随着相对速度呈类似抛物线的变化,对于本文所选模型,最小摩擦因数在相对速度为0.5m/s时得到。从图9可知,接触面上的位移和应力分布也随滑块皮带间的相对速度作相应的变化,并且接触点上的压力分布也在0.5m/s时达到最小,由此可见Stribeck模型比库仑摩擦模型更能反映接触的真实情形。
4 结论
(1)Stribeck摩擦模型反映了相对速度对摩擦因数的影响,较库仑摩擦模型更能反映接触体间相对速度对接触求解的影响。本文干摩擦滑块接触问题算例在相对速度为0.5m/s时,求解结果偏离达到最大。
(2)本文成功解决Stribeck摩擦简化条件下的干摩擦滑块接触问题。
参考文献
[1]Tan Y,Wang Y.Transient Stabilization UsingAdaptive Excitation and Dynamic Brake Control[J].Control Engineering Practice,1997,5(3):337-346.
[2]Sinou J J,Thouverez F,Jezequel L,et al.FrictionInduced Vibration for an Aircraft Brake System-Part 2:Non-linear Dynamics[J].InternationalJournal of Mechanical Sciences,2006,48(5):555-567.
[3]Fletcher D,Smith L,Kapoor A.Rail Rolling Con-tact Fatigue Dependence on Friction,PredictedUsing Fracture Mechanics with a Three-dimen-sional Boundary Element Model[J].EngineeringFracture Mechanics,2009,76(17):2612-2625.
[4]He S,Singh R,Fellow A.Dynamic TransmissionError Prediction of Helical Gear Pair Under SlidingFriction Using Floquet Theory[J].Journal ofMechanical Design,2008,130(5):052603.1-052603.9.
[5]Liu G,Parker R G.Impact of Tooth Friction andIts Bending Effect on Gear Dynamics[J].Journal ofSound and Vibration,2009,320(4/5):1039-1063.
[6]单颖春,朱梓根.复杂接触运动下非线性摩擦力的求解[J].润滑与密封,2006(3):73-77.
[7]单颖春,郝燕平,朱梓根.平面二维接触摩擦力的轨迹跟踪计算方法[J].航空工程与维修,2002(3):48-50.
[8]单颖春,朱梓根,刘献栋.凸肩结构对叶片的干摩擦减振研究——理论方法[J].航空动力学报,2006,21(1):168-173.
[9]Hamraoui M,Zouaoui Z.Modelling of Heat Trans-fer between Two Rollers in Dry Friction[J].Inter-national Journal of Thermal Sciences,2009,48(6):1243-1246.
摩擦模型 篇4
摩擦力是影响机械系统动态性能的重要因素。建立精确摩擦模型并准确确定不同工况下的模型参数, 有助于开展精密机械系统摩擦补偿、轨迹跟踪和动态预测等工程应用研究[1,2], 具有重要的理论和实际意义。迄今为止, 已提出了多种描述摩擦特性的数学模型[3]。其中, LuGre摩擦模型[4]能够较准确地描述摩擦过程中的黏滑运动、摩擦滞后、预滑动位移、变最大静摩擦力等特性, 是目前较为完善的一种摩擦模型。该模型包含6个参数, 即鬃毛刚度系数和阻尼系数2个动态参数, 库仑摩擦力、静摩擦力、黏滞摩擦系数和Stribeck速度4个静态参数。由于6个参数之间存在耦合情况, 且模型中涉及的内部状态变量在实际系统中难以测量, 使得LuGre摩擦模型的参数特别是两个动态参数的确定非常困难。Canudas等[5]提出了一种基于局部线性化的两步辨识法, 利用最小二乘法对模型中的4个静态参数和2个动态参数分两步进行了辨识。此后, 基于两步法思想, 许多学者采用不同的优化方法对模型参数进行了辨识[6]。由于系统稳态阶段的速度-摩擦力曲线即Stribeck曲线能够测量, 且摩擦力与速度之间是线性关系, 易于建立形式简单的优化目标函数, 因此, 上述两步法与传统优化方法结合的辨识方法能够有效地辨识模型4个静态参数。然而, 对于模型中的2个动态参数, 其辨识目标函数复杂, 采用上述传统优化方法辨识时, 辨识结果对参数初值依赖性较强, 辨识收敛性难以保证, 往往得到的是局部最优解。Hensen等[7]提出了一种辨识LuGre模型动态参数的频域方法。该方法以随机噪声作为系统激励, 通过测得对应的频响函数, 来获得预滑动阶段两个动态参数的值, 但辨识效果不理想。Madi等[8]提出了利用区间分析方法对LuGre模型6个参数进行辨识。
本文基于两步法思想, 提出一种最小二乘法与边界误差估计方法[9]相结合的LuGre摩擦模型参数辨识方法。
1 边界误差估计方法
边界误差估计方法假定待估系统的测量误差是有界的, 并借助该有界误差和实测数据建立待估系统的理论输出值区间向量。通过比较仿真计算值与该理论输出值之间的包含关系, 最终获得待估系统参数的辨识值。
首先引入区间理论相关的基本概念和定义:
(1) , 定义为实数域R上的区间;
(2) [X]定义为实数域R上的区间向量, 包含n个区间元素[x1]~[xn];
(3) w ([X]) 定义为区间向量[X]中最大 (长) 区间元素的长度。
假设包含待估参数p的系统模型为yM (p) 。待估系统实测输出数据y与理论输出值yM (p) 之间的误差为e (p) , 则e (p) 必在给定的误差区间内。即, 待估系统理论输出值。因此, 当参数使得待估系统仿真输出值落在区间[Y]内时, 则表示该参数为满足误差要求的待估参数辨识值。所有满足误差要求的待估参数构成一个参数集S:
显然, 参数集S的长度越小, 表明待估参数的辨识值越接近参数真实值。
式 (1) 进一步改写成:
式 (2) 表示的是一个集反演问题。该问题的求解可以利用集反演 (set inversion via interval analysis, SIVIA) 方法[10?12]来实现。
SIVIA方法是基于给定的初始参数区间来求解的, 该区间包含了满足边界误差要求的所有参数值。将参数代入系统模型进行求解, 获得仿真输出值解集。将该解集与包含实测数据的可行集进行比较, 并将比较结果作为进一步切分参数区间的判据, 经过多次切分, 最后得到满足精度要求的全局最优参数区间。
当系统模型为常微分方程时, 不能直接将已知参数代入系统模型计算得到仿真输出值。此时, 需要在已知参数的基础上求解微分方程, 然后利用所得结果间接计算系统仿真输出值。以下是利用SIVIA算法辨识常微分方程中参数的过程[13]。
设任意T>0, 令I=[0, T]为实数域R上的任意闭区间。x为任意n阶向量, x= (x1, x2, …, xn) T。考察如下动力学系统:
式中, t为时间变量;u为系统外激励;f、g为状态函数和系统输出函数;x0为初始状态。
连续函数f满足唯一性条件。求解该常微分方程, 得到的解的下边界和上边界分别为
且
对式 (4) 和 (5) 进行数值积分, 得到变量x在不同时间节点的下边界和上边界, 然后回代到式 (3) 计算系统仿真输出值yM。
2 LuGre摩擦模型参数辨识
LuGre摩擦模型综合了Dahl模型和鬃毛模型的思想, 以状态变量z表示鬃毛的平均变形, 系统摩擦力Ff表示为
式中, v为两个接触面间的相对速度;σ0、σ1分别为鬃毛的平均刚度系数和阻尼系数;σ2为系统黏滞摩擦系数;vs为Stribeck速度;uC为库仑摩擦力;uS为静摩擦力;α (v) 为描述Stribeck负斜率效应的函数。
准确地确定σ0、σ1、uC、uS、σ2和vs这6个参数是利用LuGre摩擦模型进行摩擦补偿等相关研究的关键。
为研究在黏着阶段LuGre模型中鬃毛的动态特性不失一般性, 取质量为m的物体作为研究对象。该物体放置在固定的水平面上, 受外载荷P (t) 的作用而产生水平方向的位移x, 水平面对物体的摩擦力为Ff, 如图1所示。
物体与固定水平面之间的相对速度v=dx/dt。根据牛顿运动定律, 该物体的运动微分方程为
其中, 摩擦力Ff由式 (6) ~式 (8) 表示。式 (6) 和 (7) 是连续的非线性常微分方程。令:x1=x, , x3=z, 参数向量p= (σ0, σ1, σ2, vs, uC, uS) , 则式 (6) ~式 (9) 可以转化成形如式 (10) 的微分方程组:
当系统以摩擦力Ff作为输出时, 式 (10) 所表达的方程组在形式上与式 (3) 统一。因此, SIVIA方法可以应用到LuGre摩擦模型的辨识中。
考察上述系统中6个待辨识的参数组成的参数向量p。一方面, 随着待估参数数目的增加, 辨识所耗费的时间也会增多, 特别是在辨识过程中需要求解微分方程以获得仿真输出值的情况下, 这种耗时增多会更加明显。另一方面, 当上述系统处于稳定滑动状态时, 摩擦力是滑动速度的线性表达, 其中包含了4个静态参数σ2、vs、uC、uS。而稳定滑动状态下的速度-摩擦力曲线是实验可测的。因此, 利用传统的线性辨识方法如最小二乘法即可较好地辨识LuGre摩擦模型中的4个静态参数。这极大程度上降低了SIVIA法要辨识的参数向量的维数, 有利于在保证辨识精度的同时提高辨识效率, 对于LuGre摩擦模型在摩擦补偿中的广泛应用具有重要意义。
3 辨识实例
3.1 原始系统参数
本节以文献[4]中所提系统为例来验证本文所提方法的正确性和有效性。系统模型如图1所示。滑块质量m=1kg。系统摩擦参数如表1所示。
稳态滑动阶段所施加的外力为正弦函数P (t) =1.42sinωt。该系统对应的预滑动阶段摩擦力曲线如图2所示。
3.2 静态参数辨识
对于图1所描述的单自由度质量系统, 其稳态滑动状态下的摩擦力为[5]
由式 (17) 可以看出, 稳态下的摩擦力FSS是关于速度变量v的线性函数, 而函数的参数正是待辨识的4个LuGre摩擦模型参数。另一方面, 稳态下待估系统的速度-摩擦力曲线可以通过实验获得。因此, 对于该线性参数估计问题, 我们可以借助最小二乘法来解决。辨识目标函数为
式中, FSSi为离散的第i个摩擦力实测值;为第i个时间节点上由式 (17) 得到的摩擦力理论计算值。
辨识结果如表2所示。
3.3 动态参数辨识
在预滑动阶段, 摩擦力的影响因素主要由LuGre摩擦模型中的鬃毛刚度系数σ0和阻力系数σ1来表征。为此, 考察如下LuGre模型的线性化形式[5]:
为保证上述线性化模型适用, 实验中所施加的作用力非常小, 以保证产生的预滑动速度也较小。由式 (19) 可以得到x1=x3+c。c代表积分常数, 此处则表示初始内部状态变量z0。这里, 黏滞摩擦系数σ2已在3.2节中辨识得到, 此处直接作为已知参数代入而不再对其辨识, 且在预滑动阶段, 黏滞摩擦系数σ2相比鬃毛阻尼系数σ1要小很多, 对预滑动阶段的摩擦力影响也较小。
式 (19) 在形式上与式 (3) 统一, 因此, 可以利用SIVIA方法来辨识式 (19) 中的参数σ0、σ1、z0。其初始搜索区间[P0]如表3所示。误差精度ε=0.0005。表3给出了最终辨识结果。
为了对比本文所提方法与已有辨识方法的辨识效果, 本节同时利用文献[7]的方法对本例动态参数进行了辨识, 结果如表3所示。
由表3可以看出, 利用本文方法辨识所得两个动态参数的辨识结果比文献[7]中方法的结果更为接近原始数值, 表明所提辨识方法辨识精度更高。而且, 利用本文方法还额外获得了黏着阶段鬃毛变形情况, 这也有利于对LuGre模型动态特性进行更为深入的研究。
4 结束语
本文提出了一种基于最小二乘法和区间分析法的LuGre摩擦模型参数两步辨识方法。与传统非线性参数辨识方法相比, 本文所提方法的优势在于无需建立待估参数的显性目标函数, 也无需对实测数据进行复杂的转换, 克服了传统辨识方法中辨识结果对初值依赖性较强且不是全局最优解、辨识精度及收敛性难以保证等不足。本文给出的辨识实例, 表明了本文所提方法是正确和有效的。本文所提方法为LuGre摩擦模型在液压、气动、机器人等精密机械系统的摩擦补偿以及低速往复运动系统轨迹跟踪和动态预测等领域中更广泛的应用提供了可靠的保障。
摩擦模型 篇5
关键词：搅拌摩擦焊,数值模拟,温度场,热生成,热传导
0 引言
作为一种新型固相连接技术,搅拌摩擦焊(Friction stir welding,FSW)因其高质量、高效率、低成本、无污染、易操控等特点,自20世纪90年代被英国焊接研究所发明以来,在短短数十年内获得了长足的发展。相比TIG、MIG等传统熔焊技术,FSW过程中没有飞溅和烟尘,较为绿色环保,焊后也不存在合金元素烧损、气孔等缺陷。其焊缝组织均匀细密,接头性能优异,残余应力和变形小,焊接过程易于实现自动化,因此被广泛应用于航空、航天、汽车和船舶制造等工业领域[1]。
焊接过程中的温度场是搅拌摩擦焊的一个重要研究对象[2]。研究者通常借助于热电偶来测量焊接过程中的温度[3,4],然而,采用该方法存在一定的局限,如难以获知焊核区的温度分布、热电偶自身存在热惰性,容易导致温度-时间曲线偏移[5]等。相比之下,作为一种新兴技术,通过数值模拟方法研究FSW过程中的温度场,具有高效、便捷、直观等优点,近年来已被越来越多的研究者所采用。
随着对FSW物理机制的深入了解,描述其温度场分布的模型也越来越精确。本文简要概述了当前针对FSW温度场进行数值模拟时涉及到的主要模型以及它们的大致思路。从热量输入和输出的角度将这些模型划分为产热和传热模型。在产热模型中,分别介绍了与热输入有关的热源模型、产热量计算模型和接触模型;在传热模型中,则论述了与热输出有关的热传导、对流换热和热辐射。此外,考虑到材料的某些热学和力学性能会对温度场数值模拟的结果产生影响,也在文中对其进行了简单讨论。
1 产热模型
在研究FSW过程中的温度场分布时,产热模型是影响结果准确程度的关键因素。热源模型和产热量计算模型是产热模型的重要组成部分,前者是对产热机制的解释,后者是计算产热量的数学方法。在上述模型都已确定的情况下,搅拌头与工件之间的接触状态对产热量的计算结果有决定性影响,因此,接触模型也在文中进行了重点讨论。
1.1 热源模型
早期的研究假定FSW过程中的热输入完全来自于轴肩与工件表面的摩擦热[6,7,8,9,10,11]。采用该产热机制对铝合金或镁合金薄板进行模拟时可获得较为精确的结果,因为铝镁合金的导热系数大,轴肩处产生的热量可以迅速沿工件厚度方向传导至焊缝底部,同时薄板焊接时搅拌针与工件的摩擦面较小,产热量只占总体极少部分,可以忽略不计。但在厚板焊接时,搅拌针与工件的摩擦面增大,若仍然忽略搅拌针的影响,则总体产热量偏少,模拟结果与实际存在较大偏差。Christophe Gallais等[12]针对6000系铝合金,采用“试错法”探讨了搅拌针对FSW过程中热输入的影响。发现移除搅拌针后,得到的热输入量是有搅拌针时的82.6%,即搅拌针的热输入量可达总热输入的17.4%,可见搅拌针的影响不容忽视。值得注意的是,此处搅拌针的热输入既包括搅拌针表面与工件的摩擦热,也包括因搅拌针引起的塑性变形热,两者尚未区分开来。
随后的研究中,摩擦热源的模型越来越精细。徐韦锋等[13]对2219铝合金厚板FSW过程中搅拌头插入和焊接稳定阶段的温度场进行了模拟,在热源模型中同时考虑了轴肩以及搅拌针侧面与工件的摩擦热,并将仿真结果与实验数据进行对比,发现两者基本吻合。赵俊敏等[14]在模拟TC4钛合金板焊接过程中的温度场时,进一步考虑了搅拌针端面处的摩擦热。此时焊接过程中总摩擦热Qf为:
式中:Q1、Q2和Q3分别指代搅拌头轴肩、搅拌针侧面以及搅拌针端面与工件的摩擦热。
在搅拌摩擦焊的热源模型中仅考虑摩擦产热虽然有助于简化模型,但毕竟与实际情况不符。随着对搅拌摩擦焊产热机制的了解不断深入,越来越多的研究开始考虑材料的塑性变形产热对搅拌摩擦焊温度场的影响[15,16,17,18,19,20,21,22,23,24]。如图1所示,最终焊接过程中总的热输入Qtotal等于摩擦热Qf和塑性变形热Qt之和:
1.2 产热量计算模型
有多种方法可用于计算FSW过程中的产热量,根据其计算思路可大致分为直接求解法和间接求解法。直接求解法基于上述对热源的分析,直接计算摩擦热和塑性变形热,是目前应用较为广泛的产热量计算方法。采用直接求解法计算摩擦热时,通常分别计算搅拌头轴肩与工件的摩擦热Q1、搅拌针侧面与工件的摩擦热Q2以及搅拌针端面与工件的摩擦热Q3。针对不同形状的搅拌头,Q1、Q2和Q3的计算公式略有差异。为了简化分析过程,常假定轴肩底部为无内凹角的水平面,搅拌针为无锥角、无螺纹的规则圆柱体。此时轴肩与工件的摩擦面为一圆环,如图2所示。其内圆半径为搅拌针的半径Rp,外圆半径为轴肩的半径Rs。
在该圆环上任取一微小单元,若其与圆心的距离为r,搅拌头旋转速度为ω,单位面积所受摩擦力为τcontact,则该微元面积为rdθdr,所受摩擦力为τcontactrdθdr,相对位移为ωr,故生成的摩擦热功率为[7,25,26]:
对其在该圆环上积分,可得到轴肩与工件的摩擦热Q1:
同理,当搅拌针高度为H时,搅拌针侧面与工件的摩擦热Q2为:
搅拌针端面与工件的摩擦热Q3为:
τcontact随接触状态的变化而变化。在稳定焊接阶段,τcontact主要与材料的剪切强度有关。为简化模型,可假定各接触面上的τcontact均匀分布且大小相等[27,28]。此时总的摩擦热为:
式(3)-式(7)说明了采用直接求解法计算FSW过程中摩擦热的基本步骤。根据求解思路、搅拌头形状以及采用计量单位的不同,研究者得到的计算公式在形式上往往略有差异。曹朝霞[29]考虑了轴肩带有内凹锥角的搅拌头,如图3(a)所示,当内凹锥角为α时,轴肩与工件的摩擦热Q1为:
陈婷[30]讨论了不同搅拌针形状下的产热公式,认为当搅拌针为如图3(b)所示锥角为β的圆台体时,该搅拌针侧面的产热公式为:
据此,可知对于轴肩具有内凹角为α、搅拌针锥角为β的搅拌头,焊接过程中总的摩擦热为:
当采用直接求解法求解FSW过程中的塑性变形热时,其思路与计算摩擦热时一致,核心是计算应力和塑性应变的乘积。假定σ为应力,ε·为塑性应变率。则单位体积内的塑性变形热为[29]:
式中:ηp为塑性变形产热系数,即塑性变形功转化为热量的比例,通常认为它在0.8~0.99之间取值[31]。
对式(11)在体积V和时间增量 ΔT内进行积分,则 ΔT内产生的塑性变形热为:
根据式(2),结合上述摩擦热和塑性变形热的计算公式,便可得到FSW过程中的总产热量。
除了直接求解法外,还有许多研究者采用间接方法来计算FSW过程中的产热量。该类方法不直接计算摩擦热和塑性变形热,而是基于能量守恒原理,利用焊接时外力所做功来估算焊接过程中的热输入。如王希靖等[32]认为,可通过焊接时主轴电机的轻负载(摩擦头未与工件接触)功率P0和重负载功率P′0之差,计算出搅拌摩擦加热功率Pm。其公式为:
主轴电机的轻负载功率受多种因素影响,若采用数值方法来计算,则复杂且难以保证精确[32],因此较为适宜的方法是通过试验来测得。得到加热功率Pm后,将其乘以一定的工件吸热效率η,便可计算出FSW过程中工件吸收的热量[33]。
此外,Khandkar等[10]提出了一种基于扭矩的产热量计算模型:
式中:T为扭矩,W为搅拌头的旋转频率,其单位为r/min。
采用上述间接方法来计算FSW过程中的产热量时,通常模型较为简单,计算并不复杂。然而这类方法往往需要对焊接过程中能量的传递做理想化处理,不考虑能量在电机内部和焊机传动系统上的损耗[32]。此外,由于该方法不能直接获得搅拌头上各微小单元的产热量,因此还需额外假定一个摩擦热在轴肩上线性分布的斜率[33]。
1.3 接触模型
针对FSW过程进行模拟时,搅拌头表面与工件之间的接触关系是数值模型中最关键的部分[34]。根据不同的接触模型,接触剪切应力τcontact的求解方法各有不同。在经典的库仑接触模型中,假定搅拌头表面与工件之间为刚性接触,此时有:
式中:μ和p分别表示摩擦系数和触点压力。焊接过程中的摩擦系数μ为一变量,受载荷、滑动速度、温度等因素影响。周鹏展[35]通过分析2519铝合金摩擦状态与温度的关系,认为室温时摩擦系数最大,μ= 0.5。一旦摩擦界面的瞬时高温达到铝合金中某些低熔点物相的熔点,就会形成局部液相,对摩擦界面起到润滑作用,此时摩擦系数迅速降低。在前人研究的基础上,他假定当摩擦界面温度T = 650K时,摩擦系数为0.32。随着温度趋近2519铝合金的初熔温度818K,此时的摩擦类似于流体润滑下的摩擦,因此取818K下的摩擦系数为0.05。根据以上分析,他提出了摩擦系数μ(T)与实际摩擦温度T之间的非线性关系式:
Darvazi等[36]同样考虑了温度对μ的影响,认为当温度由25℃上升到1450℃后,μ由最大值0.78减小为0,由此建立了它们的线性方程式。
上述关系式都基于较多假设,且只考虑了摩擦系数与温度的关系。杜随更等采用试验方法,针对GH2132 管状试件,通过测量摩擦焊接初始阶段的摩擦扭矩、压力、温度及相对速度,提出了一个回归方程[37]:
式中:pf是摩擦表面压应力。Vf是摩擦表面速度,对于内径和外径分别为r和R的管状试件,其大小可简单计算为:
当考虑载荷、滑动速度、温度等因素对摩擦系数的影响时,更符合实际情况,但同时也会增加建模的工作量,并提高计算代价。事实上,一些研究为简化模型,也常将μ取为定值,如μ=0.3[38]或μ=0.2[39]等。
在采用经典的库仑接触模型时,随着FSW过程的进行,接触压力p不断增大,最终计算出的接触剪应力τcontact将远远高于材料的剪切强度。张昭等[38]的研究表明,在1400r/min的高转速下,采用经典的库仑接触模型将使得惯性效应变得非常明显,从而不能合理预测FSW过程。
Zhang等[40,41]采用修正的库仑接触模型来模拟焊接接触面的行为。在该模型中,接触剪切应力存在最大值τmax,其值等于工件材料的屈服剪切应力τyield。τcontact大小等于μp和τyield中的较小值,其求解公式为:
当假定μ为常数时,τcontact的变化规律如图4所示。由于τyield通常是温度、应变率和应力状态的函数。随着温度升高,τyield不断降低,产热量也随之减小,直至最高温度稳定在某一范围。通过限定τmax=τyield,该模型有效避免了模拟温度的无限升高。
张昭等[38]建立了完全热力耦合的搅拌摩擦焊接模型,分别采用经典的库仑接触模型和修正的库仑接触模型预测FSW过程中的温度场。对结果的比较显示,采用经典的库仑接触模型会使预测温度稍高。在转速较低时,两种模型的预测结果区别不大。但在高转速下,只有修正的库仑接触模型才能获得合理的计算结果。
Schmidt等[42]系统地分析了搅拌头表面与工件之间的接触状态,认为FSW过程中存在3类接触状态:滑移接触、粘性接触,以及部分滑移/粘性接触。在不同的接触状态下,应采用不同的方法求解τcontact。
1.3.1 滑移接触
当接触剪切应力τcontact小于工件的屈服剪切应力τyield时,搅拌头周围母材仅发生弹性切变,其运动速率vmaterial= 0,这种状态称为滑移接触状态。
1.3.2 粘性接触
若摩擦剪切应力大于工件的屈服剪切应力,搅拌头周围母材将发生塑性软化并粘接在搅拌头表面。此时,这部分母材在搅拌头的带动下沿其表面加速旋转,其运动速率vmaterial逐渐增大直至与搅拌头表面速率vtool相等,最终接触剪切应力τcontact与母材的屈服剪切应力τyield达到一致,这种状态称为粘性接触状态。
1.3.3 部分滑移/粘性接触
在实际的FSW过程中,搅拌头与工件的接触面主要为粘性和滑移共存的混合状态。在该接触状态下,搅拌头周围母材加速旋转,并以低于搅拌头表面速率vtool的运动速率达到稳定。
1.3.4 不同接触状态下的接触剪切应力和产热计算模型
计算FSW过程中的产热量时,需要确定搅拌头表面与工件之间的接触剪切应力τcontact。不同的接触状态下,τcontact的计算方法各不相同。假设搅拌头表面与工件材料的接触状态为滑移接触,此时可采用库仑法则求解τcontact:
将式(20)代入式(10),可知对于搅拌针半径为Rp、锥角为β、轴肩半径为Rs、内凹角为α的搅拌头,当转速为ω 时,滑移接触下的产热量为:
若搅拌头表面与工件材料的接触状态为粘性接触,此时接触剪切应力τcontact与母材的屈服剪切应力τyield一致,即:
因此,粘性接触下的产热量为:
在滑移接触状态下,系统的热输入来自于搅拌头表面和工件的摩擦热;在粘性接触状态下,热输入来自于工件材料的塑性变形热;而对于部分滑移/粘性接触状态,热输入同时受到摩擦和塑性变形的影响[27]。其产热量为:
式中:δ为接触状态变量,是工件与搅拌头接触点上的母材运动速率vmaterial和搅拌头表面的速率vtool之比,即:
式(24)可用于预测各类接触状态下的产热量。δ=0时,为滑移接触;δ=1时,为粘性接触;而当0<δ<1时,为部分滑移/粘性接触。Darvazi等[36]认为δ是温度的线性函数,随着温度由室温(25 ℃)上升到初熔温度(1450 ℃),δ也由0变为1。当模拟的温度达到最高值时,δ=0.7。一般来说,在FSW过程中,δ是0.37~0.8之间的数值[43]。
2 传热模型
FSW过程中的温度场不仅受到热输入的影响,同时也与热输出的机制息息相关。如前所述,热输入主要来源于搅拌头与工件接触面间的摩擦热,以及搅拌头周围母材的塑性变形热。热量产生后,迅速提高了搅拌头附近工件的温度,使该区域材料与周围产生温差。随后,热量从高温区域向低温区域传递,传递的方式主要为热传导和对流换热。
2.1 FSW过程中的热传导
热传导是指两个物体之间或同一物体的不同部分之间由于温度梯度而引起的内能交换。热传导遵循傅里叶定律,其形式为:
式中:q*为热流密度(W/m2),是单位时间内通过单位面积的热量;k为导热系数(W/(m·K));dT/dx为沿向的温度梯度,负号表示热量流向温度降低的地方。
在FSW过程中,搅拌头与工件之间、工件内部以及工件与垫板之间均存在热传导。
2.1.1 搅拌头与工件间的热传导
焊接开始后,搅拌头在机械力的作用下紧压在工件表面,并与其发生相对运动产生热量。通常假定产生的热量以一定比例流向搅拌头和工件,周明智等[44]将该比例取为0.5,即热量均等地分配给搅拌头和工件。冯天涛等[45]认为85%的热量流向了工件。可建立一个简单的热流模型来估算通过搅拌头损失的热量:在搅拌头上沿轴向取两点,分别测试其温度。同时,给出该热量的一个预估值,根据预估值模拟出搅拌头内的温度分布,并将实验数据与模拟结果相对照来调整预估值。通过多次试错,最终可估算出通过搅拌头损失的热量大约占总热量的5%[26,46]。
2.1.2 工件内部的热传导
热量首先产生于搅拌头附近,随后在温度梯度的作用下向低温区传递。传递的速度取决于工件材料的导热系数k。通常来说,导热系数随温度变化而变化。温度越高,导热系数也越高。然而相对于温度的变化,导热系数的起伏并不大,在一些研究中,为简化模型,也常对工件材料的导热系数取定值[47,48,49]。
2.1.3 工件与垫板间的热传导
已有多种方法用于模拟FSW过程中工件底部与垫板之间的热传导。根据对垫板建模与否,这些方法可被大致分为两类。若对垫板进行建模,则需要对照实验数据和模拟结果,估算出工件底面和垫板间的导热系数;若不对垫板进行建模,则需要在工件底面引入一个对流换热系数,以工件与底部空气的对流换热来模拟实际中工件底面和垫板的热传导。此外,根据对工件与垫板间接触条件的不同假设,这两类方法还可进一步细分如下:(1)对垫板进行建模,并假定垫板与工件底面各处均为完美的热接触,即垫板与工件相接触的表面温度相等;(2)对垫板进行建模,且假定仅在搅拌头正下方的垫板和工件才是完美的热接触;(3)不对垫板进行建模,并假定工件底部绝热,因而对流换热系数为0;(4)不对垫板进行建模,同时假定工件底部不绝热,需要根据实验数据和模拟结果的对照来估算导热系数。
Colegrove等[50]对垫板进行了建模,并指出工件底面与垫板间的导热系数约为1000 W/(m·K)。在另一项研究中,Colegrove等[51]假定搅拌头正下方的工件与垫板为完美的热接触,而除此之外的接触面则依旧为1000 W/(m·K)。
对垫板建模将会给数值模型增添额外的计算量,因此减小垫板的模型尺寸有助于降低计算成本。考虑到在搅拌头压力作用下,其正下方的垫板与工件接触最为紧密,而其他区域随着压力减小,垫板与工件间的裂缝逐渐增大。因此,一些研究者[10,52,53]将垫板简化为位于搅拌头行进路径正下方的垫条,如图5 所示。该垫条宽度通常与搅拌头直径相等,而高度不超过垫板厚度。
Ulysse[15]假定工件底部绝热,建立了一个不包含垫板的数值模型。尽管模拟预测出的温度总是略低于实验所测,但总体而言模拟结果与实验数据较为一致。Chao等[8]针对6061-T6铝合金建立了三维有限元模型,通过对比模拟预测的温度和实验测得的数据,认为该铝合金工件底面的对流换热系数为200 W/(m2·K)。Kandkar等[9]比较了4种不同的对流换热系数,发现对流换热系数值与模拟预测的最高温度呈负相关,当对流换热系数取1000W/(m2·K)时,模拟结果与实验数据最为接近。Liu H J等[54]建立了2219铝合金的有限元模型,采用大小为500W/(m2·K)的对流换热系数来估算通过垫板的热损耗。朱文峰等[55]针对A6N01-T5铝合金进行温度场模拟,选用的对流换热系数为253 W/(m2·K)。
2.2 FSW过程中的对流换热
对流换热是指固体表面与其周围流体间由于温差而引起的热量交换。对流换热实质上是热传导和热对流共同作用的结果,当流体流经固体表面时,通过固体表面的热流会在流体内部质点的作用下向流体扩散,即发生了固体表面与流体间的热传导。同时,流体内部质点的相对运动也会导致热量在流体中扩散,从而产生热对流。对流换热可用牛顿冷却公式来描述:
式中:hf为对流换热系数;TS和TB分别为固体表面及其周围流体的温度。
FSW过程中,工件表面与其周围空气之间发生的热量交换便是典型的对流换热过程。杜岩峰等[56]采用三维实体耦合的有限元方法来分析2219 铝合金搅拌摩擦焊热过程和温度场分布,提出铝板与环境的对流换热系数为1000 N/(s·mm·℃)。朱文峰等[55]模拟了A6N01-T5铝合金的温度场,对工件和环境采用的对流换热系数为162 W/(m2·K)。
2.3 FSW过程中的热辐射
热辐射是指物体发射电磁能,并被其他物体吸收转化为热的热量交换过程。物体温度越高,单位时间内辐射的热量越多。FSW过程中,因热辐射而损失的热量可通过简化的斯蒂芬-波尔兹曼公式来计算:
式中:ε为辐射率或吸收率;σ 是波尔兹曼常数,约为5.67×10-8W/(m2·K4);TS和Tα分别为物体表面和环境的绝对温度。
Nandan R等[57]综合考虑了FSW过程中工件表面与周围空气间的对流换热和热辐射,给出热量公式为:
相对于热传导和热对流,热辐射对FSW过程中温度场的影响较小,因此在许多文献中并未对其加以考虑[49]。
3 材料模型
除了产热和传热模型外,材料自身的一些热物理性能如比热容、导热系数等也会影响到温度场数值模拟的结果。对金属材料而言,比热容和导热系数通常都是温度的函数。为了提高温度场模拟的准确性,应尽可能精确地给定不同温度下的比热容和导热系数。然而,通过实验方法来分别测得各温度下的比热容和导热系数并不现实。在大多数研究中,研究者往往利用少数几个温度下的比热容和导热系数,采用插值法得到它们与温度的函数关系[35,45]。
此外,若采用热力耦合模型研究FSW过程中的温度场,则需综合考虑温度场与应力、应变场的相互影响,此时材料的某些力学性能同样不容忽视。由前文可知,FSW过程中的产热量与搅拌头和工件之间的接触剪切应力τcontact息息相关。在修正的库仑接触模型中,τcontact等于μp和τyield中的较小值。其中p是搅拌头表面与工件的触点压力。屈服剪切应力τyield是材料的固有属性,通常是温度、应变率和应力状态的函数。为了计算出尽量准确的产热量,一方面要根据不同的温度、应变率和应力状态给出τyield的不同数值,为此往往需要建立τyield与各状态变量的函数关系;另一方面,也要求同时模拟出准确的应力场,这除了需要建立合理的应力场数值模型外,同样需要给出尽可能符合实际的相关材料属性,即对材料模型的正确性有所要求。
总而言之,搅拌摩擦焊是一个包含了产热、传热、大塑形变形和塑性流动等众多物理现象的复杂加工过程,焊接时温度、应力、应变、应变率等状态均存在较大变化。在针对温度场进行模拟时,给出的材料模型越正确,各状态下的材料属性越接近实际材料,最终算出的温度场就越精确。
4 结语
数值模拟技术具有高效、便捷、直观等一系列优点,自其被用于研究FSW过程中的温度场以来,已涌现了大量的相关数值模型。随着对FSW物理机制的深入了解、模拟计算方法的不断更新,以及计算机计算能力的持续提高,这些模型也在逐步完善:
(1)早期的热源模拟仅考虑了搅拌头轴肩与工件的摩擦热,随后搅拌针与工件的摩擦热以及焊接过程中的塑性变形热也被纳入模型。
(2)针对产热量的计算有直接方法和间接方法两种,直接方法基于对热源的分析,直接计算摩擦热和塑性变形热,是目前较为主流的计算方法;间接方法基于能量守恒原理,通过测算外力所做功,如电机的输出功率来获得焊接热输入,原理和计算都较为简单,但往往需要对焊接过程中能量的传递作理想化处理,并人为指定摩擦热在轴肩上的分布斜率。
(3)FSW过程中存在3类接触状态:滑移接触、粘性接触,以及部分滑移/粘性接触。不同的接触状态对应着接触剪切应力τcontact的不同计算方法。相比经典的库仑接触模型,修正的库仑接触模型能获得更为合理的结果。
(4)FSW过程中,温度的传递方式主要有热传导和对流换热。通常假定摩擦产生的热量以一定比例流向搅拌头和工件。根据对工件与垫板之间接触条件的不同假设,两者间的热传导有不同的模拟方法。工件表面与周围空气间的热传递方式是对流换热,对流换热系数通过反复对照试验和模拟的结果来确定。通常不考虑FSW过程中的热辐射。
摩擦模型 篇6
冷轧轧制过程是一个典型的多变量、时变、强耦合和非线性过程, 多种因素相互影响最终作用在辊缝变形区域。高精度模型设定计算是稳定、高效轧制的前提和基础, 而轧制工艺数学模型又是高精度设定计算的核心[1]。轧制过程的复杂性决定了轧制工艺数学模型往往也具有很高的复杂性, 每个模型需要包含和体现多个因素对设定结果的影响[2], 如摩擦系数工艺数学模型就是一个包含轧制速度、轧辊粗糙度、轧制长度等变量的非线性多项式方程[3,4,5]。
非线性多项式形式的摩擦系数模型方程中的参数, 一般由常数表格查表直接给定, 但是对于不同轧制生产线, 或者相同轧制生产线处于不同的轧制状况时, 如果不进行优化修正, 往往不能满足设定计算的精度要求[6]。这就产生了如何根据实际生产状况确定摩擦系数模型参数的问题。
基于此, 笔者对典型非线性多项式摩擦系数模型进行分析, 并利用现场实际轧制数据, 直接对非线性多项式摩擦系数模型的参数进行回归优化, 避免了对复杂数学模型的线性化处理过程, 只要采集到的数据真实可靠, 分析优化的结果就是更能反映现场实际轧制情况的更加优化的摩擦系数模型参数。该优化方法于2013年开始应用于国内某冷连轧生产线后, 提高了该生产线摩擦系数模型设定精度, 进而提高了轧制力设定精度。
1 典型摩擦系数模型分析
式 (1) 所示为某生产线典型的非线性多项式摩擦系数模型, 包含实际轧制速度v, 轧辊表面粗糙度R和实际轧制长度L 3个自变量以及7个模型参数。
式中, u为摩擦系数;u0为基本摩擦系数参数;为速度变化影响参数;v0为参考轧制速度;CR为轧辊粗糙度影响参数;R0为轧辊表面参考粗糙度;cW为轧制长度影响参数;L0为轧辊轧制带钢的基准轧制长度。其中, 摩擦系数u是模型因变量, 需要根据模型自变量v, R和L以及模型参数计算得到。模型自变量中的实际轧制速度v可以通过现场速度测量仪表以200 ms的测量周期获得;轧辊表面粗糙度R和实际轧制长度L来自于生产线在线轧辊数据, 可以认为一个换辊周期内R不变, L随着轧制过程不断累计, 每一卷带钢对应一个轧制长度。u0, , v0, CR, R0, cW和L0是摩擦系数模型方程中反映轧制工况的参数, 对摩擦系数计算结果有直接影响, 一般通过查表获得。
优化摩擦系数模型的目的就是得到更能反映现场实际轧制工况的模型参数, 即在获得海量轧制力、压下量、张应力、轧制速度、轧辊累计轧制长度、轧辊粗糙度等过程数据的情况下, 先通过轧制理论公式反算得到实际摩擦系数, 然后通过优化算法得到模型方程中的参数。
针对上述非线性摩擦系数模型, 基于以下考虑, 本系统只将模型参数中的u0, 和cW作为优化对象: (1) v0, R0和L0是3个基本参数, 由生产线基本状况决定, 按照常数处理, 不进行优化; (2) 轧辊表面粗糙度数值基本稳定, 因此不对CR参数进行优化; (3) u0是基本摩擦系数, 该值的大小直接决定了最终的计算结果, 反映了轧制速度对摩擦系数的影响, cW反映了轧制长度对摩擦系数的影响, 因此对其进行优化。
2 模型优化方法
回归分析是最常用的参数优化方法, 通常使用的回归分析算法[3], 例如一元线性回归、多元线性回归、线性逐步回归算法都不适用于上述非线性多项式摩擦系数模型, 同时也无法通过变量变换的方法将其转化为多元线性回归。因此本文考虑将回归问题转化为多元非线性优化问题, 也就是寻找最优摩擦系数方程参数, 使得摩擦系数模型计算结果与实际摩擦系数最接近。
优化问题的核心是选择搜索方向和确定步长因子。梯度下降法是传统的优化问题解决方法, 其利用迭代点的负梯度方向就是函数值下降最快的方向这一特点, 将负梯度方向作为迭代的搜索方向。但是负梯度的特点决定了梯度法在远离极小点的时候逼近速度较快, 而接近极小点的时候逼近速度较慢 (只有线性的收敛速度) 。另一优化问题解决算法牛顿法将函数展开成Taylor级数, 利用函数负梯度和二阶导数矩阵构造搜索方向, 在靠近最优点附近的时候, 能够产生理想的搜索方向, 但是迭代发散问题是牛顿法的一个障碍[7]。Levenberg-Marquardt优化算法是梯度下降法和牛顿法的结合, 综合了上述两种方法各自的优点, 具有良好的迭代速度和收敛特性。它利用了二阶梯度的信息, 具有很快的收敛速度。当初始点远离最优点时, 负梯度方向是最快速下降的方向;当靠近最优点附近时, 在牛顿法迭代过程中引入步长因子和一维搜索, 保证迭代点的严格下降性, 这样就产生了一个理想的搜索方向。
Levenberg-Marquardt算法自身是一种寻找函数最小值的优化算法, 与其他优化算法一样, Levenberg-Marquardt算法寻找函数最小值的过程也是一个迭代过程:首先给定待回归参数的初始向量β, 然后在每次迭代过程中计算β的步长因子增量向量δ, 获得新的待回归参数向量β=β+δ, 计算模型方程的残差平方和Q, 不断迭代, 最终获得使Q最小的待回归向量β。LevenbergMarquardt算法用于计算δ的模型如下:
式中, J为包含线性化拟合矩阵和极小化函数梯度向量信息的雅可比矩阵;λ为阻尼因子;I为搜索方向矩阵;y为目标值矩阵;f为当前β对应的模型计算值矩阵。
将Levenberg-Marquardt算法应用于轧制过程摩擦系数模型参数的优化时, 必须结合轧制过程工艺特点和要求。具体步骤如下:
(1) 根据工艺数据相关性分析确定对摩擦系数模型设定精度起关键作用的数据对象, 构造待分析的数据结构, 结构中包括轧制速度、轧制长度、轧辊表面粗糙度和反算得到的摩擦系数;
(2) 给定待回归参数u0, 和cW的初始值;, duv
(3) 利用摩擦系数模型, 计算在当前参数u0, 和cW情况下的各组轧制速度、轧制长度和轧辊表面粗糙度对应的设定摩擦系数;
(4) 计算设定摩擦系数和样本中反算的摩擦系数之间的残差平方和, 执行第9步;
(5) 构造参数u0, 和cW的优化迭代的阻尼因子λ和搜索方向I;
(6) 实例化待优化多项式非线性方程, 并进行数值求导, 构造线性化拟合矩阵;
(7) 求解线性化拟合矩阵和极小化函数梯度向量, 获得待优化参数u0, 和cW的步长因子增量向量;
(8) 计算新的待优化参数u0, 和cW, 执行第3步;
(9) 判断残差平方和是否满足允许的最小偏差0.001, 不满足则重复第5~8步;满足则退出, 将当前使用的优化参数u0, 和cW作为优化结果。
Levenberg-Marquardt优化算法流程如图1所示。
3 实际应用效果
基于上述摩擦系数模型优化方法和流程, 在Visual Studio 2010环境下, 采用C++语言实现了Levenberg-Marquardt优化算法, 以及与之相关的工艺数据读取和处理、工艺数据相关性分析、优化结果存储和应用等功能, 并将其应用于国内某冷连轧生产线模型系统, 对摩擦系数模型参数进行优化, 验证其优化效果。
在轧制过程中, 摩擦系数是一个无法精确测量的物理量, 但是摩擦系数模型的计算精度直接影响最终轧制力模型的设定精度。通过对累计10个钢种的6 715卷带钢 (对比钢卷2 808卷) 的轧制力设定精度的变化情况, 验证基于Levenberg-Marquardt算法的摩擦系数模型参数的优化结果, 如表1所示。从表1可以看出, 在现场多个钢种的实际轧制过程中, 应用本文所述的优化方法, 轧制力模型精度都得到了较大的提高, 最高可以改善3.42%。
图2显示了在同钢种、同规格情况下, 钢种M3 A33在优化前后的轧制力精度 (偏差百分比) 变化情况, 可以看出优化后的轧制力精度要明显高于优化前的轧制力精度, 并且更加稳定。
4 结论
面对工况不断变化的轧制过程, 本文根据轧制过程摩擦系数模型的特点, 在分析轧制工艺的基础上, 充分利用海量的实际工艺数据, 通过Levenberg-Marquardt优化算法对非线性多项式摩擦系数模型参数进行回归分析, 避免了对非线性多项式模型回归分析过程中为满足线性化要求而删除重要变量的缺点, 获得的摩擦系数模型参数具有更高的可靠性和可信性, 从而提高了轧制过程摩擦系数模型设定精度, 进而提高了轧制力模型的设定精度, 最终实现高精度轧制。
参考文献
[1]Olof Wiklund.Process control of cold rolling mills[J].Scandinavian Journal of metallurgy, 1998, 27 (2) :82-87.
[2]郭立伟, 高雷, 陈丹.数据挖掘技术在冷连轧机板形控制系统中的应用研究[J].冶金自动化, 2012, 36 (2) :30-33
[3]杨节.轧制过程数学模型[M].北京:冶金工业出版社, 2002, 36-42.
[4]高雷, 郭立伟, 陈丹.冷连轧轧制力模型中变形抗力和摩擦系数的分析[J].轧钢, 2013, 30 (4) :12-15.
[5]Venkata Reddy N, Suryanarayana G.A set-up model for tandem cold rolling mills[J].Journal of Materials Processing Technology, 2001, 116 (2) :269-277.
[6]白振华, 周庆田, 窦爱民, 等.冷连轧高速轧制过程中摩擦因数机理模型的研究[J].钢铁, 2007, 42 (5) :47-50.
摩擦模型 篇7
三、“静静”突变模型;四、“动动”突变模型.
一、“静动”突变模型
“静动”突变是指物体在摩擦力和其他力作用下处于静止状态,当其他力发生变化时,如果物体不能保持静止状态,则物体受到的静摩擦力将“突变”成滑动摩擦力.
例1如图1所示,在水平桌面上放一木块,用从零开始逐渐增大的水平拉力拉木块直到沿桌面运动,在此过程中,木块所受到的摩擦力f的大小随拉力F大小变化的图象正确的是()
解析:1.开始时,当F=0时,桌面对木块没有摩擦力,f=0.当木块受到水平拉力F较小时,木块静止不动,但有向右运动的趋势,产生静摩擦力f静,根据水平方向上二力平衡可知:f静=F,体现在图像上,直线倾角α=45°.
图2
2. 随着水平拉力F不断增大,木块向右运动的趋势也逐渐增强,,桌面对木块的静摩擦力也相应的增大,直到水平拉力F足够大时,木块开始滑动,此时桌面对木块的静摩擦力达到最大值fm.
3. 当静摩擦力f达到最大静摩擦力之后,木块开始滑动,木块运动状态将发生了变化,此时静摩擦力“突变”为滑动摩擦力,F动=μFN=μmg保持不变,且F动<fm.
基于上述分析推证:(B)正确.
二、“动静”突变模型
“动静”突变是指在摩擦力和其他力作用下,做减速运动的物体突然停止滑行时,物体将不受摩擦力作用,或所受的滑动摩擦力“突变”成静摩擦力.
例2把一重为G的物体,用一个水平的推力F=kt(k为恒量,t为时间)压在竖直的足够高的平整墙上,如图3所示,从t=0开始物体所受的摩擦力随t的变化关系是图4中的()
图4
解析:(1)由于物体受的水平推力为F=kt,从t=0开始水平推力F=0,物体所受滑动摩擦力为0.当F较小时,物体所受的摩擦力Ff小于物体的重力G,物体将沿墙壁下滑,此时物体与墙壁间的摩擦力为滑动摩擦力.
(2)当压力F不断增大,物体受到滑动摩擦力Ff不断增大.当滑动摩擦力Ff大小等于重力G时,由于惯性作用,物体不能立刻停止运动,物体受到的摩擦力仍为滑动摩擦力.(3)随着时间t的增加,推力F继续增大,滑动摩擦力Ff>G时,物体做减速运动直至停止运动,此时滑动摩擦力“突变”为静摩擦力,静摩擦力与正压力大小无关,物体处于静止状态,即使推力再增大,也不会影响物体的静摩擦力大小,此时静摩擦力的大小等于重力G.基于上述分析推证:(B)正确.
三、“静静”突变模型
“静静”突变是指物体在摩擦力和其他力的作用下处于静止状态,当作用在物体上的其他力的合力发生变化时,如果物体仍然保持静止状态,则物体受到的静摩擦力的大小或方向将发生“突变”.
例3一木块放在水平桌面上,在水平方向共受到三个力即F1、F2和摩擦力的作用,木块处于静止状态,如图5所示,其中F1=10 N,F2=2 N,若撤去F1,则木块受到的摩擦力为()
(A)10 N方向向左
(B)6 N,方向向右
(C)2N,方向向右
(D)零
解析:(1)当物体受F1、F2及摩擦力三个力的作用而处于平衡状态时,由平衡条件可知物体所受的摩擦的作用大小为8 N,方向与F2的方向相同,水平向左,由此可知最大静摩擦力fmax≥8 N.
2. 当撤去力F1后,F2=2 N<fmax,物体仍处于静止状态,由平衡条件可知物体所受的静摩擦力大小和方向均发生“突变”,且与作用在物体上的F2等大反向.基于上述分析推证:(C)正确.
四、“动动”突变模型
“动动”突变是指在摩擦力和其他力作用下处于运动状态,而滑动摩擦力存在于发生相对运动的物体之间,当两物体的速度达到相同时,物体所受的滑动摩擦力的大小或方向将发生“突变”.
例4传送带以恒定的速率v=10 m/s逆时针方向运动,已知它与水平面成θ=37°,如图6所示,PQ=16 m,将一个小物体无初速度地放在P点,小物体与传送带间的动摩擦因数为μ=0.5,问当传送带逆时针转动时,小物体运动到Q点的时间为多少?(g=10 m/s2,sin37°=0.6,cos37°=0.8)
基于上述分析推证:小物体从P点运动到Q点的时间t=t1+t2=2 s.
参考文献
[1]叶秀清.摩擦力的突变问题归类分析[J].物理教师,2014,35(3):82-84.
摩擦模型 篇8
光电稳定平台的扰动因素主要来源于载体对平台的扰动[3]。该扰动通过轴系间的摩擦力逐个框架耦合至传感器视轴[[4,5],间接的带动视轴一同运动,从而影响平台的视轴稳定精度。由此可知轴系间摩擦力是影响其视轴稳定精度的一个重要因素。因此应用合适的摩擦模型从最大程度上抑制摩擦对平台的影响,已经成为提高光电稳定平台视轴精度的关键问题。
目前工程上应用最多的摩擦模型为Stribeck模型,但该模型在速度为零时摩擦力有突变,没有完全拟合实际的摩擦曲线[[6,7,8]]。而Lu Gre模型能够很好地解决这些问题。该模型引入了刚鬓模型,采用6个参数能够精确地描述了摩擦力的静态和动态特性,是一个比较完善的摩擦模型[[9,10]。
文中针对摩擦影响光电稳定平台视轴精度问题,通过对某型两轴四框架航空光电稳定平台内方位速度环进行分析,建立了基于Lu Gre模型的摩擦补偿模型,并对与摩擦力矩相等效的控制电压的进行分析。通过实验数据和遗传算法对模型的摩擦参数进行辨识。通过所辨识的结果对系统的摩擦力矩进行补偿并进行实验验证。结果表明采用Lu Gre摩擦模型补偿后能够更好地提高平台的抗干扰性能,平台的扰动隔离度明显提升。为将该摩擦补偿方法应用到实际工程项目中提供了基础。
1 航空光电稳定平台系统建模
文中所研究对象为某型两轴四框架航空光电稳定平台。由于方位和俯仰轴相互垂直,因此两者之间的相互影响很小。文中以内方位速度环为例对摩擦力补偿方法进行讨论。在考虑电流环但忽略电感的作用下,得到平台系统的速度环开环模型结构图如图1所示。
由图1得含有摩擦力矩的系统动力学方程为
式(1)中,Tf为摩擦力的扰动力矩;J为折合到电机上的总的转动惯量;w为负载相对于平台的角速度;u为控制电压;k伺服放大器和直流电机等效模型的比例系数。
经实验测试,内方位速度环的传递函数为
式(2)中,f(s)代表与系统中摩擦力矩相等效的控制电压。在对摩擦力进行补偿时只需考虑对与摩擦力相等效的控制电压进行补偿即可,因为在实际应用中,由于被控系统的复杂性和保密性等因素常常限制了对实际被控系统内部特性的深入了解,这种方法能够有效地回避平台参数未知的问题。
根据公式(2),对图1所示的速度开环模型进行简化,如图2所示。由此可知,只要能够准确辨识出与摩擦力矩相等效的控制电压f(s)的模型,对系统中摩擦力矩进行前馈补偿,便能够达到抑制摩擦力的目的。
2 摩擦力矩模型
Lugre模型是1995年由C.Canudas De Wit等通过大量实验提出的[11,12]。该模型融合了鬃毛模型的思想,并在此基础上做了进一步的完善[[13,14]。在实际运行情况下平台其他因素的影响会导致摩擦特性随着运动方向和所处位置的不同而变化,其表达式如下。
式中,θ为负载相对平台的角位移;w为负载相对平台的角速度;σ0+、σ0-为正反向鬃毛刚性系数;σ1+、σ1-为正反向阻尼系数;σ2+、σ2-为正反向黏性系数;z代表鬃毛的平均形变量;Ts+、Ts-为正反向最大静摩擦力矩;Tc+、Tc-为正反向库伦摩擦力矩;ws+、ws-为正反向切换速度。对平台摩擦力矩进行补偿需要对摩擦力模型的参数进行准确的辨识。
3 平台模型参数辨识
上述模型的参数辨识过程分为两步。第一步为静态参数辨识,静态参数主要是用来描述系统处于稳态时的摩擦特性。其辨识方法主要是基于实验和数据拟合对静态模型参数进行辨识;第二步为动态参数辨识,动态参数是用来描述系统处于临界状态时摩擦特性。其辨识方法主要是在静态模型的基础上,通过实验和遗传算法对动态模型参数辨识。
3.1 静态参数辨识
当系统处于稳态时,摩擦力矩和速度之间的关系为
与上述摩擦力矩相等效的控制电压为
通过等效电压的方法把对摩擦力模型中较为难辨识的参数转换为对公式(7)中参数A+,A-,B+,B-,C+,C-,ws+,ws-进行辨识,这样转换可以使参数辨识较为容易。
为了对静态参数进行辨识,构造速度闭环系统,使系统以不同速率做恒速运动,且在低速运动时使速度以小间隔变化,从而得到摩擦力的负斜率特性。图3为测得速度和对应控制电压关系图。由于系统做匀速运动,所以对应的控制电压刚好为系统克服摩擦力矩相对应的电压。
图3中纵坐标u为电压所对应的码值,利用Matlab对该曲线进行分段拟合,辨识得到的参数结果如表1所示。
3.2 动态参数辨识
在静态模型参数已知的基础上,对平台的动态参数进行辨识。由公式(2)~式(7)和静态模型的基础上可以得到与平台摩擦力矩相等效的控制电压模型。
为了准确得到系统的补偿模型需要对公式(8)中的参数a+,a-,b+,b-进行辨识。给系统速度环一个正弦速度,采集平台运动速度记为w0(t),并通过加入带有摩擦的系统模型得到含有参数a,b且速度记为w[(a,b),t],以二者的方差作为评价函数。
利用遗传算法,并根据实际经验对变量范围进行限制,对动态模型参数进行辨识,参数辨识结果如表2所示。
为了验证所辨识出模型的正确性,在上述辨识模型的基础上通过仿真实验来模拟自由减速运动,在电机开路的情况下给定一个正向初速度,让光电稳定平台在只有摩擦力的作用下做自由减速运动,由于摩擦的滞后特性使得仅在摩擦力作用下引起速度过零振荡,图4为辨识模型所产生的自由减速曲线与原始曲线的比较,从图中可知辨识模型所获得的自由减速曲线与原始曲线基本重合。
4 补偿试验与结果分析
4.1 平台摩擦力矩的补偿方案
由于摩擦力矩对系统的影响可以等效为一个控制电压对系统的影响,因此只需在控制电压的输入端加入一个与摩擦力矩相等效的大小相等方向相反控制电压。其补偿基本原理为以系统的实际输出角速度为输入量,利用Lu Gre模型计算出平台每个时刻的摩擦力矩相对应的控制电压,在速度环控制器输出端叠加该电压即可。施加摩擦补偿的系统结构图如图5所示。
4.2 平台低速性能试验
载体处于静止状态,让平台以1 Hz,1°/s的正弦运动,对加入摩擦力矩补偿前和补偿后的系统闭环速度跟踪性能进行对比试验,采集陀螺在两种不同控制方法下的速度,实验结果如图6所示。比较可知,未加入摩擦力矩补偿前在速度过零时平台会出现爬行现象,而加入摩擦力矩补偿后,由于对摩擦给系统的影响进行了实时补偿,低速性能有了显著提高。
4.3 视轴稳定精度对比试验
为了验证加入摩擦力矩补偿后载体对平台速度扰动的抑制能力,将光电稳定平台安装在能够模拟载体运动的飞行模拟转台上并进行实验,图7为实验时用到模拟载体运动的摇摆台和光电稳定平台。
平台给定速度为零,飞行模拟转台以2.5 Hz以内的任意频率做正弦运动的情况下,通过采集陀螺信号来观察验证加入摩擦补偿后平台速度环对扰动的抑制能力,表3给出了模拟飞行转台以幅值为1°、频率在0.1~2.5 Hz之间每隔0.5 Hz的步进频率进行摇摆时,验证平台加入摩擦补偿后对各个频率扰动隔离度的提高程度。由表3可知加入摩擦力补偿后系统速度环对各个频率的扰动抑制能力明显增强,最优情况达到16.90 d B。
为了更直观比较加入摩擦补偿前后光电稳定平台速度对扰动的抑制能力,图8给出了当平台期望转速为零时,模拟飞行转台以幅值为1°、频率为2.5Hz做正弦运动时加入摩擦补偿前后速度稳定情况。图9是对图8数据进行频谱分析。
从图8可知,采用Lu Gre模型对摩擦力进行补偿,在2.5 Hz的扰动下平台速度的最大幅值明显减小,加入摩擦补偿后其最大峰值为0.2°/s与补偿前最大峰值为1.0°/s相比,减小了80%。由此可知采用Lu Gre模型对摩擦力进行补偿后,平台的抗干扰能力明显提高。
由图9可得,加入摩擦力补偿后扰动残余量明显小于未加摩擦补偿时扰动残余量,在频率2.5 Hz下加入摩擦力补偿后扰动残余量大约是未补偿时的0.142 8倍,即其扰动隔离度提高了16.90 d B。而且针对由其他因素引起其他频率干扰也有一定抑制作用。
5 结论
(1)利用Lu Gre模型准确建立了航空光电稳定平台摩擦模型,确定与该平台摩擦力矩相等效的控制电压,通过补偿该等效电压实现对摩擦力的补偿。
(2)通过实验和遗传算法对平台Lu Gre摩擦模型进行参数辨识,实现了Lu Gre模型在航空光电稳定平台上的应用。
(3)加入摩擦补偿后,平台在低速运行时爬行现象明显减弱。