李鹤, 韩萍, 张玲利
东北大学 机械工程与自动化学院, 辽宁 沈阳 110819
收稿日期: 2015-02-16
基金项目: 国家自然科学基金资助项目(51175071); 中央高校基本科研业务费专项资金资助项目(N120203001, N130803001).
作者简介: 李鹤(1975-), 男, 河南方城人, 东北大学教授, 博士生导师.

摘要: 利用切比雪夫分割法离散延时微分方程中的时间周期项，然后采用二阶Newton差商方法求出每个时间元端点处铣刀的运动向量，进而得到铣削系统的传递矩阵，并用Floquet理论判断铣削系统的稳定性.利用以上方法，研究了一个两自由度铣削系统稳定性问题，得到了系统稳定性叶瓣图.结果表明：在同等计算精度的要求下，切比雪夫分割法比平均分割法所取的点数要少，计算效率较高；将时间周期分割成相同的份数，切比雪夫分割法计算精度高.
关键词: 铣削切比雪夫分割法传递矩阵二阶Newton差商法稳定性叶瓣图
Milling Stability Analysis Based on Chebyshev Segmentation
LI He, HAN Ping, ZHANG Ling-li
School of Mechanical Engineering & Automation, Northeastern University, Shenyang 110819, China
Corresponding author: LI He, professor, E-mail: hli@mail.neu.edu.cn
Abstract: Using the Chebyshev segmentation method to discrete the time period contained in delay differential equation, the Newton second-order difference quotient method was used to calculate the motion vector of cutter at each time endpoint, then the transfer matrix of milling system was obtained, and the system stability was determined by the Floquet theory. Using the above methods, a two-degrees of freedom milling system stability issues was investigated to get system stability lobe diagrams. The results showed that under the same calculation accuracy, the points needed to cut the time period by the Chebyshev segmentation are less than those by the average segmentation, and the computational efficiency of the Chebyshev segmentation is higher. When the time period is divided into the same parts, the stability lobe diagrams gotten by the Chebyshev segmentation method are more accurate than that of the other.
Key words: millingChebyshev segmentationtransfer matrixNewton second-order difference quotient methodstability lobe diagram 在模具、电子、航空、航天等领域，常常要求高效、高精度地从毛坯中铣削去除大量材料.实践表明，颤振是影响高效、高精度铣削的重要因素之一^[1].颤振使得刀具与工件的相对位置发生改变，降低铣削效率和加工表面的精度，还会加剧刀具的磨损甚至破坏，缩短机床和刀具的使用寿命.因此，铣削过程稳定性分析和铣削颤振的控制理论研究已经成为提高铣削加工效率的重要手段之一，具有重要的研究价值.
铣削稳定性分析是合理选取主轴转速和切深等加工参数以达到优化工艺参数，避免发生颤振，提高生产效率的重要途径.铣削过程中多齿刀具旋转而引起的切削力周期变化和切削过程间断等因素，使得铣削稳定性分析比连续正交切削稳定性分析困难得多.
目前国内外学者针对这类铣削稳定性问题，已经提出许多近似计算方法预测由临界加工参数构成的铣削过程稳定性边界曲线即稳定性叶瓣图(stability lobe diagram，SLD)^[2].文献[3]提出了用于判断延时系统稳定性的半离散法(semi - discretization method，SDM).文献[4]针对单自由度断续切削的情况，提出了时间有限元分析法(temporal finite element analysis，TFEA).文献[5]基于铣削过程动态响应的直接积分格式提出了适用于大(小)径向切深、大(小)轴向切深等多种工况半解析的全离散法(full-discretization method，FDM).文献[6]证明，全离散法与零阶半离散法具有同阶局部截断误差.
以上的铣削稳定性预测方法，不论是半离散法、时间有限元法还是全离散法都是将铣削延时微分方程的时间周期平均分割.本文提出的切比雪夫分割方法是将时间周期T不均匀分割成m个时间元，刀具任意时刻的状态项和时间周期项用时间元端点的状态项线性表示，延时项采用二阶Newton差商的方法用时间元端点的状态项插值表示.首先求出铣削系统在每个时间元上的传递矩阵，再求铣削系统在一个时间周期内的传递矩阵，再用Floquet理论^[7]来判断铣削系统的稳定性.通过时域内仿真结果可以看出：在相同精度要求下，切比雪夫不平均分割时间周期方法的计算效率比平均分割时间周期方法的计算效率高.
1 切比雪夫分割法逼近法近似求函数解析表达式是在函数定义域区间段上取n+1个互异的节点，给定这些节点的函数值，采用Lagrange或者Newton插值方法^[8]构造出一个唯一的n次插值多项式来逼近原函数.采用均匀取点法构造出来的n+1个节点的插值型积分公式一般具有n次代数精度.如果适当地选取节点来构造求积公式，就可以得到代数精度更高的求积公式.切比雪夫分割法^[9]就是这样一种取点分割方法，它的取点分割过程如下.
如图 1所示，将半圆弧按照角θ分割，圆弧上的分割点在x轴上的投影点即是该区间段[-1,1]上的分割点.
图1(Fig. 1)

图1 切比雪夫分割图(N0=4) Fig. 1Image of Chebyshev segmentation

在区间[a,b]上，x轴上取点的位置表达式为

2 全离散法只考虑再生颤振的铣削延时微分方程为

式中：A为铣削系统特性的常数矩阵；B为再生效应动态铣削力决定的周期矩阵B(t)=B(t+T)；T为时间周期等于延时量.
切比雪夫分割方法将时间周期T不均匀分割为m个时间元：

在第一个时间元内t_k≤t≤t_k+1(k=0，1，…，m-1)，式(4)的响应可以写成直接积分格式：

X(t_k)表示在t_k时刻刀具的状态向量^[6]，则在t_k+1时刻刀具的状态向量X(t_k+1)为

此时式(10)中的积分时间区间为[0,Δt_k].对周期系数矩阵B(t_k+1-ξ)，X(t_k+1-ξ)在时间元t_k≤t≤t_k+1上做线性逼近^{[10, 11, 12]}；对于延时项X(t_k+1-ξ-T)可以用时间元[t_k-T,t_k+1-T]两端点时刻铣刀的状态向量X(t_k+1-T)和X(t_k-T)采用二阶Newton差商的方法线性表示，见式(11)；将式(11)带入式(10)得式(12).


式(12)中的变量F_k+1，F_k-1，F_0,k分别为

设[I-F_k+1]^-1存在，则铣削系统每个时间元上的传递矩阵Ψ_k可以表示为

铣削系统在一个时间周期T上的传递矩阵Ψ可以用这个周期内时间元的传递矩阵序列Ψ_k (k=0，1，…，m-1)表示，即

最后，根据Floquet理论，铣削系统的稳定性可由系统传递矩阵Ψ的特征值决定：若传递矩阵Ψ的所有特征值的模长均小于1，则系统稳定；否则系统不稳定.
3 切比雪夫分割法下的稳定性分析3.1 铣削系统参数两自由度铣削系统如下^[12]：两个刀齿的槽式直铣刀；铣刀齿数N=2；质量m_tx=0.03993kg，m_ty=0.03993kg；阻尼比ζ=0.011；系统的固有频率f_n=ω_n/(2π)=922Hz；铣刀转速Ω变化范围为0.5×10⁴~2.5×10⁴r/min；切深w的变化范围为0~0.01m；K_t=6×10⁸N/m²；K_n=2×10⁸N/m²；铣削系统特性的常数矩阵A：

再生效应动态铣削力决定的周期矩阵B：

K_t是线性化切向切削力系数，K_n是线性化法向切削力系数，a/D表示铣刀的径向侵入比.
3.2 稳定性分析结果将转速-轴向切深平面平均划分为200×100网格点，即转速步长100r/min，轴向切深步长为0.1mm.分别计算径向侵入比等于0.1和0.2两种工况下的铣削稳定性极限.根据前面所述的理论基础和算法原理，编写计算程序，利用计算机绘制稳定性叶瓣图,如表 1所示.
表1(Table 1)

表1 两自由度铣削模型的算法比较Table 1 Algorithm comparison of two degrees of freedom milling model 切比雪夫分割法下的全离散(m=30) 平均分割法下的全离散法(m=40) a/D=0.1 计算时间/s 242.819 501.195 a/D=0.2 计算时间/s 251.978 509.060

表1 两自由度铣削模型的算法比较 Table 1 Algorithm comparison of two degrees of freedom milling model

由表 1可以看出在用叶瓣图来预测铣削系统的稳定性时，切比雪夫分割法取的点数比平均分割法取的点数少的情况下，叶瓣图仍能达到相同的精度.由两种分割方法的计算时间的对比可知，获得相同精度叶瓣图时，切比雪夫分割法的计算效率比平均分割法计算效率高.
图 2表示a/D=0.1，m =40时切比雪夫分割法和平均分割法得到的稳定性叶瓣对比图，这两种分割方法得到的铣削系统的稳定性范围有所不同.在图 2中取4个工作点：点1，点2，点3，点4，这4个工作点对应的铣削参数和两种分割方法预测的4个工作点的稳定性如表 2所示.对表 2中所列工作点进行时域仿真，仿真结果如图 3所示.工作点的仿真结果如图 3所示.由图 3可知，工作点1，点2，点3对应的铣削过程是稳定的，点4对应的铣削过程是不稳定的.图 3表明4个工作点的时域仿真结果与表 2中切比雪夫分割法预测的4个工作点的稳定性完全吻合，而与平均分割法预测的4个工作点的稳定性相差较多.因此，取相同分割点时，切比雪夫分割法对铣削系统的稳定性预测结果比平均分割法对铣削系统的稳定性预测结果更准确.
图2(Fig. 2)

图2 两种分割方法下的稳定性叶瓣图对比(顺铣N=2，a/D=0.1 ) Fig. 2 Comparison of two different segmentation’s stability lobes diagram

表2(Table 2)

表2 两种取点分割法得到的部分工作点的稳定性Table 2 Stability of sectional working points for two way to take segmentation 工作点 Ω w 切比雪 夫分割法 平均分 割法 r·min-1 mm 点1 15000 1.5 稳定 不稳定 点2 20000 1 稳定 稳定 点3 20000 2 稳定 不稳定 点4 20000 3 不稳定 不稳定

表2 两种取点分割法得到的部分工作点的稳定性 Table 2 Stability of sectional working points for two way to take segmentation

图3(Fig. 3)

图3 部分工作点的时域仿真结果图 Fig. 3 Time domain simulation’s results for part of working points (a)—点1； (b)—点2; (c)—点3； (d)—点4.

4 结 论本文介绍了切比雪夫分割方法和Newton二阶差商法对离散项进行插值积分的方法.分析两自由度铣削系统示例，由绘制的系统稳定性叶瓣图可知：在达到相同预测精度的条件下，切比雪夫不平均分割法需要取的点数比平均分割法所需的点数少，并且切比雪夫不平均分割法的计算效率比平均分割法计算效率高；在分割点数相同的情况下，切比雪夫不平均分割法预测的系统稳定性比平均分割法预测的系统稳定性更加准确.

参考文献

[1]	Tlusty J,Polacek M.The stability of machine tools against self-excited vibrations in machining [J].International Research in Production Engineering,ASME, 1963,1(8):465-474.(1)
[2]	Merrit H E.Theory of self-excited machine tool chatter [J].ASME Journal of Engineering for Industry,1965,8(7):447-454.(1)
[3]	Insperger T,Stépán G.Semi-discretization method for delayed systems [J].International Journal for Numerical Methods in Engineering,2002,55(5):503-518.(1)
[4]	Bayly P V,Halley J E,Mann B P,et al.Stability of interrupted cutting by temporal finite element analysis [J].Journal of Manufacturing Science and Engineering,Transactions of the ASME,2003,125(2):220-225.(1)
[5]	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.(1)
[6]	Insperger T.Full-discretization and semi-discretization for milling stability prediction:some comments [J].International Journal of Machine Tools and Manufacture, 2010,50(7):658-662.(2)
[7]	Farkas M.Periodic motions [M].New York:Springer-Verlag,1994.(1)
[8]	张铁,阎家斌.数值分析 [M].北京:冶金工业出版社,2007:141-187.(Zhang tie,Yan Jia-bin.Numerical analysis [M].Beijing:Metallurgical Industry Press,2007:141-187.)(1)
[9]	Yang W Y,Cao W W,Chung T S,et al.Applied numerical methods using MATLAB[M].New York :John Wiley & Sons,2005:132-135.(1)
[10]	Tan S J,Zhong W X.Precise integration method for duhamel terms arising from non-homogenous dynamic systems [J]. Chinese Journal of Theoretical and Applied Mechanics,2007,39(3):374-381.(1)
[11]	Zhong W X,Williams F W.A precise time step integration method [J].Proceedings of the Institution of Mechanical Engineers.Part C:Mechanical Engineering Science,1994,208(6):427-430.(1)
[12]	Ding Y,Zhu L M,Zhang X J.Second-order full-discrimination method for milling stability prediction [J].International Journal of Machine Tools & Manufacture,2010,50(10):926-932.(2)

1

本文献在全文中的定位：
... 去除大量材料.实践表明，颤振是影响高效、高精度铣削的重要因素之一^[1].颤振使得刀具与工件的相对位置发生改变，降低铣削效率和加工表面 ...



1

本文献在全文中的定位：
... 构成的铣削过程稳定性边界曲线即稳定性叶瓣图(stability lobe diagram，SLD)^[2].文献[3]提出了用于判断延时 ...



1

本文献在全文中的定位：
... stability lobe diagram，SLD)^[2].文献[3]提出了用于判断延时系统稳定性的半离散法(semi - discretization method，SDM) ...



1

本文献在全文中的定位：
... 了用于判断延时系统稳定性的半离散法(semi - discretization method，SDM).文献[4]针对单自由度断续切削的情况，提出了时间有限元分析法(temporal finite e ...



1

本文献在全文中的定位：
... 的情况，提出了时间有限元分析法(temporal finite element analysis，TFEA).文献[5]基于铣削过程动态响应的直接积分格式提出了适用于大(小)径向切深、 ...



2

本文献在全文中的定位：
... )轴向切深等多种工况半解析的全离散法(full-discretization method，FDM).文献[6]证明，全离散法与零阶半离散法具有同阶局部截断误差. ...

... X(t_k)表示在t_k时刻刀具的状态向量^[6]，则在t_k+1时刻刀具的状态向量X(t ...



1

本文献在全文中的定位：
... 传递矩阵，再求铣削系统在一个时间周期内的传递矩阵，再用Floquet理论^[7]来判断铣削系统的稳定性.通过时域内仿真结果可以看出：在相同精度 ...



1

本文献在全文中的定位：
... +1个互异的节点，给定这些节点的函数值，采用Lagrange或者Newton插值方法^[8]构造出一个唯一的n次插值多项式来逼近原函数.采用均匀取点法构 ...



1

本文献在全文中的定位：
... 来构造求积公式，就可以得到代数精度更高的求积公式.切比雪夫分割法^[9]就是这样一种取点分割方法，它的取点分割过程如下. ...



1

本文献在全文中的定位：
... )在时间元t_k≤t≤t_k+1上做线性逼近^{[10, 11, 12]} ...



1

本文献在全文中的定位：
... _k+1上做线性逼近^{[10, 11, 12]}；对于延时 ...



2

本文献在全文中的定位：
... ^{[10, 11, 12]}；对于延时项X(t_k+1-ξ-T)可 ...

... 两自由度铣削系统如下^[12]：两个刀齿的槽式直铣刀；铣刀齿数N=2；质量m ...