刘 宇, 王振宇, 杨慧刚, 张义民

(1.东北大学 机械工程与自动化学院， 沈阳 110819； 2.菏泽学院 蒋震机电工程学院，山东 菏泽 274015)

随着微铣刀加工技术的不断提高，微铣削的应用领域也越来越广，如芯片上的微通道总线，微型燃料电池，微喷头，光纤微孔等。然而铣削加工过程中的再生颤振问题仍严重制约着机床的加工效率和工件的加工精度。

目前颤振预测的方法主要有四种：Altintas等[1]提出的定向法，该方法根据定向系数将二维铣削动力学参数定向到切削合力方向，又被称为准单自由度法；Altintas等[2]提出的零阶近似法，该方法给出了精确计算时变切削力的解析表达式，为精确近似时变的切削力提供了新的思路；Elbeyli等[3]提出的半离散法和Ding等[4]提出的全离散方法，将切削过程均分成若干等份，通过求解延时微分方程获得铣削加工的颤振稳定性叶瓣图。此后，研究发现由于刀具的磨损、材料的分布不均匀、测量误差、噪声干扰、加工环境的变化、高主轴转速下系统的非线性问题以及是否使用冷却液等，致使实际加工过程中切削力系数、固有频率等都是不确定参数[5-6]，这些不确定参数又可根据是否与时间有关分为随机参数和时变参数。稳定性叶瓣图对模态参数的变化特别敏感，因此这些参数的变动必然对颤振稳定性预测产生很大的影响。在上述研究中，所有的方法都将加工系统的结构参数、材料切削性能参数等看作确定性量。

针对加工过程中参数的随机性问题，有些学者提出了包含随机参数的颤振稳定性预测方法。Duncan等首次提出了参数的不确定性对稳定性叶瓣图的影响，指出加工过程中的频响函数和切削力系数都是随机参数。Park等[7]将鲁棒稳定性引入到颤振预测模型中，并将固有频率和切削力系数当作不确定参数，通过对这些参数取最大最小值并结合鲁棒稳定性原理得出颤振稳定性叶瓣图。Sims等[8]将模糊算法引入到颤振预测中，用隶属度表示参数的变化范围，并得到模糊稳定性叶瓣图。Totis[9]将概率算法引入车削过程中，建立了铣削动态不确定模型，分析了系统的颤振鲁棒稳定性。Liu[10]将动态结构的可靠性分析引入车削系统的颤振预测中，得到一系列不同可靠度下的颤振稳定性叶瓣图。上述颤振预测方法虽然考虑到参数的随机特性，但在其预测中都将随机参数都看作时不变的，并没有考虑因刀具磨损导致的刀刃半径和切削力系数等参数的时变特征。

刀具磨损是一个时变的过程，因此在其影响下的刀刃半径和切削力系数也是时变的。Li等[11]研究发现微铣削中刀具的磨损使得刀刃半径和X，Y，Z方向上的切削力均逐渐增。Karandikar等[12]研究发现切向切削力和径向切削力随刀具磨损量的增加而增加。

研究表明微铣削中刀刃半径的大小对切削力和颤振稳定性都有显著的影响。Bierman等[13]研究表明刀刃半径的大小对系统的颤振稳定性有显著的影响。Afazov等[14-15]发现微铣削中刀刃半径越大切削力越大，并用非线性方程描述了刀刃半径与切削力系数之间的关系，绘制出了不同刀刃半径下的稳定性叶瓣图。Jin等[16]建立了包含刀刃半径的切削力系数模型，考虑略刀刃半径对过程阻尼的影响进而得到对应的颤振稳定性叶瓣图。然而他们只是将刀具的大小当作定值忽略了刀具磨损导致的刀刃半径的随机性和时变问题。

针对加工参数的随机特性和刀具磨损导致的参数的时变问题，本文将动态结构的时变可靠性分析思想引入到铣削加工系统的颤振稳定性分析中，用伽马过程描述刀刃半径随切削时间的变化关系，进而得到了不同切削刻下切削力系数的大小。建立了颤振时变稳定性和颤振时变可靠性模型，分析了微铣削加工中的时变稳定性和时变可靠性并绘制出系统的时变稳定性和时变可靠度曲线。

本文的结构如下：第一部分为微铣削的力学建模，将切削力系数用含刀刃半径的方程表示；第二部分为刀具的磨损量建模，用Gamma过程描述了刀刃半径随时间的变化关系；第三部分为系统颤振时变稳定性建模与分析与模型简化，给出了包含时间信息的极限切深与颤振可靠度的数学表达式；第四部分建立了系统给的时变可靠性模型；第五部分为算例研究；第六部分为结论。

1 刀刃半径时变规律建模

Gamma过程是一个具有独立、非减且时间和状态都是连续的随机过程。其增量服从Gamma 分布。由Li等[17-18]的论文可知刀具磨损是一个典型的连续时间、 连续状态的随机过程，并且由于其不可自我修复，因此也是一个增量非减的过程，故本文用Gamma过程来描述刀具的磨损量与时间之间的关系。定义r(t)(其中t的单位为小时)为t时刻刀刃半径的增量。由Gamma过程的定义可知其概率密度函数为

(1)

式中：Ga为Gamma分布函数；Γ为Gamma函数；μ，v(t)分别为伽马分布中的尺度参数的倒数和形状参数；I(0, ∞)为当x∈[0, ∞]时其值为1, 否则为0的函数。R(t)的均值和方差分别表示为

(2)

Moriwaki等[19]研究表明，期望的劣化值与能量规律成正比，所以刀刃半径增量的期望值可以表示为

(3)

式中：μ、b和c都是大于零的实数，可由曲线拟合得到。所以Gamma过程中的形状参数可以表示为

v(t)=ctb

(4)

所以在t时刻是刀刃半径的大小可表示为

(5)

式中：r0为加工初始时刻刀刃半径的大小；t时刻是刀刃半径的均值可表示为

(6)

微铣刀磨损前后的形状如图1所示。

图1(a)为微铣刀整体尺寸，图1(b)为刀具磨损前的刀刃半径大小，图1(c)为磨损后刀刃半径的大小。

2 包含刀刃半径时变切削力系数建模

微铣削中由于切深较小一般在25 μm以下，微铣刀的刀刃半径一般在3～20 μm[19]，因而刀刃半径的变化直接影响切削力的大小。故，如何精确的建立刀刃半径与切削力系数之间的关系成为了研究刀刃半径对加工过程稳定性影响的关键。本文结合了 Jin等论文中的切削力系数模型，研究了加工时间与切削力系数之间的关系具体表达式如下

(a)(b) (c)

图1 微铣刀磨损前后刀刃半径的变化对比

Fig.1 Comparison of cutting edge radius between new and wear tool in micro-milling

Kt(h,r(t))=αthdt+βthptr(t)qt=

(7)

Kr(h,r(t))=αrhdr+βrhprr(t)qr=

(8)

式中：Kt，Kr分别为切向切削力系数和径向切削力系数；α,β,d,p,q为常数。本文Deform-2D仿真技术获取微切削材料AL6061的切削力，进而求取其切削力系数，刀具前角为5°，后角为7°，刀刃半径r在[1,20] μm内变化，切厚h在[0.5, 20]μm内变动，进给速度为4.5 μm每刃，其余参数与Jin等均相同，网格单元尺寸为0.05μm。各参数拟合后的大小如表1所示。

表1 切削力系数的大小

3 微铣削中颤振时变稳定性建模

3.1 颤振稳定性建模

如图2所示，我们将mx,my,kx,ky,cx,cy定义为两自由度铣削模型在X方向和Y方向上的模态质量、模态刚度和模态阻尼，h表示曲面法向方向的动态切削厚度，h=x(t)-x(t-T)，则系统的动力学方程如下

(9)

(10)

式中：N为铣刀刀刃数；φi为第i个齿的回转角，Ft,i,Fr,i分别为第i个齿在切向和径向方向上的分力，其大小可以表示为

图2 微铣削力学模型

(11)

Fr,j=KrFt,j

(12)

式中：a为轴向切深；将上述公式写成矩阵的形式如下

(13)

式中：Kt为切向切削力系数；Δx,Δy分别为X和Y方向上当前和上一周期的位移差，即切厚。αxx,αxy,αyx,αyy为平均定向因子，其表达式分别为

(14)

(15)

(16)

(17)

式中：φj为浸入角；φst,φex分别为铣刀的切入角和切出角，槽铣时其大小分别为0和π。Kn为径向切削力系数与切向切削力系数的比值，具体表达式如下

Kn(h,r(t))=Kt(h,r(t))/Kr(h,r(t))

(18)

将上述时变平均定向因子代入式(20)和式(21)，然后将代入后的值代入式(19)时变特征方程，如下所示

a0(r(t))Δ2+a1(r(t))Δ+1=0

(19)

a0(r(t))=φxx(iω)φyy(iω)×

(αxx(r(t))αyy(r(t))-αxy(r(t))αyx(r(t)))

(20)

a1(r(t))=αxx(r(t))φxx(iω)+αyy(r(t))φyy(iω)

(21)

式中：φxx,φyy分别为铣削系统在X和Y方向上的频响函数。为绘制系统的稳定性叶瓣图我们用切削时间t处刀刃半径的均值代替该时刻的刀刃半径大小。求解式(19)可得系统的时变特征值λ(E[r(t)])，及其对应的时变特征值的实部Re(λ(E[r(t)]))和虚部Im(λ(E[r(t)]))。则系统在t时刻的极限切深可表示如下

(22)

(23)

(24)

3.2 颤振稳定性简化模型

为便于计算微铣削颤振时变可靠度，我们将Kn看作常数，大小为加工初始时刻切向切削力和径向切削力系数的比值。由式(14)～式(17)可知铣削系统的平均定向因子αxx,αxy,αyx,αyy均为常数，则同一主轴转速在不同切削时刻下对应的特征值也是固定的，所以系统的极限切深和主轴转速可简化为

(25)

极限切深对应的主轴转速可表示为

(26)

(27)

为检验简化后的模型精度，本文令Kn的值为0.357 7,并在同一加工条件下选取主轴转速分别为40 000和60 000时比较了颤振稳定性模型化简前后其极限切深随加工时间的变化，如图3所示。其不同切削时刻下极限切深的波动大小如表2所示。

(a)(b)

图3 简化前后的极限切深对比

Fig.3 Comparison between model and simplified model of the limit cut depth

由图3和表2可知颤振时变稳定性模型简化后对其极限切深的影响很小，最大波动率为4.57%。所以本文在进行时变颤振稳定性预测时将径向切削力系数和切向切削力系数的比值看作定值。

表2 模型简化前后极限切深的波动率

4 微铣削颤振时变可靠性建模与计算

4.1 颤振时变稳定性建模

定义铣削加工过程中当轴向切深小于对应主轴转速下的极限切深时系统可靠，否则不可靠，则系统的功能函数可以表示为

g(X)=blim-b0

(28)

式中：b0为铣削过程中加工的轴向切深；g(X)为系统失效的概率密度函数，系统的时变可靠度可以表示为

(29)

式中：XR为积分区间。定义加工过程中考虑时变参数的颤振可靠性分析为颤振时变可靠性分析。

4.2 颤振时变稳定性计算

则将式(25)代入式 (28)和式(29)可得包含时变切削力系数的功能函数和颤振时变可靠度，如式(30)、式(31)所示

(30)

(31)

由式(31)可知，在加工参数和主轴转速一定时，刀刃半径的大小是影响系统颤振可靠度大小的唯一因素。由式(7)和式(8)知在加工参数和主轴转速一定时系统的切削力系数随刀刃半径的增加而单调递增，因此，在加工参数和主轴转速一定时，刀刃半径的概率分布直接决定了系统的颤振可靠度，故系统的颤振时变可靠度可简化成如下形式

(32)

式中：0～rtb为积分区间；XR，rtb为给定主轴转速下的极限切深为b0时对应的刀刃半径大小。可由式(33)求得。

(33)

由于式(33)不存在解析解，故本文采用二分法计算给定切深和主轴转速下对应的颤振临界点的刀刃半径rtb大小。当rtb有解时，对式(32)进行定积分可得该时刻下的颤振可靠度

(34)

式中：γ为不完全Gamma函数具体如式(35)所示

(35)

式中：a,z,η为已知量。当式(33)无解时，若在二分法的求解范围[ra,rb]的端点值大于零则Pr(X)=0，否则Pr(X)=1。

4.3 给定切深和主轴转速下的极限切深和颤振可靠度曲线

颤振时变可靠度是用来描述加工过程中不发生颤振的概率。计算流程如图4所示。

图4 颤振时变可靠度流程图

在铣削过程中轴向切深b0主轴转速Ω和进给速度υ一定时，将式(6)计算得到的不同时刻下的刀刃半径的均值代入式(9)和式(10)得时变切削力系数，然后将其先后代入式(14)～式(21)的系统在不同切削时间下的特征值。然后将其代入式(25)可得t时刻下的极限切深，由式(26)和式(27)可得对应极限切深下的主轴转速，进而可得系统在t时刻的稳定性叶瓣图。将计算得到的不同切削时间下的特征值代入式(33)得到积分上限rtb由式(34)可得在给定加工条件下切削时间为t时系统的颤振可靠度值。然后以时间t为横坐标，以t时刻所对应的可靠度为纵坐标画图，得系统在给定切深和主轴转速下的颤振时变可靠度曲线。

5 算例研究

5.1 拟合Gamma过程参数

Gamma过程中参数的具体数值可由数据拟合得到。实验中采用WC/Co材料的不带涂层的两刃铣刀，刀刃半径为3 μm，切削材料为AL6061。主轴转速为65 000 r/min，进给速度为4.5 μm每齿，切深为45 μm。 实验中刀刃半径与切削时间的关系[20]如图5(a)中的圆点所示。我们通过数据拟合方式获得刀刃半径的Gamma过程的模型参数，并绘制了刀刃半径在不同切削时刻下的Gamma分布和均值，如图5(a)中各时刻对应的曲线和方块所示。其中“方块”在“圆”和“三角形”中间，如图5(a)中局部放大的部分。图5(b) 为Gamma过程预测的刀刃半径在不同切削时刻下的标准差。

(a)(b)

图5 刀刃半径的拟合结果

Fig.5 Fitting results of cutting edge radius

Gamma过程的拟合参数大小如表3所示。

表3 Gamma过程拟合参数

5.2 微铣削加工过程中颤振时变稳定性预测

本文选用超精密立式数控铣床(Kern Micro 2255) 进行铣削加工的颤振稳定性研究，该铣削加工系统的刀尖点模态参数[21]如表4所示。

由式(25)～式(27)可得系统在不同切削时间下的颤振稳定性叶瓣图如图6所示。图6显示随着切削时间的增加其叶瓣图逐渐降低。

表4 铣削系统的刀尖点模态参数

5.3 计算给定切深和主轴转速下的系统颤振时变可靠度

图7分别比较了颤振时变可靠度与刀刃半径和极限切深随切削时间的变化关系，图7(a)显示随着切削时间的增加系统的刀刃半径逐渐增加，系统的颤振可靠度逐渐降低。图7(b)中显示随着切削时间的增加系统的极限切深和颤振可靠度均逐渐降低。

图6 不同加工时刻下的稳定性叶瓣 Fig.6 Time-varyingchatterstabilitylobediagram(a)(b) 图7 颤振可靠度与刀刃半径和极限切深随切削时间的对比关系 Fig.7 Therelationshipbetweenreliabilityandtooledgeradius,criticaldepthwithcuttingtime

6 结 论

微铣削加工过程由于刀具的尺寸小，主轴转速高使得加工过程中道具磨损较快。而刀具磨损使得刀刃半径和切削力系数随加工时间明显增加，进而对系统的颤振稳定性预测结果产生很大的影响，使得传统的颤振稳定性预测方法很快失效。本文提出的颤振时变稳定性预测和时变可靠性分析方法将刀具磨损的影响考虑到颤振稳定性预测中。在不同的加工时刻上获得了不同的叶瓣图。并在给定切深和主轴转速的情况下给出了系统的极限切深和颤振可靠度随切削时间的变化关系，使得颤振预测的结果能够准确的反应不同加工时间上的颤振稳定性。

[2] ALTINTA Y, BUDAK E. Analytical prediction of stability lobes in milling[J]. CIRP Annals-Manufacturing Technology, 1995, 44(1): 357-362.

[3] ELBEYLI O, SUN J Q, ÜNAL G. A semi-discretization method for delayed stochastic systems[J]. Communications in Nonlinear Science and Numerical Simulation, 2005, 10(1): 85-94.

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

[5] DUNCAN G, KURDI M, SCHMITZ T L, et al. Uncertainty propagation for selected analytical milling stability limit analyses[C]//34th North American Manufacturing Research Conference.[S.l.]:[s.n.], 2006.

[6] ZHANG X, ZHU L, ZHANG D, et al. Numerical robust optimization of spindle speed for milling process with uncertainties[J]. International Journal of Machine Tools and Manufacture, 2010，61(1): 9-19.

[7] PARK S S, RAHNAMA R. Robust chatter stability in micro-milling operations[J]. CIRP Annals-Manufacturing Technology, 2010, 59(1): 391-394.

[8] SIMS N D, MANSON G, MANN B. Fuzzy stability analysis of regenerative chatter in milling[J]. Journal of Sound and Vibration, 2010, 329(8): 1025-1041.

[9] TOTIS G. RCPM—a new method for robust chatter prediction in milling[J]. International Journal of Machine Tools and Manufacture, 2009, 49(3/4): 273-284.

[10] LIU Yu. Mechanical systems and signal processing[J]. Mechanical Systems and Signal Processing, 2016，66/67(3): 232-247.

[11] LI H Z, ZENG H, CHEN X Q. An experimental study of tool wear and cutting force variation in the end milling of Inconel 718 with coated carbide inserts[J]. Journal of Materials Processing Technology, 2006, 180(1/2/3): 296-304.

[12] KARANDIKAR J M, ZAPATA R E, SCHMITZ T L. Combining tool wear and stability in high-speed machining performance prediction[D]. Gainesville: University of Florida, 2010: 7-10.

[13] BIERMANN D, BASCHIN A. Influence of cutting edge geometry and cutting edge radius on the stability of micromilling processes[J]. Production Engineering, 2009, 3(4/5): 375-380.

[14] AFAZOV S M, RATCHEV S M, SEGAL J. Modelling and simulation of micro-milling cutting forces[J]. Journal of Materials Processing Technology, 2010, 210(15): 2154-2162.

[15] AFAZOV S M, ZDEBSKI D, RATCHEV S M, et al. Effects of micro-milling conditions on the cutting forces and process stability[J]. Journal of Materials Processing Technology, 2013,213(5):671-684.

[16] JIN X, ALTINTAS Y. Chatter stability model of micro-milling with process damping[J]. Journal of Manufacturing Science & Engineering, 2013, 135(3): 1-9.

[17] 李常有, 张义民, 王跃武. 恒定加工条件及定期补偿下的刀具渐变可靠性灵敏度分析方法[J]. 机械工程学报, 2012, 48(12): 162-168.

LI Changyou, ZHANG Yimin, WANG Yuewu. Gradual reliability and its sensitivity analysis approach of cutting tool in invariant machining condition and periodical compensation[J]. Journal of Mechanical Engineering, 2012, 48(12): 162-168.

[18] NOORTWIJK J M V, WEIDE J A M V D, KALLEN M J, et al. Gamma processes and peaks-over-threshold distributions for time-dependent reliability[J]. Reliability Engineering & System Safety, 2007, 92(12): 1651-1658.

[19] MORIWAKI T, TOBITO M. A new approach to automatic detection of life of coated tool based on acoustic emission measurement[J]. Journal of Engineering for Industry, 1990，112(3)：212-218.

[20] MALEKIAN M, PARK S S, JUN M B G. Modeling of dynamic micro-milling cutting forces[J]. International Journal of Machine Tools & Manufacture, 2009， 49(7): 586-598.

[21] THEPSONTHI T, ÖZELT T. 3D Finite element process simulation of micro-end milling Ti-6Al-4V titanium alloy: experimental validations on chip flow and tool wear[J]. Journal of Materials Processing Technology, 2015，221：128-145.