吴　石，边立健，刘献礼，宋盛罡，姜彦翠

(哈尔滨理工大学 机械动力工程学院，哈尔滨　150080)



薄板件铣削颤振稳定性的非线性判据实验研究

吴石，边立健，刘献礼，宋盛罡，姜彦翠

(哈尔滨理工大学 机械动力工程学院，哈尔滨150080)

铣削过程中非线性动力学行为一直伴随整个切削过程，为准确地判定和预测加工过程的颤振稳定性，基于实验方法，研究了两端固定薄板件铣削颤振稳定性的非线性判据。实验中以薄板件振动信号为研究对象，基于相平面法、庞加莱法和频谱分析了不同加工参数时的振动信号，绘制并讨论了最大Lyapunov指数与主轴转速和铣削深度的变化关系。最后以最大Lyapunov指数作为判据，通过等高线法确定铣削颤振稳定域，并和基于全离散法得出的铣削颤振稳定域进行比较分析，实验得出了航空铝合金7075-T6薄板件颤振稳定域的非线性判据。

薄板件；铣削颤振；非线性判据；李雅谱诺夫指数

现代航空航天对零件的重量、强度的要求逐渐提高，使薄壁结构件被广泛应用于飞机的大梁、壁板等众多主承力结构件上，有的成品零件壁厚仅仅1 mm。薄壁件虽然具有比强度高、相对重量较轻等众多优点，但是薄壁件的加工却存在诸多问题，其中切削加工过程中的振动失稳(如颤振)就严重制约着薄壁工件的加工质量，为实现薄板件铣削加工参数的优化，确保铣削过程的稳定性，必须采用方便有效的方法对加工过程的颤振稳定性进行准确的判定、预测。

在实际加工中主轴-刀具-工件，以及夹具系统是复杂的非线性系统[1]，传统线性理论对于切削厚度较小时的非线性系统采用线性理论近似方法，可以满足低速要求，但是随着高速铣削的发展，采用传统线性理论已经不能保证所建模型的精度，难以准确预测颤振临界值和工件的表面位置误差、粗糙度等。于是GRADISEK等[2-3]分析了不同切削深度对加工中非线性振动的影响，并基于分岔模型进行论证，提出颤振发生时振动信号中具有低维混沌振动现象。DAVID等[4]基于非线性动态切削力建立了二阶时滞的再生型颤振模型，分析了系统振动由Hopf分岔进入混沌振动状态，得到了Hopf分岔稳定性边界条件以及表征值。SZALAI等[5]将刀具对工件切削简化为碰撞行为，建立了高速铣削的非线性动力学模型，认为高速铣削时存在混沌振动，并对稳定性边界进行了亚临界分岔分析。STEFANSKI等[6-7]针对非线性铣削振动指出，动态系统的李雅普诺夫指数(Lyapunov Exponent，LE)是分析混沌运动的一种有效方法。这是由于动态切削力的非光滑特性等导致动力学模型的相空间是多维的，虽然计算两个轨迹间距离及相邻轨迹间发散度是比较困难的，但是可以采用最大Lyapunov指数来进行分析铣削系统的振动信号。

李忠群等[8]通过数值方法求解非线性铣削动力学方程，并根据稳定性判据得到了时域颤振稳定域。孔繁森等[9]基于Lyapunov指数和Kolmogorov熵分析了切削系统从无颤振状态到颤振状态这一转变过程中的非线性行为的变化特征。李西彬等[10]研究中发现加工过程振动信号的非线性特征明显，在相同工艺参数条件下却不一定获得完全相同的加工效果。同时，在铣削过程实验中，发现薄壁工件振动信号的非线性特征比基于理论模型的振动信号要更为明显和复杂[11]，目前薄壁工件的多场耦合振动的分岔与混沌运动的研究，也取得了一些成果，YEH等[12-14]对热弹耦合矩形板振动的分岔与混沌进行了研究。

虽然非线性动力学模型在颤振预测方面占有巨大优势，但是如何基于非线性理论判定薄壁件切削的稳定性，并在实验的基础上进行准确的颤振稳定域预测，建立一个针对具体工件的铣削颤振稳定性判据仍是个挑战。本文将以航空铝合金7075-T6两端固定的薄板件为研究对象，探讨以工件振动信号的最大Lyapunov指数作为阈值，通过实验判定铣削加工过程颤振稳定域，确定两端固定的薄板件加工的颤振稳定域的非线性判据，并为从理论上研究热、力耦合作用下的薄板件切削振动的非线性稳定性奠定基础。

1　颤振过程的分析方法

由于铣削系统的非线性(动态铣削力的非光滑性、系统的非线性阻尼和刚度等)，使切削振动具有明显的非线性特性，而最大李雅谱诺夫(Maximum Lyapunov Exponent)指数、庞加莱(Poincare)映射及相平面图法等分析方法，是有效地分析铣削颤振非线性动力学行为的有效手段。

1.1最大Lyapunov指数

最大Lyapunov指数是描述初始时刻两个无限靠近的点随时间演化而分离的特征，作为混沌运动的特征参量，表示相轨迹的最大发散程度，或对初值的最大敏感程度[15]。

设变切深铣削振动信号x1，x2，… xN(单变量时间序列)，式中N是时间序列的总数。根据PACKED等[16]提出的时间延迟思想，可以重构出所观察到的动力学系统的相空间，基于这个思想，对时间序列进行相空间重构，得到重构的轨迹X，表示为

X=[X1,X2,…,XM]T

(1)

式中，M为相空间重构后轨迹点的个数，Xi为铣削振动系统在间断时间点i的状态，它可以表示为

Xi=[xi,xi+τ,…,xi+(m-1)·τ]

(2)

式中，τ为时间延迟，m为嵌入维数，其中M=N-(m-1)τ。对重构相空间X分为n段：[X1,X2,…,XT]，[XT+1,XT+2,…X2T]，…，[X(n-1)T+1,X(n-1)T+2,…,XnT]。其中每段的长度T=M/n，称为演化时间。

(3)

式中:Δt是采样时间间隔。

在不同的相空间维数下，对铣削时薄壁件的振动信号进行分析计算，进而研究不同加工参数下的铣削颤振稳定性。

1.2Poincare映射和相平面图

Poincare映射的分析方法是横坐标为系统响应的位移值，纵坐标为系统响应的速度值，表现为一个运动点，其计算方式为间隔每一个时间段进行数据提取。若铣削系统为周期振动，则对应于Poincare图上的一个孤立点，若振动信号的周期为N时，在图中出现N个独立点，且周期数同孤立点数相同。如果振动信号变化剧烈发生颤振时，则在图中显示为离散点堆积图[16]。

相平面图是振动系统的解轨迹投影在相空间的运动曲线，对于确定性的运动系统，有下列几种形式：振动信号周期为1或N时，对应相图中出现1条或N条自封闭的曲线，振动信号具有强颤振特性时封闭的曲线会逐渐拓展为杂乱的曲线群。

2　薄板件铣削加工实验

实验在VDL-1000E型机床上进行，工件模型为两端固定的薄板件，该件长120 mm，宽80 mm，壁厚为5 mm，实验装置与传感器布置如图1所示。测力仪型号为KISTLER9257B，振动加速度传感器为PCB356A25，置于加工面的背侧，其灵敏度为10.42 mv/g。刀具为山特维克(R216.64-08030-AO09G 1610)直径8 mm、螺旋角30°的四刃硬质合金球头铣刀。数采分析系统为东华DH5922的信号采集系统，采样频率为5 kHz。加工参数：顺铣，径向铣削宽度ae=0.2 mm，进给速度fz=0.1 mm/tooth，铣削实验中工件与刀具的相对倾角设为30°，走刀路径为沿水平方向来回进给。

表1　薄板件模态参数

假设刀具为刚性、工件为柔性，工件前二阶模态参数如表1所示。首先基于刀具-工件模型双自由度系统，根据全离散法[17]获得的两端固定薄板件的铣削颤振稳定域如图3所示，根据图中所标示的铣削加工参数点进行铣削加工实验，即每一个参数点对应一定的切削深度和主轴转速，这样总共进行121个加工参数的切削振动测试。如图2所示的两端固定的薄壁样件一共加工了二个，针对第一个样件铣削实验时将每一转速作为一组实验，一共进行13组实验，转速从3 000 r/min增加到6 000 r/min。针对第二个样件铣削实验时将每一转速作为一组实验，也一共进行13组实验，转速从6 000 r/min减少到3 000 r/min。

加工区域在工件上如图2所示。每组包括10个或9个不同的切深。同组实验过程中忽略上一加工参数实验后由于工件材料去除对下一参数实验所产生的影响。同时进行下一组实验前，先对工件进行加工修模，使每组实验工件薄壁部分的厚度及高度等尺寸相同，这样可以尽量确保各组实验工件的模态特性相对一致性。这两个试验样件在相同参数切削时，薄壁件振动情况基本相同，总体误差小于9%。然后针对每一种工况，基于相平面法、庞加莱法和频谱分析了其实测的工件振动信号，判定在此加工参数下是否发生了颤振。并计算该振动信号的最大Lyapunov指数。

图1　薄板件铣削加工实验Fig.1　The milling experiments of thin wall

图2　两端固定薄板件的13组加工区域分布Fig.2　The13 groups at both ends of the fixed thin wall processing regional distribution

同时，在实验中发现，在薄壁件加工中瞬时铣削力信号中的混沌(不确定性)特性，没有对振动信号进行混沌(不确定性)分析明显，故薄板件铣削加工实验以薄板件的振动信号为研究对象。

图4分别是两端固定薄板件在不同加工参数时(A(3 500 r/min，0.3 mm)、B(3 500 r/min，0.5 mm)、C(4 500 r/min，0.3 mm)及D(5 000 r/min，0.7 mm))的工件加速度振动信号以及其相平面图、Poincare截面图及频谱分析结果。

图3　两端固定薄板件的铣削颤振稳定域(○为铣削加工参数点)Fig.3　Milling chatter stability region of fixed at both ends of the thin wall(○ is the parameter point of milling)

由图4可知，当以加工参数A(3 500 r/min，0.3 mm)进行铣削实验时，其工件振动的加速度时域信号幅值为5 m/s2左右，研究中取振动信号的4 096个采样点，通过滤波、离散傅里叶变换后发现能量在刀齿切削频率(3 500×4/60=233 Hz)及谐波频率(467 Hz)附近，且能量分布较均匀，如图4(d)所示。从振动信号的相平面图4(b)可知，工件的振动响应逐渐收敛为围绕中心点的封闭曲线，而Poincare截面图4(c)中的点也较集中且仅有数个点，此加工参数下振动信号的最大Lyapunov指数(嵌入维数等于2)虽然为0.277 3，但经测试其表面粗糙度、以及相平面、Poincare截面图、频谱分析后，综合判定该加工参数时两端固定的薄板件铣削处于无颤振的稳定切削状态。

当以加工参数B(3 500 r/min，0.5 mm)进行铣削实验时，工件振动加速度时域信号如图4(e)所示，振动幅值达到11.2 m/s2左右。从振动信号的相平面图4(f)可知，工件的振动响应逐渐发散为多条曲线且没有规律，而Poincare截面图4(g)中的点也变得分布不均，此加工参数下振动信号的最大Lyapunov指数(嵌入维数等于2)为0.613 1。经滤波、傅里叶变换后得到的频谱如图4(h)所示，发现峰值不仅仅处于刀齿切削频率(3 500×4/60=233 Hz)及谐波频率(467 Hz)附近，而且能量有向刀具结构频率(1 278 Hz)处积聚的趋势，经测试其表面粗糙度、以及相平面、Poincare截面图、频谱分析后，综合判定该加工参数时两端固定的薄板件铣削处于颤振切削状态。

(e) 振动信号(3500r/min,0.5mm)(f) 相平面图(3500r/min,0.5mm)(g) Poincare截面图(3500r/min,0.5mm)(h) 频谱分析(3500r/min,0.5mm)

(i) 振动信号(4500r/min,0.3mm)(j) 相平面图(4500r/min,0.3mm)(k) Poincare截面图(4500r/min,0.3mm)(l) 频谱分析(4500r/min,0.3mm)

(m) 振动信号(5000r/min,0.7mm)(n) 相平面图(5000r/min,0.7mm)(o) Poincare截面图(5000r/min,0.7mm)(p) 频谱分析(5000r/min,0.3mm)

图4不同加工参数下工件振动信号及其相平面图,Poincare截面图,频谱分析
Fig.4 Workpiece vibration signal and its phase plane，Poincare-sectional view， spectrum analysis under different processing parameters

当以加工参数C(4 500 r/min，0.3 mm)进行实验时，振动加速度时域信号如图4(i)所示，工件出现了明显振动，最大振幅甚至达到13.3 m/s2左右。其振动信号的相平面图4(j)杂乱无章，Poincare截面图4(k)中的点也变得发散，且此参数下振动信号的最大Lyapunov指数(嵌入维数等于2)为0.611。将振动加速度信号通过滤波、傅里叶变换后发现能量向刀具结构频率(1 278 Hz)处积聚，并且能量增大非常多，如图4(l)所示，符合非稳定切削颤振情况。

同理，当以加工参数D(5 000 r/min，0.7 mm)进行实验时，振动加速度时域振动信号如图4(m)所示，工件发生轻微颤振，加速度的振动幅值在9.7 m/s2左右。从相平面图4(n)可知，工件的振动响应逐渐发散为多条曲线，振动信号的整个相平面图比较混乱，而图4(o)中Poincare截面图中的点也分布不均，该加工参数下对振动加速度信号进行最大Lyapunov指数(嵌入维数等于2)计算，其值为0.609。经傅里叶变换后得到的频谱如图4(p)所示，发现峰值不仅仅处于刀齿切削频率(5 000×4/60=333 Hz)及谐波频率(666 Hz)附近，而且能量向刀具结构频率(1 278 Hz)处积聚，因此在该铣削参数加工时发生颤振，处于非稳定切削状态。

通过振动信号的相平面图、Poincare截面图及频谱分析，发现当振动信号的最大Lyapunov指数超过0.59以后铣削处于非稳定切削状态。

3　非线性稳定域预测

基于相平面法、庞加莱法和频谱分析判定在该加工参数下是否发生了颤振后，对所有实验工况的工件振动加速度信号进行最大Lyapunov指数(嵌入维数等于2)计算，结果如图5所示。由图5可知，在每一转速下随着铣削深度的变化，其最大Lyapunov指数值也在随之变化，而在最大Lyapunov指数(嵌入维数等于2)为0.3～0.6附近发生突变。但是，当转速为5 250 r/min时折线在铣削深度0.9 mm处突然降低，其相平面图、Poincare截面图，频谱分析也说明此加工参数处于稳定切削状态，该点处于hof分叉叶瓣内，颤振稳定域边界附近，切削稳定稳定性边界存在一定的不确定性。

由于刀具转速和铣削深度是影响铣削稳定性的两个较主要的因素，那么通过插值方法可获得最大Lyapunov指数(嵌入维数等于2)与主轴转速和铣削深度的变化关系如图6所示，这样可以较直观看到主轴转速和铣削深度的变化对振动信号的非线性特征的影响。

图5　不同铣削深度时的振动加速度信号的最大Lyapunov指数Fig.5　The maximum Lyapunov exponent of vibration acceleration signals of different milling depth

图6　振动加速度信号的最大Lyapunov指数与主轴转速、铣削深度的变化关系Fig.6　Changes in the relationship between the maximum Lyapunov exponent of vibration acceleration signals with the speed of spindle and milling deep

由图6可知，在主轴转速4 000 r/min、5 500 r/min左右分别存在一个Lyapunov指数较小的区域，选用此处加工参数易获得较好的加工质量(该加工参数下的两端固定薄壁件的表面粗糙度小于其它工况下的表面粗糙度)。根据图6分别选定最大Lyapunov指数为0.6、0.59、0.61时作为阈值，画出等高线。以此作为稳定域边界，来分析铣削加工的稳态临界铣削深度与主轴转速的关系，图7为以不同最大Lyapunov指数值阈值时铣削稳定预测图与基于全离散法的铣削颤振稳定域的比较。

在图7中，其中标注“×”的表示在该切削参数实验中会导致颤振发生，而标注“○”则表示在该切削参数实验中切削处于稳定切削状态。其中，曲线1为采用全离散法得到的稳定性临界曲线；曲线4为振动信号最大Lyapunov指数值为0.61时所获得的稳定性临界曲线，这两条曲线的重合度较好，整体误差度小于9%，这也同时证明了此方法的有效性。曲线2、3为不同最大Lyapunov指数值阈值时铣削稳定预测图。本文所获的稳定性临界曲线在极限切深上低于全离散法曲线，这是由于加工系统非常复杂，存在着各种非线性因素所导致。

图7　铣削加工颤振稳定域临界曲线及实验结果(×非稳定切削参数；○稳定切削参数)Fig.7　Milling chatter stability region critical curve and experimental results (×The non stable cutting parameters；○ The stable cutting parameters)

本文基于实验所获得的稳定性曲线与一些实验结果有一定的误差，出现误差的原因可能由于以下的因素：① 由于动力学参数的变化铣削稳定域边界是在一定范围内变化的。② 利用插值技术所获得的最大Lyapunov指数与主轴转速及铣削深度的关系图本身具有一定的误差。③ 还需大量实验数据来进一步精确此颤振稳定性临界曲线。④ 实际加工中铣削稳定域的边界存在一定的不确定性。

4　结　论

本文以航空铝合金7075-T6两端固定薄板件为研究对象，研究了以最大李雅普诺夫指数作为阈值确定铣削颤振稳定性，并以此作为稳定判据。具体结论如下：

(1)以铝合金7075-T6两端固定薄板件振动信号为研究对象，基于相平面法、庞加莱法和频谱分析了不同加工参数时的振动信号，以此确定铣削过程是否发生颤振。

(2)Lyapunov指数是衡量系统动力学不确定性特征的重要判据，不规则振动中含有混沌，混沌是由铣削系统的非线性引起的，会导致不规则振动。稳定切削时与颤振切削时的不规则振动的程度不同。

(3)提取不同加工参数情况下两端固定薄壁件振动信号的最大Lyapunov指数，绘制并讨论了最大Lyapunov指数与主轴转速和铣削深度的变化关系。最大Lyapunov指数随铣削深度的增大而增大，实验中发现该指数会在某一切深时发生数值突变，由此可以将振动信号的最大Lyapunov指数作为稳定性的非线性判据。

(4)针对两端固定薄板件加工实验，以振动信号的最大Lyapunov指数值为0.61时作为铣削颤振的非线性判据，通过等高线法得到铣削颤振稳定性临界曲线。虽然本文所获得的稳定性临界曲线与基于全离散法所获取的临界曲线有一定的误差，但这种误差是可以接受的。另外，加工时可以以工件振动信号的最大Lyapunov指数作为加工参数优化的一种手段。

[1]GAO S H,LONG X H,MENG G.Noniinear response and nonsmooth bifurcations of an unbalanced machine-tool spindle-bearing system[J].Nonlinear Dynamics,2008,54:365-377.

[2]GRADISED J,GOVEKAR E,GRABED I.A chaotic cutting process and determining optimal cutting parameter values using neural networks[J].Institution Journal Machine Tools and Manufacture,1996,36(10):1161-1172.

[3]GRADISED J,GOVERKAR E,GRABED I.Time series analysis in metal cutting:Chatter versus chatter-free cutting[J].Mechannical Systems and Signal Processing,1998,12(6):839-854.

[4]DAVID E,GILSINN D E.Estimating critical hopf bifurcation parameters for a second-order delay differential equation with application to machine tool chatter[J].Nonlinear Dynamics,2002，30(2)：103-154.

[5]SZALAI R,STEPAN G,HOGAN S J.Global dynamics of low immersion high-speed milling[J].Chaos,2004,14(4):1069-1077.

[6]STEFANSKI A,KAPITANIAK T.Using chaos synchronization to estimate the largest Lyapunov exponent of non-smooth systems[J].Discrete Dynamics in Nature Sociences,2000,4(3):207-215.

[7]STEFANSKI A,KAPITANIAK T,DABROWSKI A.The largest lyapunov exponent of dynamical systems with time delay[C]//Proceedings of IUTAM Symposium on Chaotic Dynamics and Control of Systems and Processes in Mechanics,Rome,Italy,2003.6.

[8]LI Z Q,LIU Q.Solution and analysis of chatter tability for end milling in the time-domain[J].Chinese Journal of Aeronautics,2008,21(2):169-178.

[9]孔繁森，刘鹏，王晓明.切削振动加速度时间历程演化过程的动力学特征[J].振动与冲击,2011,30(7):10-15.

KONG Fansen,LIU Peng,WANG Xiaoming.Dynamic characteristics of evolution process of cutting vibration acceleration time history[J].Journal of Vibration and Shock,2013,30(7):10-15.

[10]刘志兵，王西彬. 微细铣削振动信号非线性特征的试验研究[J]. 兵工学报,2010,31(1):84-87.

LIU Zhibing,WANG Xibin.An expermental study on nonlinear character istic of vibration signal in Micro-end-milling Process[J].Acta Armamentar,2010,31(1):84-87.

[11]吴石，刘献礼，肖飞.铣削颤振过程中的振动非线性特征试验[J].振动测试与诊断,2012,32(6):935-940.

WU Shi,LIU Xianli,XIAO Fei.Experimental study of the nonlinear characteristics of vibration in milling chatter[J].Journal of Vibrat ion,Measurement and Diagnosis,2012,32(6):935-940.

[12]王平，陈蜀梅，王知人.大挠度简支矩形薄板受热力磁耦合作用分岔与混沌[J].振动与冲击,2013,32(7):129-134.

WANG Ping,CHEN Shumei,WANG Zhiren.Bifurcation and chaos of a thin rectangular plate simply upported with large deflection in a coupled environment of heating,force and magnetic field[J].Journal of Vibration and Shock,2013,32(7):129-134.

[13]YEH Y,CHEN C K,LAI H.Chaotic and bifurcation dynamics for a simply supported rectangular plate of themo-mechanical coupling in large deflection[J].Chaos,Solitons Fractals,2002,13(7):1493-1506.

[14]YEH Y.The effect of thermo-mechanical coupling for a simply supported orthotropic rectangular plate on nonlinear dynamics[J].Thin-Walled Structures,2005,43(8)：1277-1295.

[15]杨智春,张蕊丽.基于最大李雅普诺夫指数的壁板热颤振特性分析[J].西北工业大学学报,2009,27(6):770-776.

YANG Zhichun,ZHANG Ruili.Analysis of panel therm al flutter using max imum Lyapuno v ex ponent[J].Journal of Northwestern Polytechnical University,2009,27(6):770-776.

[16]闻邦椿，李以农，徐培民,等.工程非线性振动[M].北京：科学出版社，2007.

[17]DING Y,ZHU L M,ZHANG X J,et al.A full-discretization method for prediction of milling stability[J].International Journal of Machine Tools and Manufacture,2010,50(5):502-509.

Tests for milling chatter stability nonlinear criterion of thin parts

WU Shi,BIAN lijian,LIU Xianli,SONG Shenggang,JIANG Yancui

(School of Mechanical Engineering,Harbin University of Science and Technology,Harbin 150080,China)

The nonlinear dynamic behavior of milling process accompanies the entire cutting process.In order to accurately determine and predict the chatter stability of machining process,a nonlinear criterion for milling shatter stability of a thin part fixed at its two ends was studied with tests.Vibration signals of the thin part were taken as the study object in test,based on the phase plane method,Poincare method and spectra the vibration signals under different processing parameters were analyzed and plotted.The relationships between the maximum Lyapunov exponent and the spindle rotating speed,the former and the milling depth were discussed.Finally,taking the maximum Lyapunov exponent as the criterion,the milling chatter stability domain of the thin part was determined with the contour method,compared with the milling chatter stability domain obtained with the full discrete method,the nonlinear criterion for the milling chatter stability domain of a thin part made of aviation aluminum alloy 7075-T6 was obtained with tests.

thin part; milling chatter; nonlinear criterion; lyapunov exponent

国家自然科学基金资助项目(51275139);哈尔滨市科技创新人才研究基金项目(2011RFLXG019)

2015-04-07修改稿收到日期：2015-07-17

吴石 男，博士,教授，1971年8月生

TB122;TG506

A DOI:10.13465/j.cnki.jvs.2016.17.032