五轴侧铣可用于加工含有复杂曲面的零件，目前已经在航空工业、汽车、纺织机械和模具行业中得到了广泛应用[1]。在侧铣加工过程中，铣削力可能会引起刀具和工件系统的振动，导致被加工零件表面质量差[2-3]、尺寸精度低[4-5]、零件废品率高等问题。因此，有必要对零件精加工过程中铣削力进行预测。铣削力精确预测能够为抑制工艺系统颤振、提高零件加工质量、优化工艺参数提供可靠的理论依据。

国内外许多学者对铣削力建模进行了大量研究工作，当前的铣削力预测模型依据建模方法可分为微元铣削力模型、多元回归分析模型、有限元模型和人工神经网络模型[6]。其中，微元铣削力模型应用最广泛，学者们基于该模型提出了不同的铣削力计算方法。Aydin等[7]基于微元铣削刃力分布提出了一种侧铣加工铣削力计算方法，该方法仅适用于三轴切削。Li等[8]建立了一种平头铣刀五轴侧铣铣削力模型，该模型考虑了刀具底部铣削刃对铣削力的影响。在五轴侧铣微元铣削力模型中，铣削瞬时未变形切屑厚度(instantaneous undeformed cutting thickness，IUCT)计算受到了重点关注。Li等[9]建立了切削厚度与切削点坐标方程，把切削厚度作为自变量求解一元二次方程获得切削厚度。Sun等[10]建立了刀齿扫掠面方程，采用迭代法求解切削厚度。该方法需要求解扫掠面的Jacobi矩阵，切削厚度计算耗时较长。Zhang等[1]提出了一种解析算法计算IUCT，该方法比Sun等[10]的模型计算效率更高。Guo等[11]基于切削厚度误差补偿提出了五轴侧铣IUCT计算方法。Wang等[12]提出了一种简化的IUCT计算方法，该方法计算耗时比Guo等[11]的模型更短。Ozturk等[13]提出了一种矢量投影方法计算切削厚度，该方法需要频繁变换微元切削刃点坐标来构造矢量。Zhu等[14]基于向量法提出了一种通用铣刀IUCT计算模型，若采用该模型计算平头铣刀的侧铣铣削力，也存在多次变换坐标导致的计算效率较低的问题。

针对平头铣刀五轴侧铣IUCT计算效率低的问题，本文提出了一种基于拆分思想的五轴侧铣加工IUCT计算方法，结合微元铣削力模型，建立了五轴侧铣加工的铣削力预测模型。

刀刃微元坐标系的建立过程如下：以选定的微元切削刃圆周与当前切削刃交点为原点o，过原点o平行于刀具中心线指向刀柄方向为a1轴，过原点o指向微元切削刃所在圆周的圆心方向为a2轴，根据右手准则确定a3轴[15]。

刀具坐标系c-xyz的建立过程如下：圆柱平底铣刀底面圆心为原点c，由铣刀中心线指向刀柄方向为z轴正方向，铣刀进给方向为x轴正方向，根据右手准则确定y轴正方向，如图 1a所示。切削时刀具姿态发生变化，刀具坐标系随之变化，如图 1b所示。

工件坐标系w-x0y0z0是一个固定的坐标系，其x0、y0轴分别平行于机床工作台的长度、宽度方向，由右手准则确定z0轴方向。原点w位置是任意的，对于圆台工件，工件坐标系原点可以选择在圆台下表面中心，如图 1b所示。

螺旋平底立铣刀单齿铣削瞬时受力情况如图 1a所示。将铣刀沿轴向离散为一系列厚度为dz的微元铣削刃，图中所示的微元铣削刃当前所受的铣削力可以分解到微元刀刃坐标系上沿a1、a2、a3 3个方向的力，即轴向分力dFa、径向分力dFr、切向分力dFt。假设图中微元铣削刃为第i个，其轴向坐标为z，φ为该微元坐标系a2轴负方向与刀具坐标系y轴夹角。设θ为刀具底端微元铣削刃的旋转角，则第j个刀齿上第i个微元铣削刃旋转角φ可以表示如下：

其中：R为铣刀半径，β为铣刀螺旋角，n为切削时主轴转速，t为时间，$ \phi_{\mathrm{p}}$为齿距角，N为刀具齿数。

根据Altintas[16]提出的铣削力模型，分别计算微元铣削刃的切向、径向、轴向铣削力为

其中：Ktc为切向作用力系数，Krc为径向作用力系数，Kac为轴向作用力系数，Kte为切向刃口力系数，Kre为径向刃口力系数，Kae为轴向刃口力系数，h(φ)为IUCT，dz为轴向微元铣削刃厚度，g(φ)为微元铣削刃t时刻是否切削材料的判断函数，$ \phi $st为切入角，$ \phi $ ex为切出角，计算示意图如图 2所示。

顺铣时，刀具切入角为

其中ae为径向切深。

顺铣时，刀具切出角为

逆铣时，刀具切入角为

逆铣时，刀具切出角为

微元铣削刃厚度为

其中：ap为轴向切深，Ms为刀具轴向离散的微元数。

将式(4)中微元铣削力转化到刀具坐标系c-xyz下，x、y、z各方向微元铣削力dFx、dFy、dFz如下：

式(11)中坐标变换矩阵T表示如下：

对式(11)中刀具各方向微元铣削力求和，得到刀具坐标系中x、y、z各方向铣削力为

五轴铣削加工时，需要通过坐标变换把刀具坐标系c-xyz中的铣削力转换到工件坐标系中，表示如下[11]：

其中：α(t)是铣削t时刻主轴绕x轴摆角，β(t)是铣削t时刻主轴绕y轴摆角。旋转矩阵A(α(t))、B(β(t))的表达式为

对于节1中h(φ)，即IUCT的计算，若直接采用已有的解析法或向量法计算，则采用微元铣削力模型预测平头铣刀五轴侧铣铣削力的计算效率较低。铣削力计算效率提高的本质是IUCT模型复杂度的降低，因此本文提出了一种复杂度较低的针对平头铣刀五轴侧铣IUCT计算模型，如图 3所示。cj为t时刻刀具位置，zj为t时刻刀轴矢量，zj*为假设刀轴姿态不变时t时刻刀轴矢量，zj-1为t-Δt时刻刀轴矢量，则zj-1与zj*平行，Sj-1为t-Δt时刻刀具扫掠面，γ为刀轴矢量zj-1与zj的夹角，Sj为t时刻刀具扫掠面，Sj*为假设刀具姿态不变时t时刻刀具扫掠面，cj-1为t-Δt时刻刀具位置，其中Δt为单齿铣削时间，其计算公式为

PQ的物理意义为t时刻在刀具坐标系中轴向高度为z的微元铣削刃IUCT。五轴侧铣加工微元铣削刃IUCT根据几何关系被拆分为2部分，第1部分为NQ，NQ可以通过MN求得，MN等效于五轴机床侧铣刀具姿态不变时铣削加工的切削厚度。第2部分为五轴机床侧铣加工t时刻与t-Δt时刻刀轴矢量由于姿态旋转而增加或减少的厚度NP。通过使用该模型计算IUCT，不用迭代求解非线性方程组，减小了计算量，从而提高了五轴侧铣加工铣削力计算效率。

五轴机床刀具姿态不变切削时的厚度MN等效于三轴机床直线切削的切削厚度，MN[16]计算如下：

其中fz为每齿进给量。

如图 3局部放大图所示，由几何关系可得

对刀具轴向坐标为z的铣削刃，时间间隔为Δt的刀轴矢量由于姿态旋转而增加或减少的厚度为

刀具在不同位置tanγ计算方法如下：运用CAM软件提取待加工工件的刀位点信息，提取的信息包括刀具的位置坐标x、y、z及刀具的姿态矢量i、j、k[17]。假设x2、y2、z2为提取当前时刻t刀具位置点坐标，x1、y1、z1为前一个刀具位置点坐标，相邻刀位点间的距离公式为

则相邻刀位点间的坐标插值个数为

插值后获得的刀具在每个刀位点的位置及姿态坐标分别为$x^{\prime}、y^{\prime}、z^{\prime}、i^{\prime}, j^{\prime}、k^{\prime} $。若t时刻刀具姿态矢量为$\boldsymbol{z}_{2}=\left(i_{2}^{\prime}, j_{2}^{\prime}, k_{2}^{\prime}\right), t-\Delta t $时刻刀具姿态矢量为$ \boldsymbol{z}_{1}=\left(i_{1}^{\prime}, j_{1}^{\prime}, k_{1}^{\prime}\right)$，则t与t-Δt时刻刀轴矢量夹角γ的正切值计算公式为

图 3只是五轴侧铣加工刀具姿态变化的一种情况，下面分情况讨论IUCT的计算方法。取刀具坐标系中一轴向矢量$ \boldsymbol{z}_{0}=\left[\begin{array}{lll} 0 & 0 & 1 \end{array}\right]^{\mathrm{T}}$，若t时刻z0在工件坐标系中的z向坐标$z_{w 2} $小于t-Δt时刻刀轴矢量z0在工件坐标系中的z向坐标zw1，此时$ S_{j}^{*}$位于Sj-1与Sj之间(见图 3)，则IUCT为NQ与NP相加；其中，$ z_{w 2}、z_{w 1}$的表达式为

若$ z_{w 2}$与$ z_{w 1}$相等，此时Sj*与Sj重合，IUCT为NQ；否则IUCT为NQ与NP相减。综上所述，五轴侧铣加工IUCT计算公式为

将式(27)代入式(4)，按照节1微元铣削力模型的步骤，即可计算五轴侧铣加工铣削力。

为了验证提出的铣削力模型的准确性，本节开展了五轴侧铣加工实验，实验设备是基于3PRRU并联机构的五轴混联机床，使用Kistler Type 9255C测力仪采集力信号，采样频率为10 kHz。测力仪数据首先传输到电荷放大器，再经过数据采集系统，最后传输到计算机，实验测量平台简图如图 4所示。

毛坯材料为航空铝7075-T7451，待加工工件形状如图 5所示。圆台上底半径r1=50.0 mm，下底半径r2=52.6 mm，高h=30 mm。铣刀选用直径为12 mm的整体硬质合金平底立铣刀，铣刀齿数N=2，铣刀螺旋角β=25°，最大切深apmax=16 mm。主轴转速n=2 000 r/min，每齿进给量fz=0.18 mm，径向切深ae=1 mm，轴向切深ap=10 mm。实验中分3层铣削圆台侧面，铣削方式为逆铣。侧铣加工实验装置如图 6所示。

为了预测铣削力，首先需要标定铣削力系数，采用Altintas[16]提出的平均铣削力法确定铣削力系数，满刀铣槽切削参数及平均铣削力计算结果如表 1所示。根据表 1平均铣削力计算结果，求得铣削力系数为：Ktc=1 141.7 N/mm2，Krc=455.9 N/mm2，Kte=21.3 N/mm2，Kre=21.7 N/mm2。

为了分析本文提出的铣削力预测模型计算效率，将本文铣削力与文[14]所提出的矢量投影法计算切削厚度得到的铣削力进行对比。使用的计算机型号为Dell G3，处理器为Intel(R) Core(TM) i5-8300H @2.30 GHz，已安装的内存为16 GB，操作系统为64位Windows 10专业版，计算软件为MATLAB 2014。为了排除该文提出的模型的计算效率的偶然性，设计了2组仿真铣削参数与矢量投影法进行对比。第一组仿真铣削参数为时间积分步数Ns=300，轴向积分步数Ms=100∶100∶500，计算时长为一个周期，其他参数与节3.1所述一致；第二组仿真铣削参数为轴向积分步数Ms=300，时间积分步数Ns=60∶60∶300，计算时长为一个周期，其他参数与节3.1所述一致。

为了便于比较计算效率，定义计算效率提升比为

其中：tv为矢量投影法计算耗时，tm为本文提出的拆分法计算耗时。对比结果如图 7所示。

根据图 7a和7b可以得出矢量投影法计算耗时跟轴向积分步数大约成线性关系，拆分法计算耗时与轴向积分步数线性度更高，计算效率提升比随轴向积分步数在41%上下波动。根据图 7c和7d可知计算耗时与时间积分步数呈非线性关系，计算效率提升比随时间积分步数波动，随着时间积分步数增加，效率提升比有增大趋势。综合来看，拆分法计算效率较高，效率提升比在40%左右。

为了验证提出模型的准确性，取轴向积分步数为300，时间积分步数为300，将铣削力预测结果与第2层侧铣加工实验结果进行比较，对比结果如图 8和9所示。

根据图 8可知，五轴侧铣过程中x向与y向铣削力预测的变化趋势与实际测量的铣削力变化趋势一致，铣削力峰值轮廓线零点、拐点出现的时间基本吻合。对比预测值与实测值的峰值，可以看出预测值比实测值偏大。原因是实际加工过程中，刀具受力变形，刀具工件接触区域减小，进而实际铣削力减小[12]。根据图 9可知，五轴侧铣过程中x向、y向铣削力不是周期性变化，相邻切削周期铣削力峰值有波动，铣削力预测结果与实验结果基本一致，验证了本文提出的铣削力模型具有良好的准确性。

通过对五轴侧铣加工刀具姿态变化的研究和分析，本文提出了一种五轴侧铣加工IUCT模型。基于微元铣削力模型，建立了五轴侧铣加工铣削力预测模型。与其他铣削力计算方法相比，本文提出的铣削力预测模型计算效率较高。通过五轴侧铣加工实验，验证了提出的铣削力预测模型的准确性。结果表明：该模型可用于预测五轴侧铣加工的铣削力，对确定加工工艺参数、分析零件加工过程中刀具工件受力情况具有一定的参考价值。