CN109840380B - 一种考虑多模态振动与工件加工响应的稳定性预测方法 - Google Patents

Abstract

本发明涉及机械加工中的金属切削领域，具体涉及一种考虑刀具‑工件多阶模态振动与工件实时加工频率响应的稳定性预测方法。本发明所述加工测量系统包括工件、力锤、电容传感器、阻尼杆、磁力表座、第一夹具、第二夹具、电荷放大器、PC机和刀具；所述第一夹具和第二夹具设置在工件下端，所述第一夹具和第二夹具上均设有标尺，所述电容传感器设置在阻尼杆上，所述阻尼杆设置在磁力表座上；本发明构建了刀具‑工件四自由度铣削模型同时考虑到刀具螺旋角滞后效应，并基于此引入便于计算侵入体积的过程阻尼模型，本发明利用等效三角截面积能快速确定过程阻尼，进一步提升稳定性边界预测的准确性。

Description

一种考虑多模态振动与工件加工响应的稳定性预测方法

技术领域

本发明涉及机械加工中的金属切削领域，具体涉及一种考虑刀具-工件多阶模态振动与工件实时加工频率响应的稳定性预测方法。

背景技术

由于轻量化设计的日益普及，薄壁工件在航空航天、通信等领域得到了广泛的应用。薄壁工件具有复杂几何拓扑结构和刚性差的特点，给铣削加工带来了很大的挑战。其切削过程涉及到机床、铣刀、工件和夹具的协同配合，最大限度地提高机床的性能，提高加工效率，使整个过程稳定可控具有重要意义。国际著名期刊《Journal of ManufacturingScience and Engineering》于2012年发表了“Discrete-Time Prediction of ChatterStability,Cutting Forces and Surface Location Errors in Flexible MillingSystems”首次提出采用刀具-工件微元模态阵型的方式分别对每一个轴向切深进行判定，进而得出稳定性图。但是文中只考虑了刀具-工件系统的二自由度，并未深入考虑材料去除过程中的工件实时频率响应，即只能得到静态的叶瓣图，不能得到反映加工全过程的稳定性叶瓣图。在实际铣削过程中，机床操作人员主要依赖切削经验数据以及试切来确定加工参数，且材料在去除过程中工件的频率、模态阵型以及刀具切削位置的变化会造成稳定性判断过程中动态参数变化，这些均会对稳定性图造成影响，当忽略以上因素时，会使加工参数选择具有盲目性。尤其针对连续切削过程，材料去除效应引起的工件模态参数变化已不能忽略，且随着材料去除，工件的模态振动会主导整个加工过程。使用传统单一因素的静态稳定性图不能准确预测加工状态，更不能获取工件加工过程中的实时频率响应，这给加工过程引入大量未知因素，尤其是对弱刚性系统的铣削稳定性的预测造成极大困扰，对工件的加工质量造成严重威胁，极易出现加工结果不满足实际要求，甚至工件报废。

发明内容

本发明所要解决的技术问题是提供一种同时考虑刀具和工件多阶模态振动的模型、准确预测工件加工频率响应的模型和工件动力学参数辨识辅助装置，以便确定柔性系统的稳定性加工参数，实时获知材料去除过程中工件模态参数变化的稳定性预测方法。

为解决以上技术问题，本发明采用如下技术方案：

①构建机床直角坐标系X-Y-Z，并建立同时考虑系统多模态振动的四自由度加工测量系统；所述加工测量系统包括工件、力锤、电容传感器、阻尼杆、磁力表座、第一夹具、第二夹具、电荷放大器、PC机和刀具；

所述第一夹具和第二夹具设置在工件下端，所述第一夹具和第二夹具上均设有标尺，所述电容传感器设置在阻尼杆上，所述阻尼杆设置在磁力表座上，所述磁力表座设置在第二夹具上，所述力锤和电容传感器输出端经电荷放大器与PC机连通；

②建立柔性系统切削力模型与加工动力学方程；

③确定工件在材料去除过程中，加工频率响应与模态阵型；

④基于步骤①和步骤②绘制综合三维稳定性图并选取参数；

401沿着刀具加工路径定义刀位点为steps(i＝1,2,…steps)，并沿刀轴方向离散刀具微元，即确定刀触位置为sty(k＝1,2,…sty)；

402将主轴转速离散为stx份，并将相邻刀齿间的迟滞时间离散为m份,并存储至元胞数组 {i,sty,stx,m}；

403用有限元软件分析工件的初始振动模态，获取单元质量矩阵和单元刚度矩阵；

404在第i个刀位点第k个轴向切深处确定要去除的质量单元，将其组装成质量和刚度变化矩阵：ΔM_w和ΔK_w；

405在刀具和工件接触区内，被离散的刀具微元节点会映射在加工工件上，在工件上确定距离刀具映射点处最近的节点，并提取该节点位置处的模态阵型和频率；

406迭代元胞数组{i,sty,stx,m}中的每一个元素，并利用改进的全离散方法确定综合考虑多模态振动和过程阻尼的稳定性叶瓣图。

所述步骤②包括：201将刀具按照轴向切深最大值离散成sty个高度为Δb的切削微元，确定第j齿的相对瞬时未切厚度表达式：

h_j(t)＝f_zsin(φ_j)+[((η_t(t)-η_t(t-τ))-((η_w(t)-η_w(t-τ))] (一)

其中

式(一)中f_z为每齿进给量，t为时间，φ为齿位角，j为刀齿序数，x_t(t),y_t(t)为刀具位移坐标，x_w(t),y_w(t)为工件位移坐标，τ为相邻两刀齿之间的迟滞，进而式(一)演变为：

h_j(t)＝f_zsin(φ_j)+[sin(φ_j),cos(φ_j)]((q_t(t)-q_t(t-τ))-(q_w(t)-q_w(t-τ)))(二)

刀具和工件的多阶振动模态定义为：

202确定考虑刀具后刀面侵入已加工表面的过程阻尼模型为：

其中K_pd为侵入系数,μ为摩擦系数，dV_j为刀具第j^th齿的侵入体积增量；

计算刀具后刀面和工件之间的有效正交后角γ_j为:

式(四)中侵入体积增量为：

其中v_c为刀具的切削速度：

v_t,j和v_w,j分别为刀具和工件第j齿对应的振动速度， l_w为后刀面破损长度，A_j为第j齿侵入三角横截面积；

由式(一)确定刀具和工件的相对振动速度为：

203.综合考虑步骤201-202，构建整体切削力在切向和法向的投影，整体切削力包括剪切项F_ts,j和F_rs,j,犁削项F_tp,j和F_rp,j以及过程阻尼项F_tpd,j和F_rpd,j，具体表达式为：

其中

Τ_j＝[sin(φ_j),cos(φ_j)]，F_t,j和F_r,j分别为作用在刀具的第j齿上切向力和法向力，其中K_t为切向切削力系数，K_r为法向切削力系数，K_te为切向刃口力系数，K_re法向刃口力系数，g(φ_j(t))为选择函数：

其中φ_s为切入角，φ_e为切出角；

204将刀具每个微元对应的力投影到坐标轴X-Y上，并确定铣削力模型的动态项表达式为：

其中

205建立系统动力学方程：

其中M为系统的模态质量，C为系统的模态阻尼，K为系统的模态刚度，F(t)为系统的动态切削力；

将式(十一)转化到模态空间，基于模态坐标变换Q(t)＝UΓ(t)＝[U_t,-U_w]Γ(t)，其中U_t和U_w分别为刀具(10)和工件(1)的模态阵型。

所述步骤③包括：301对于工件的初始状态，确定工件无阻尼齐次方程为：

302设定材料去除的起始点S和当前点C作为判定几何条件，由起始点S和C构成六面体 S1S2S3S4-C1C2C3C4，当工件(1)有限元节点坐标(x_nw,y_nw,z_nw)位于六面体内，即满足：

时，将工件位于坐标范围内的节点保存，并依次叠加构成ΔM_w和ΔK_w；

303随着材料的去除，工件的质量和刚度矩阵会发生变化，工件的无阻尼其次方程变为：

其中ΔM_w和ΔK_w分别为工件质量变化矩阵和刚度变化矩阵；

304利用二范数确定工件上距离刀具微元映射节点(x_tp,y_tp,z_tp)最近的节点，其表达式如下：

min(||(x_nw,y_nw,z_nw)-(x_tp,y_tp,z_tp)||) (十五)

305综合考虑式(十二)和式(十四)，得到工件(1)结构动力修改后的模态阵型为：U_w＝φ_wφ_wm，

其中φ_w和φ_wm分别为式(十二)和式(十四)的模态阵型，则对于在e处激励，在u处响应的频率响应函数H_ue(ω)表示为：

其中m_w为工件被激发的模态数量，r为模态阶数，ω_r为第r阶固有频率，ζ_r为第r阶阻尼比。

所述步骤201中的齿位角将螺旋角引发的滞后效应引入，其表达式为：

其中N为刀齿总数，β为刀具的螺旋角，k为刀具微元序数，d为刀具直径，Ω为主轴转速。

本发明的积极效果如下：本发明构建了同时考虑刀具和工件的多阶振动模型，将刀具的模态振型和工件的模态振型同时考虑到系统的动力学方程中，这有效的避免了单模态状态下预测的失准，提高加工参数选取的可靠性；本发明构建了刀具-工件四自由度铣削模型同时考虑到刀具螺旋角滞后效应，并基于此引入便于计算侵入体积的过程阻尼模型，本发明利用等效三角截面积能快速确定过程阻尼，进一步提升稳定性边界预测的准确性；本发明采用刀具切削微元递增法，逐次增加切削微元的侵入深度，在模态空间坐标下，分别在每个切削节点处判断系统是否稳定，相对于物理坐标，模态坐标下的运算能极大简化系统之间的耦合，降低计算的复杂度；本发明将有限元分析软件Ansys与结构动力学修改方法相结合，利用节点坐标条件筛选出切除材料的节点，进而快速组装成质量和刚度变化矩阵，不需要反复模拟工件某一加工状态下的频率，利用判断条件去除材料也能避免人工组装单元刚度矩阵；本发明将刀具切削路径方向上的模态变化引入到稳定性叶瓣图的绘制中，得到沿刀具路径方向的稳定性叶瓣图，这有效的保证了加工参数选择的准确性；本发明第一夹具和第二夹具上均设有标尺,这能有效保证力锤和电容传感器在系统的坐标中迅速准确定位；本发明磁力表座与阻尼杆相连，阻尼杆可以完成自由伸缩，这就能时电容传感器任意定位在工件表面，能快速测出工件在某一位置的加工响应，便于与预测值比较。

附图说明

图1为本发明加工测量系统示意图；

图2为本发明四自由度铣削模型示意图；

图3为本发明材料去除几何条件判定示意图；

图4为本发明稳定性叶瓣图绘制流程图；

图5为本发明在N1点的频率响应实验值和预测值；

图6为本发明在N2点的频率响应实验值和预测值；

图7为本发明在N3点的频率响应实验值和预测值；

图8为本发明三维多模态综合稳定性叶瓣图；

图9为本发明不同情况下二维稳定性叶瓣图。

在图中：1工件、1-1刀具侵入工件体积、2力锤、3电容传感器、4阻尼杆、5磁力表座、6 标尺、7-1第一夹具7-2第二夹具、8电荷放大器、9PC机、10刀具，A为刀具节点，B为刀具映射在工件上的节点，C为工件上距离刀具节点最近的节点，D为测试频率响应的点的位置，N1第一测试点，N2第二测试点，N3第三测试点，b为轴向切深，a_e为径向切深。

具体实施方式

下面结合附图和具体实例对本发明进行详细说明。

①构建机床直角坐标系X-Y-Z，并建立同时考虑系统多模态振动的四自由度加工测量系统，如图1所示；所述加工测量系统包括工件1、力锤2、电容传感器3、阻尼杆4、磁力表座5、第一夹具7-1、第二夹具7-2、电荷放大器8、PC机9和刀具10；所述第一夹具7-1 和第二夹具7-2设置在工件1下端，所述第一夹具7-1和第二夹具7-2上均设有标尺6，所述电容传感器3设置在阻尼杆4上，所述阻尼杆4设置在磁力表座5上，所述磁力表座5设置在第二夹具7-2上，所述力锤2和电容传感器3输出端经电荷放大器8与PC机9连通；所述阻尼杆4为气动弹簧阻尼杆，阻尼杆4能根据测量需要，在坐标系X-Y-Z的标尺6上任意定位。

②建立柔性系统切削力模型与加工动力学方程；

201如图1和2所示，将刀具10按照轴向切深最大值离散成sty个高度为Δb的切削微元，确定第j齿的相对瞬时未切厚度表达式：

h_j(t)＝f_zsin(φ_j)+[((η_t(t)-η_t(t-τ))-((η_w(t)-η_w(t-τ))] (一)

其中

式(一)中f_z为每齿进给量，t为时间，φ为齿位角，j为刀齿序数，x_t(t),y_t(t)为刀具10 位移坐标，x_w(t),y_w(t)为工件1位移坐标，τ为相邻两刀齿之间的迟滞，进而式(一)演变为：

h_j(t)＝f_zsin(φ_j)+[sin(φ_j),cos(φ_j)]((q_t(t)-q_t(t-τ))-(q_w(t)-q_w(t-τ))) (二)

其中刀具10和工件1的多阶振动模态定义为：

所述齿位角φ将螺旋角引发的滞后效应引入，表达式为：

其中N为刀齿总数，β为刀具10的螺旋角，k为刀具10微元序数，d为刀具10直径，Ω为主轴转速；

202确定考虑刀具10后刀面侵入已加工表面的过程阻尼模型为：

其中K_pd为侵入系数,μ为摩擦系数，dV_j为刀具第j齿的侵入体积增量；

K_pd通过压痕实验确定，其表达式为：

其中ρ为弹塑性变形区域内的距离，ν为泊松比，E为杨氏模量；

计算刀具10后刀面和工件1之间的有效正交后角γ_j为:

式(四)中侵入体积增量为：

其中v_c为刀具10的切削速度：

v_t,j和v_w,j分别为刀具10和工件1第j齿对应的振动速度，l_w为后刀面破损长度，A_j为第j齿侵入三角横截面积；

由式(一)确定刀具10和工件1的相对振动速度为：

203.综合考虑步骤201-202，构建整体切削力在切向和法向的投影，整体切削力包括剪切项F_ts,j和F_rs,j，犁削项F_tp,j和F_rp,j以及过程阻尼项F_tpd,j和F_rpd,j，具体表达式为：

其中

Τ_j＝[sin(φ_j),cos(φ_j)]，F_t,j和F_r,j分别为作用在刀具10 的第j齿上切向力和法向力，其中K_t为切向切削力系数，K_r为法向切削力系数，K_te为切向刃口力系数，K_re法向刃口力系数，g(φ_j(t))为选择函数：

其中φ_s为切入角，φ_e为切出角；当铣削方式为顺铣时，满足

当铣削方式为逆铣时，满足

其中a_e为径向切深；

204将刀具10每个微元对应的力投影到坐标轴X-Y上，并确定铣削力模型的动态项表达式为：

其中

205由步骤204得到每一个刀具10微元节点处的动态铣削力，所有节点构成整体切削力：

且作用在刀具10和工件1上的切削力满足作用力和反作用力，故作用在工件1上的整体切削力满足：F_w,xy＝-F_t,xy；建立系统动力学方程：

其中M为系统的模态质量，C为系统的模态阻尼，K为系统的模态刚度，F(t)为系统的动态切削力，Q(t)为物理坐标；将式(十一)转化到模态空间，基于模态坐标变换 Q(t)＝UΓ(t)＝[U_t,-U_w]Γ(t)，其中U_t和U_w分别为刀具10和工件1的模态阵型，其满足：

和

其中m_t为刀具10被激发的模态数量，m_w为工件1被激发的模态数量；因此式(十一)分别转化为刀具10和工件1在模态坐标下对应的动力学方程：

将此方程利用国际著名期刊《InternationalJournal of Advanced Manufacturing Technology》公开的Third-order updated full-discretization method for milling stability prediction进行求解。

③确定工件1在材料去除过程中，加工频率响应与模态阵型；

301对于工件1的初始状态，确定工件1无阻尼齐次方程为：

302如图3所示，设定材料去除的起始点S和当前点C作为判定几何条件，由起始点S和C 构成六面体S1S2S3S4-C1C2C3C4，当工件1有限元节点坐标(x_nw,y_nw,z_nw)位于六面体内，即满足：

时，将工件1位于坐标范围内的节点保存，并依次叠加构成ΔM_w和ΔK_w；

303随着材料的去除，工件1的质量和刚度矩阵会发生变化，工件1的无阻尼其次方程变为：

其中ΔM_w和ΔK_w分别为工件1质量变化矩阵和刚度变化矩阵；

304利用二范数确定工件1上距离刀具10微元映射节点(x_tp,y_tp,z_tp)最近的节点，其表达式如下：

min(||(x_nw,y_nw,z_nw)-(x_tp,y_tp,z_tp)||) (十五)

305综合考虑式(十二)和式(十四)，得到工件1结构动力修改后的模态阵型为：U_w＝φ_wφ_wm，

其中φ_w和φ_wm分别为式(十二)和式(十四)的模态阵型，即式(十二)和式(十四)特征方程所对应的特征向量；则对于在e处激励在u处响应的频率响应函数H_ue(ω)表示为：

其中m_w为工件1被激发的模态数量，r为模态阶数，ω_r为第r阶固有频率，ζ_r为第r阶阻尼比。

④基于步骤①和步骤②绘制综合三维稳定性图并选取参数，其流程图如图4所示；

401沿着刀具10加工路径定义刀位点为steps(i＝1,2,…steps)，并沿刀轴方向离散刀具10微元，即确定刀触位置为sty(k＝1,2,…sty)；

402将主轴转速离散为stx份，并将相邻刀齿间的迟滞离散为m份,并存储至元胞数组 {i,sty,stx,m}；

403用有限元软件Ansys分析工件1的初始振动模态，获取单元质量矩阵和单元刚度矩阵；

404在第i个刀位点第k个轴向切深处确定要去除的质量单元，将其组装成质量和刚度变化矩阵：ΔM_w和ΔK_w；

405在刀具10和工件1接触区内，被离散的刀具10微元节点会映射在加工工件1上，在工件1上确定距离刀具10映射点处最近的节点，并提取该节点位置处的模态阵型和频率；

406迭代元胞数组{i,sty,stx,m}中的每一个元素，并利用改进的全离散方法确定综合考虑多模态振动和过程阻尼的稳定性叶瓣图。

实施例

如图1所示，本实施例中工件1的材料为Al7075，几何参数为长100mm，高60mm，厚4mm，密度为2770kg/m³，杨氏模量为71GPa，泊松比为0.33，第一夹具7-1和第二夹具 7-2夹持工件1的高度为20mm，刀具10的模态参数如表1所示，刀具10的刀齿数为4，螺旋角为30°，工件1的有限元分析值如表2所示。当刀具10和工件1的模态振型U_t和U_w需要多个节点时，使用插补法获得。

表1

表2

在绘制稳定性叶瓣图之前，通过槽铣来标定切削力系数，每齿进给量从0.04mm开始至0.20mm，增量为0.04mm。得到切向切削系数为K_t＝8.9×10⁸Mpa，法向切削力系数 K_r＝2.4×10⁸Mpa，切向刃口力系数K_te＝2.1×10³N/m，法向刃口力系数 K_re＝3.05×10³N/m，刀具后刀面磨损宽度的测量值为80μm，压痕力系数为 K_sp＝1×10¹³N/m³。

如图1所示，在第二夹具7-2上滑动磁力表座5，确定磁力表座5在标尺6上的坐标，调节阻尼杆4的高度，确定电容传感器3的最终位置，力锤2的敲击位置需要由第一夹具7-1上的标尺6确定，基于此力锤2和确定电容传感器3分别在N1点、N2点、N3点进行Y方向的频率响应测试，由式(十六)得到

和

预测值，即考虑材料去除的加工频率响应。为了获得较大的测试面积，测试时b取20mm，a_e取1mm，预测结果和测试结果如图 5、图6和图7所示，工件在此状态下的测量值和预测值在一个数量级上且基本一致，频率预测值得误差在6％以内，这说明本发明所提出的工件频率响应预测方法具有较高的预测性能，也说明了本模型的准确性。

按照图4流程图的步骤，在Ansys中对工件的模态振型进行提取，在Matlab中进行程序编制，得到考虑多模态效应的综合三维稳定性图，如图8所示。从图8中截取主轴转速为5000rpm时的稳定性叶瓣图，并增加不考虑过程阻尼和不考虑材料去除的叶瓣图，如图9所示；并进行切削实验，实验结果表明本发明所提出的模型能更为准确的预测加工的稳定性。

本实施例中所选用的Ansys软件，通过参数化设计语言用于提取工件1模态振型和固有频率，所用到的主要代码如下:

FILE,file,full,'D:\..\..\FEM analysis\Aluminum\NewAl7075_files\dp0\SYS\MECH'

hbmat,'D:\..\..\FEM analysis\Aluminum\Aluminum',inert,,ascii,mass,no,no

hbmat,'D:\..\..\FEM analysis\Aluminum\Aluminum',stiff,,ascii,stiff,no,yes

/output,'D:\..\..\FEM analysis\Aluminum\ELEMMAT','emat',,

fileaux2,'D:\..\..\FEM analysis\Aluminum\NewAl7075_files\dp0\SYS\MECH\file.emat'

form,long

dump,all

/output

EMATWRITE,yes

本实施例中所选用的Matlab程序代码，主要用于实现判断工件1去除材料节点坐标以及确定工件1上距离刀具10微元最近节点，所用的主要代码如下：

以上所述实施方式仅为本发明的优选实施例，而并非本发明可行实施的穷举。对于本领域一般技术人员而言，在不背离本发明原理和精神的前提下对其所作出的任何显而易见的改动，都应当被认为包含在本发明的权利要求保护范围之内。

Claims (4)

1.一种考虑多模态振动与工件加工响应的稳定性预测方法，其特征在于步骤如下：

①构建机床直角坐标系X-Y-Z，并建立同时考虑系统多模态振动的四自由度加工测量系统；所述加工测量系统包括工件(1)、力锤(2)、电容传感器(3)、阻尼杆(4)、磁力表座(5)、第一夹具(7-1)、第二夹具(7-2)、电荷放大器(8)、PC机(9)和刀具(10)；

所述第一夹具(7-1)和第二夹具(7-2)设置在工件(1)下端，所述第一夹具(7-1)和第二夹具(7-2)上均设有标尺(6)，所述电容传感器(3)设置在阻尼杆(4)上，所述阻尼杆(4)设置在磁力表座(5)上，所述磁力表座(5)设置在第二夹具(7-2)上，所述力锤(2)和电容传感器(3)输出端经电荷放大器(8)与PC机(9)连通；

②建立柔性系统切削力模型与加工动力学方程；

③确定工件(1)在材料去除过程中，加工频率响应与模态阵型；

④基于步骤①和步骤②绘制综合三维稳定性图并选取参数；

401沿着刀具(10)加工路径定义刀位点为steps(i＝1,2,…steps)，并沿刀轴方向离散刀具(10)微元，即确定刀触位置为sty(k＝1,2,…sty)；

402将主轴转速离散为stx份，并将相邻刀齿间的迟滞时间离散为m份,并存储至元胞数组{i,sty,stx,m}；

403用有限元软件分析工件(1)的初始振动模态，获取单元质量矩阵和单元刚度矩阵；

404在第i个刀位点第k个轴向切深处确定要去除的质量单元，将其组装成质量和刚度变化矩阵：ΔM_w和ΔK_w；

405在刀具(10)和工件(1)接触区内，被离散的刀具(10)微元节点会映射在加工工件(1)上，在工件(1)上确定距离刀具(10)映射点处最近的节点，并提取该节点位置处的模态阵型和频率；

406迭代元胞数组{i,sty,stx,m}中的每一个元素，并利用改进的全离散方法确定综合考虑多模态振动和过程阻尼的稳定性叶瓣图。

2.根据权利要求1所述的一种考虑多模态振动与工件加工响应的稳定性预测方法，其特征在于，所述步骤②包括：

201将刀具(10)按照轴向切深最大值离散成sty个高度为Δb的切削微元，确定第j齿的相对瞬时未切厚度表达式：

h_j(t)＝f_zsin(φ_j)+[((η_t(t)-η_t(t-τ))-((η_w(t)-η_w(t-τ))] (一)

其中

式(一)中f_z为每齿进给量，t为时间，φ为齿位角，j为刀齿序数，x_t(t),y_t(t)为刀具(10)位移坐标，x_w(t),y_w(t)为工件(1)位移坐标，τ为相邻两刀齿之间的迟滞，进而式(一)演变为：

h_j(t)＝f_zsin(φ_j)+[sin(φ_j),cos(φ_j)]((q_t(t)-q_t(t-τ))-(q_w(t)-q_w(t-τ))) (二)

其中刀具(10)和工件(1)的多阶振动模态定义为：

202确定考虑刀具(10)后刀面侵入已加工表面的过程阻尼模型为：

其中K_pd为侵入系数,μ为摩擦系数，dV_j为刀具第j^th齿的侵入体积增量；

计算刀具(10)后刀面和工件(1)之间的有效正交后角γ_j为:

式(四)中侵入体积增量为：

其中v_c为刀具(10)的切削速度：

v_t,j和v_w,j分别为刀具(10)和工件(1)第j齿对应的振动速度，l_w为后刀面破损长度，A_j为第j齿侵入三角横截面积；

由式(一)确定刀具(10)和工件(1)的相对振动速度为：

203综合考虑步骤201-202，构建整体切削力在切向和法向的投影，整体切削力包括剪切项F_ts,j和F_rs,j,犁削项F_tp,j和F_rp,j以及过程阻尼项F_tpd,j和F_rpd,j，具体表达式为：

其中

Τ_j＝[sin(φ_j),cos(φ_j)]，F_t,j和F_r,j分别为作用在刀具(10)的第j齿上切向力和法向力，其中K_t为切向切削力系数，K_r为法向切削力系数，K_te为切向刃口力系数，K_re法向刃口力系数，g(φ_j(t))为选择函数：

其中φ_s为切入角，φ_e为切出角；

204将刀具(10)每个微元对应的力投影到坐标轴X-Y上，并确定铣削力模型的动态项表达式为：

其中

205建立系统动力学方程：

其中M为系统的模态质量，C为系统的模态阻尼，K为系统的模态刚度，F(t)为系统的动态切削力；

将式(十一)转化到模态空间，基于模态坐标变换Q(t)＝UΓ(t)＝[U_t,-U_w]Γ(t)，其中U_t和U_w分别为刀具(10)和工件(1)的模态阵型。

3.根据权利要求1所述的一种考虑多模态振动与工件加工响应的稳定性预测方法，其特征在于，所述步骤③包括：

301对于工件(1)的初始状态，确定工件(1)无阻尼齐次方程为：

302设定材料去除的起始点S和当前点C作为判定几何条件，由起始点S和C构成六面体S1S2S3S4-C1C2C3C4，当工件(1)有限元节点坐标(x_nw,y_nw,z_nw)位于六面体内，即满足：

时，将工件(1)位于坐标范围内的节点保存，并依次叠加构成ΔM_w和ΔK_w；

303随着材料的去除，工件(1)的质量和刚度矩阵会发生变化，工件(1)的无阻尼其次方程变为：

其中ΔM_w和ΔK_w分别为工件(1)质量变化矩阵和刚度变化矩阵；

304利用二范数确定工件(1)上距离刀具(10)微元映射节点(x_tp,y_tp,z_tp)最近的节点，其表达式如下：

min(||(x_nw,y_nw,z_nw)-(x_tp,y_tp,z_tp)||) (十五)

305综合考虑式(十二)和式(十四)，得到工件(1)结构动力修改后的模态阵型为：U_w＝φ_wφ_wm，其中φ_w和φ_wm分别为式(十二)和式(十四)的模态阵型，则对于在e处激励，在u处响应的频率响应函数H_ue(ω)表示为：

其中m_w为工件(1)被激发的模态数量，r为模态阶数，ω_r为第r阶固有频率，ζ_r为第r阶阻尼比。

4.根据权利要求2所述的一种考虑多模态振动与工件加工响应的稳定性预测方法，其特征在于：所述步骤201中的齿位角将螺旋角引发的滞后效应引入，其表达式为：

其中N为刀齿总数，β为刀具(10)的螺旋角，k为刀具(10)微元序数，d为刀具(10)直径，Ω为主轴转速。

Images