CN105058166A - 基于铣刀精确建模的刀尖点频响函数预测方法 - Google Patents

Abstract

本发明公开了基于铣刀精确建模的刀尖点频响函数预测方法，该方法将机床-主轴-刀柄-铣刀划分为机床-主轴-刀柄-部分刀杆、剩余刀杆、铣刀过渡段以及多段刀齿；考虑铣刀的实际建模，根据实际情况将铣刀分为对称型铣刀和非对称型铣刀，以2刃铣刀和4刃铣刀为例分别详细计算了铣刀各部分的截面惯性矩，并提出一种铣刀过渡段的通用建模方法，提高铣刀刀尖点频响函数的预测精度。解决铣刀精确建模，准确预测刀尖点频响函数的问题。

Description

基于铣刀精确建模的刀尖点频响函数预测方法

技术领域

本发明属于振动测试技术领域，具体涉及基于铣刀精确建模的刀尖点频响函数预测方法。

背景技术

机床的切削颤振会降低被加工工件的表面质量并影响机床的使用寿命。目前避免切削颤振最有效的方法是通过加工***的稳定性叶瓣图选取稳定的切削加工参数。稳定性叶瓣图由加工***中刀尖点的频响函数计算得到，刀尖点频响函数可通过实验模态测试确定，但在实际加工中经常更换铣刀导致的重复模态测试不仅耗时，而且会引入人为误差。针对这个问题，不少文献或者专利从结合面建模及参数识别的角度提出了铣刀刀尖点频响函数的预测方法，但是铣刀模型的准确建立是影响整机频响函数预测的重要方面。

从目前的研究来看，关于刀尖点频响函数的预测方法依然存在以下问题。第一，大部分文献虽然考虑了铣刀的建模，但是依然采用不同的等效方法，具有局限性且在没有做大量的实验的前提下，无法将其推广；第二，部分考虑铣刀实际建模的文献或专利，或者讨论某一类铣刀的建模，或者简化铣刀建模，没有考虑铣刀过渡段的精确建模。事实上，从铣刀截面惯性矩的计算角度来看，铣刀可以分为对称型铣刀和非对称型铣刀；而铣刀过渡段的简化建模也不符合实际。

发明内容

针对现有技术的不足，本发明专利的目的在于提供基于铣刀精确建模的刀尖点频响函数预测方法，解决铣刀精确建模，准确预测刀尖点频响函数的问题。

本发明专利解决的技术问题可以采用以下技术方案来实现：基于铣刀精确建模的刀尖点频响函数预测方法，包括以下步骤：

步骤1：将铣刀铣刀一端安装到刀柄中，将刀柄带有锥度的一端安装在与刀柄尾部锥度相配合的机床-主轴中，从而组成机床-主轴-刀柄-铣刀***；

步骤2：将机床-主轴-刀柄-铣刀***划分为机床-主轴-刀柄-部分刀杆、剩余刀杆、铣刀过渡段、刀齿四个部分，并将刀齿划分为多段，其中子结构I-IV为刀齿部分，子结构V为铣刀过渡段部分，子结构VI为剩余刀杆部分，子结构VII为机床-主轴-刀柄-部分刀杆部分；

步骤3：利用RCSA方法耦合子结构I-VII获得刀尖点频响函数矩阵；

步骤4：计算R_7b7b，R_7b7b是子结构VII末端的频响函数矩阵；

步骤5：采用Timoshenko梁模型得到单个子结构I-VI在自由状态下的频响函数矩阵R_ij的表达式；

步骤6：考虑铣刀的精确建模，计算子结构I-VI的截面惯性矩，结合步骤5，得到子结构I-VI的频响函数矩阵；

步骤7：将步骤6中计算得到的子结构I-VI的频响函数矩阵代入步骤3中的RCSA耦合公式，得到RC_7a7a、RC_7a1、RC_17a、RC₁₁，然后将RC_7a7a、RC_7a1、RC_17a、RC₁₁和R_7b7b代入步骤3中的刀尖点频响函数矩阵计算公式，最终得到刀尖点频响函数矩阵G₁₁，即可获得刀尖点频响函数HC₁₁；

步骤8：对铣刀进行模态锤击试验获得铣刀在x、y向上的实测刀尖点频响函数，并对比步骤7中预测得到的铣刀刀尖点频响函数。

所述的预测方法，其中：

步骤3具体为：利用RCSA方法耦合子结构I-VI获得刀尖点频响函数矩阵，两端自由的子结构的频响函数矩阵表示为如下通式：

$R_{ij}=\left[\begin{array}{cc}{H}_{ij}& L_{ij}\\ N_{ij}& P_{ij}\end{array}\right]=\left[\begin{array}{cc}{X}_i/F_j& X_i/M_j\\ {\mathrm{\θ}}_i/F_j& {\mathrm{\θ}}_i/M_j\end{array}\right]$

式中，i，j＝1，2a或i，j＝2b，3a或i，j＝3b，4a或i，j＝4b，5a或i，j＝5b，6a或i，j＝6b，7a；

H_ij、L_ij、N_ij和P_ij依次表示单个两端自由的子结构为在j点激励，在i点获得响应的位移/力、位移/力矩、转角/力以及转角/力矩的频响函数。F_j、M_j依次表示为j点处所受到的外力和外力矩；X_i、Θ_i依次表示为在外力和外力矩的作用下，i点处的平动位移和转角；1，2a，2b，3a，3b，4a，4b，5a，5b，6a，6b，7a分别表示子结构I-VII从铣刀刀尖到机床的两端端点；

其中，自由状态下相邻两个子结构耦合后的频响函数矩阵，可表示为以下通式：

${\mathrm{RC}}_{ij}=\left[\begin{array}{cc}{\mathrm{HC}}_{ij}& {\mathrm{LC}}_{ij}\\ {\mathrm{NC}}_{ij}& {\mathrm{PC}}_{ij}\end{array}\right]$

式中，RC_ij表示自由状态下相邻两个子结构耦合后的频响函数矩阵，i，j＝1，3a，或i，j＝1，4a，或i，j＝1，5a，或i，j＝1，6a，或i，j＝1，7a，，HC_ij、LC_ij、NC_ij和PC_ij依次表示子结构I-II在j点激励，在i点获得响应的位移/力、位移/力矩、转角/力以及转角/力矩的频响函数。

根据自由状态下相邻两个子结构在结合处的兼容条件和平衡方程，最终计算得到自由状态下相邻两个子结构的频响函数矩阵与自由状态下相邻两个子结构频响函数矩阵之间的关系，以子结构I和子结构II为例说明。

$\left\{\begin{array}{c}{\mathrm{RC}}_{11}=R_{11}-R_{12a}{(R_{2a2a}+R_{2b2b})}^{-1}{R}_{2a1}\\ {\mathrm{RC}}_{13a}=R_{12a}{(R_{2a2a}+R_{2b2b})}^{-1}{R}_{2b3a}\\ {\mathrm{RC}}_{3a1}=R_{3a2b}{(R_{2a2a}+R_{2b2b})}^{-1}{R}_{2a1}\\ {\mathrm{RC}}_{3a3a}=R_{3a3a}-R_{3a2b}{(R_{2a2a}+R_{2b2b})}^{-1}{R}_{2b3a}\end{array}\right.$

同理可通过上式顺序耦合子结构I-VII，最终得到刀尖点频响函数矩阵G₁₁：

G₁₁＝RC₁₁-RC_17a(RC_7a7a+R_7b7b)^-1RC_7a1

式中，G₁₁还表示为：

$G_{11}=\left[\begin{array}{cc}{\mathrm{HC}}_{11}& {\mathrm{LC}}_{11}\\ {\mathrm{NC}}_{11}& {\mathrm{PC}}_{11}\end{array}\right]$

式中，HC₁₁、LC₁₁、NC₁₁和PC₁₁依次为刀尖点的位移/力、位移/力矩、转角/力以及转角/力矩频响函数，其中HC₁₁简称为刀尖点频响函数。

所述的预测方法，其中：

步骤4具体为：计算R_7b7b，R_7b7b是子结构VII末端的频响函数矩阵，表示为：

$R_{7b7b}=\left[\begin{array}{cc}{H}_{7b7b}& L_{7b7b}\\ N_{7b7b}& P_{7b7b}\end{array}\right]$

式中，H_7b7b通过对子结构VII进行锤击试验得到，并采用一阶有限差分法计算R_7b7b中其余参数L_7b7b，N_7b7b，P_7b7b：

$R_{7b7b}=\left[\begin{array}{cc}{H}_{7b7b}& \frac{H_{77}-H_{78}}{S}\\ \frac{H_{77}-H_{78}}{S}& \frac{{N_{7b7b}}^2}{H_{7b7b}}\end{array}\right]$

式中，S表示刀杆上接近刀柄末端的点7和刀杆上接近刀刃的点8之间的距离，其中点7是刀杆上靠近刀柄末端的点，点8是刀杆上靠近刀齿部分的点，H₇₇表示在点7的原点平动位移频响函数，H₇₈表示在点8的跨点平动位移频响函数，H₇₇、H₇₈通过锤击试验得到，即在点7和点8处分别进行模态锤击实验，使用力锤分别对点7和点8施加激振力，并用安装在铣刀刀杆处与点7正对的加速度传感器记录振动加速度响应，利用数据采集***同步采集力信号和加速度信号，经过多次测量，使用计算机中的ModelVIEW软件对信号进行处理，得到频响函数H₇₇和H₇₈。

所述的预测方法，其中：

步骤5具体为：采用Timoshenko梁模型计算单个子结构I-VI在自由状态下的频响函数矩阵R_ij，其中

$H_{ij}=\underset{r=0}{\overset{\mathrm{\∞}}{\Σ}}\frac{{\mathrm{\Φ}}_r\left(x_i\right){\mathrm{\Φ}}_r\left(x_j\right)}{(1+i\mathrm{\γ}){\mathrm{\ω}}_r^2-{\mathrm{\ω}}^2}$

$N_{ij}=\underset{r=0}{\overset{\mathrm{\∞}}{\Σ}}\frac{{\mathrm{\Φ}}_r^{\′}\left(x_i\right){\mathrm{\Φ}}_r\left(x_j\right)}{(1+i\mathrm{\γ}){\mathrm{\ω}}_r^2-{\mathrm{\ω}}^2}$

$L_{ij}=\underset{r=0}{\overset{\mathrm{\∞}}{\Σ}}\frac{{\mathrm{\Φ}}_r\left(x_i\right){\mathrm{\Φ}}_r^{\′}\left(x_j\right)}{(1+i\mathrm{\γ}){\mathrm{\ω}}_r^2-{\mathrm{\ω}}^2}$

$P_{ij}=\underset{r=0}{\overset{\mathrm{\∞}}{\Σ}}\frac{{\mathrm{\Φ}}_r^{\′}\left(x_i\right){\mathrm{\Φ}}_r^{\′}\left(x_j\right)}{(1+i\mathrm{\γ}){\mathrm{\ω}}_r^2-{\mathrm{\ω}}^2}$

式中， ${\mathrm{\Φ}}_r\left(x\right)=A_r\left(C_{1}sin\right(\frac{{\mathrm{\α}}_r}{L}x)+C_{2}cos(\frac{{\mathrm{\α}}_r}{L}x)+C_3\sinh (\frac{{\mathrm{\β}}_r}{L}x)+C_4\cosh (\frac{{\mathrm{\β}}_r}{L}x\left)\right)$

$C_1=L,C_2=-\frac{(\mathrm{\α}-\mathrm{\λ})(cos\mathrm{\α}-\cosh \mathrm{\β})}{(\mathrm{\λ}-\mathrm{\α})sin\mathrm{\α}+\frac{\mathrm{\λ}\mathrm{\α}}{\mathrm{\β}\mathrm{\δ}}(\mathrm{\β}-\mathrm{\δ})sinh\mathrm{\β}}L,C_3=\frac{{\mathrm{\α}}_r-{\mathrm{\λ}}_r}{{\mathrm{\δ}}_r-{\mathrm{\β}}_r}L,$

$C_4=-\frac{{\mathrm{\λ}}_r{\mathrm{\α}}_r}{{\mathrm{\β}}_r{\mathrm{\δ}}_r}\frac{(\mathrm{\α}-\mathrm{\λ})(cos\mathrm{\α}-\cosh \mathrm{\β})}{(\mathrm{\λ}-\mathrm{\α})sin\mathrm{\α}+\frac{\mathrm{\λ}\mathrm{\α}}{\mathrm{\β}\mathrm{\δ}}(\mathrm{\β}-\mathrm{\δ})sinh\mathrm{\β}}L$

$\mathrm{\α}=\sqrt{\mathrm{\Ω}+\mathrm{\ϵ}},\mathrm{\β}=\sqrt{-\mathrm{\Ω}+\mathrm{\ϵ}},\mathrm{\Ω}=\frac{b^2(s^2+R^2)}{2}$

$\mathrm{\ϵ}=b\sqrt{\frac{1}{4}{b}^2{(s^2+R^2)}^2-(b^2{s}^2{R}^2-1)},b^2=\frac{{\mathrm{\ρA\ω}}^2{L}^4}{EI}$

$s^2=\frac{EI}{k^{\′}{\mathrm{AGL}}^2},R^2=\frac{1}{{\mathrm{AL}}^2},\mathrm{\λ}=\mathrm{\α}-\frac{b^2{s}^2}{\mathrm{\α}},k^{\′}=\frac{6(1+\mathrm{\υ})}{7+6\mathrm{\υ}}$

根据振型函数的正交性可得

${\∫}_{x=0}^L{\left\{{\mathrm{\Φ}}_r\left(x\right)\right\}}^T\mathrm{\ρ}A\left\{{\mathrm{\Φ}}_r\left(x\right)\right\}dx=\left\{\begin{array}{c}1,r=s\\ 0,r\≠s\end{array}\right.$

令r，s＝1,2,3…n可得A_r(r＝1,2,3…n)

式中，Φ,Φ'分别表示两端自由的Timoshenko梁的平动振型函数和平动振型函数对x的导数；α_r,β_r均是是无量纲频率参数，A_r是振型函数的模态系数，γ表示阻尼因子，ω表示角频率，r表示模态阶数；ρ、E和G依次表示梁材料的密度、杨氏模量和剪切模量，A、I、L、k’和υ依次表示梁截面的面积、截面惯性矩、长度、剪切系数和梁材料的泊松比。

所述的预测方法，其中：

步骤6具体为：考虑铣刀的精确建模，计算子结构I-VI的截面惯性矩，结合步骤5，得到子结构I-VI的频响函数矩阵：

根据铣刀截面惯性矩的计算将铣刀分为对称型铣刀和非对称型铣刀，即在铣刀截面建立x-y笛卡尔二维坐标系，如果铣刀刀齿截面对x轴的惯性矩和对y轴的惯性矩不等，则该铣刀属于非对称型铣刀，否则属于对称型铣刀；以下计算以2刃铣刀和4刃铣刀各子结构截面惯性矩，其中2刃铣刀属于非对称型铣刀，4刃铣刀属于对称型铣刀；①铣刀刀齿子结构I-IV截面惯性矩的计算

2刃和4刃铣刀截面具有中心对称的特点，先求解1/2的2刃铣刀和1/4的4刃铣刀端面截面惯性矩的计算公式，在铣刀端面处建立2刃铣刀和4刃铣刀刀齿截面的x-y笛卡尔二维坐标系并计算铣刀刀齿的截面惯性矩：

$\left\{\begin{array}{c}{I}_{x,2}={\∫}_0^{\mathrm{\π}}{\∫}_0^{R_{eq2}}{\mathrm{\ρ}}^3{\sin }^2\mathrm{\θ}d\mathrm{\ρ}d\mathrm{\θ}-\frac{{{\mathrm{\πf}}_d}^4}{128}\\ I_{y,2}={\∫}_0^{\mathrm{\π}}{\∫}_0^{R_{eq2}}{\mathrm{\ρ}}^3{\cos }^2\mathrm{\θ}d\mathrm{\ρ}d\mathrm{\θ}-\[\frac{{{\mathrm{\πf}}_d}^4}{128}+\frac{{{\mathrm{\πf}}_d}^2}{8}{\left(\frac{D-f_d}{2}\right)}^2\]\end{array}\right.$

$\left\{\begin{array}{c}{I}_{x,4}={\∫}_0^{\mathrm{\π}/2}{\∫}_0^{R_{eq4}}{\mathrm{\ρ}}^3{\sin }^2\mathrm{\θ}d\mathrm{\ρ}d\mathrm{\θ}-\[\frac{{{\mathrm{\πf}}_d}^4}{128}+\frac{{{\mathrm{\πf}}_d}^2}{8}{\left(\frac{D-f_d}{2}\right)}^2\]\\ I_{y,4}={\∫}_0^{\mathrm{\π}/2}{\∫}_0^{R_{eq4}}{\mathrm{\ρ}}^3{\cos }^2\mathrm{\θ}d\mathrm{\ρ}d\mathrm{\θ}-\frac{{{\mathrm{\πf}}_d}^4}{128}\end{array}\right.$

$\begin{array}{cc}{R}_{eq2}=-acos\mathrm{\θ}+\sqrt{(r^2-a^2)+a^2{\cos }^2\mathrm{\θ}}& 0\≤\mathrm{\θ}\≤\mathrm{\π}\end{array}$

$\begin{array}{cc}{R}_{eq4}=asin\mathrm{\θ}+\sqrt{r^2-a^2{\cos }^2\mathrm{\θ}}& 0\≤\mathrm{\θ}\≤\mathrm{\π}/2\end{array}$

式中，r表示铣刀切削刃圆弧半径，a表示铣刀切削刃圆弧中心到铣刀中心的距离，f_d表示铣刀螺旋槽直径，R_eqi表示i刃铣刀刀齿截面的等效半径，i＝2,4，I_x,j表示j刃铣刀刀齿截面对x轴的惯性矩，I_y,j表示j刃铣刀刀齿截面对y轴的惯性矩，j＝2,4；

根据上式将I_x,2、I_y,2展开，由于I_x,2、I_y,2的展开式较为复杂，因此通过中间量I₂₁、I₂₂、I₂₃和矩阵A、B、C表示：

$\left\{\begin{array}{c}{I}_{x,2}=\frac{AB}{4}-\frac{{\mathrm{\πf}}_d^4}{128}\\ I_{y,2}=\frac{AC}{4}-\[\frac{{\mathrm{\πf}}_d^4}{128}+\frac{{\mathrm{\πf}}_d^2}{8}{\left(\frac{D-f_d}{2}\right)}^2\]\end{array}\right.$

$\left\{\begin{array}{c}{I}_{21}=\frac{2r^4}{a^3}(\frac{t}{8}-\frac{\cos t\sin t}{8}+\frac{\cos {\mathrm{tsin}}^{3}t}{4}){|}_0^{\arcsin \frac{a}{r}}\\ I_{22}=\frac{4r^6}{a^5}(\frac{t}{16}-\frac{\sin 2t}{24}+\frac{\sin 4t}{192}+\frac{\cos t{\sin }^{5}t}{6}){|}_0^{\arcsin \frac{a}{r}}\\ I_{23}=\frac{2r^2}{a}(\frac{t}{2}+\frac{\sin }{4}){|}_0^{\arcsin \frac{a}{r}}\end{array}\right.$

$\left\{\begin{array}{c}A=\[r^4,a^4,2r^2{a}^2,4a^2{r}^2,-4a^4,-4{\mathrm{ar}}^2,-4a^3\]\\ B={\[\frac{\mathrm{\π}}{2},\frac{\mathrm{\π}}{4},\frac{-\mathrm{\π}}{4},\frac{\mathrm{\π}}{8},\frac{\mathrm{\π}}{16},I_{21},I_{21}-I_{22}\]}^T\\ C={\[\frac{\mathrm{\π}}{2},\frac{\mathrm{\π}}{4},\frac{\mathrm{\π}}{4},\frac{3\mathrm{\π}}{8}\frac{\mathrm{\π}}{16},I_{23}-I_{21},I_{23}-3I_{21}+I_{22}\]}^T\end{array}\right.$

同理将I_x,4、I_y,4展开，由于I_x,4、I_y,4的展开式较为复杂，因此通过中间量I₄₁、I₄₂、I₄₃和矩阵D、E、F表示：

$\left\{\begin{array}{c}{I}_{x,4}=\frac{DE}{4}-\[\frac{{\mathrm{\πf}}_d^4}{128}+\frac{{\mathrm{\πf}}_d^2}{8}{\left(\frac{D-f_d}{2}\right)}^2\]\\ I_{y,4}=\frac{DF}{4}-\frac{{\mathrm{\πf}}_d^4}{128}\end{array}\right.$

$\left\{\begin{array}{c}{I}_{41}=\frac{r^2}{a}(\frac{t}{2}+\frac{\sin 2t}{4}){|}_0^{\arcsin \frac{a}{r}}\\ I_{42}=\frac{r^4}{a^3}(\frac{t}{8}-\frac{\cos t\sin t}{8}+\frac{\cos {\mathrm{tsin}}^{3}t}{4}){|}_0^{\arcsin \frac{a}{r}}\\ I_{43}=\frac{2r^6}{a^5}(\frac{t}{16}-\frac{\sin 2t}{24}+\frac{\sin 4t}{192}+\frac{\cos t{\sin }^{5}t}{6}){|}_0^{\arcsin \frac{a}{r}}\end{array}\right.$

$\left\{\begin{array}{c}D=\[r^4,a^4,-2r^2{a}^2,4a^2{r}^2,-4a^4,4{\mathrm{ar}}^2,-4a^3\]\\ E={\[\frac{\mathrm{\π}}{4},\frac{\mathrm{\π}}{8},\frac{-\mathrm{\π}}{8},\frac{3\mathrm{\π}}{16},\frac{\mathrm{\π}}{32},I_{41}-I_{42},-2I_{41}+4I_{42}-I_{43}\]}^T\\ F={\[\frac{\mathrm{\π}}{4},\frac{\mathrm{\π}}{8},\frac{\mathrm{\π}}{8},\frac{\mathrm{\π}}{8},\frac{\mathrm{\π}}{32},I_{42},I_{43}-I_{42}\]}^T\end{array}\right.$

式中通过铣刀截面的几何边界条件确定r，a，对于2刃铣刀通过求解得到：

$\left\{\begin{array}{c}r=(D-f_d)/2\\ a=f_d/2\end{array}\right.$

同理可得4刃铣刀刀齿的参数r，a：

$\left\{\begin{array}{c}r=\frac{D^2-2{\mathrm{Df}}_d+2{f_d}^2}{2D}\\ a=\frac{{\mathrm{Df}}_d-{f_d}^2}{D}\end{array}\right.$

然后根据中心对称可得到2刃铣刀和4刃铣刀刀齿端面剩余部分的截面惯性矩，经过叠加最终得到2刃和4刃铣刀刀齿端面总的截面惯性矩：

$\begin{array}{cc}\left\{\begin{array}{c}{I}_{x,2tot}=2I_{x,2}\\ I_{y,2tot}=2I_{y,2}\end{array}\right.& \left\{\begin{array}{c}{I}_{y,4tot}=2(I_{y,4}+I_{x,4})\\ I_{x,4tot}=I_{y,4tot}\end{array}\right.\end{array}$

②计算距离铣刀刀齿端面不同距离处的截面惯性矩

距离铣刀刀齿端面不同距离处的截面和刀齿端面的形状相同，不同的是其在位置上相对刀齿端面转动了一个滞后角ψ，为便于方便，在距离铣刀刀齿端面不同距离处的截面建立局部坐标系u-v，在刀齿端面建立全局坐标系x-y，并将局部坐标系u-v上的截面惯性矩转换到全局坐标系x-y上，转换公式如下所示：

$\left\{\begin{array}{c}{I}_{xtot}={\cos }^2{\mathrm{\ψI}}_u+{\sin }^2{\mathrm{\ψI}}_v-I_{uv}\sin 2\mathrm{\ψ}\\ I_{ytot}={\cos }^2{\mathrm{\ψI}}_u+{\sin }^2{\mathrm{\ψI}}_v-I_{uv}\sin 2\mathrm{\ψ}\end{array}\right.$

$\mathrm{\ψ}=\frac{2ztan\mathrm{\β}}{D}$

式中，z表示刀齿任意截面距离铣刀刀齿端面的距离；β表示铣刀螺旋角；D表示铣刀刀齿圆柱基体直径。对于2刃和4刃铣刀，I_uv＝0，将①中的I_x,2tot和I_y,2tot、I_x,4tot和I_y,4tot替换I_u和I_v，可分别得到2刃铣刀和4刃铣刀刀齿部分距离铣刀刀齿端面不同距离处的截面惯性矩；

③铣刀过渡段子结构截面惯性矩的计算

铣刀过渡段截面惯性矩的计算采用分块计算，然后叠加的方法。根据铣刀截面的特点，提出了一种通用的划分方法，即在铣刀过渡段截面建立以其形心为坐标原点的笛卡尔二维坐标系，然后以其原点为中心，划分出若干块三角形区域和扇形区域，通过分别计算每块区域的截面惯性矩，利用坐标转换公式，经过叠加，最后得到整体的截面惯性矩。

I_tx表示三角形截面对x轴的惯性矩；I_ty表示三角形截面对y轴的惯性矩；I_sx表示扇形截面对x轴的惯性矩；I_sy表示扇形截面对y轴的惯性矩；

根据坐标转化，经过叠加，得到2刃铣刀过渡段总的截面惯性矩的计算公式：

$\left\{\begin{array}{c}{I}_{x-2tot}=2(I_{tx}+I_{sx})\\ I_{y-2tot}=2(I_{ty}+I_{sy})\end{array}\right.$

同理可得4刃铣刀过渡段总的截面惯性矩的计算公式：

$\left\{\begin{array}{c}{I}_{x-4tot}=2(I_{tx}+I_{sx}+I_{ty}+I_{sy})\\ I_{y-4tot}=I_{x-4tot}\end{array}\right.$

④铣刀刀杆子结构截面惯性矩的计算

铣刀刀杆部分是一段圆柱段，其截面为圆形截面，故其截面惯性矩的计算公式为

$I_x=I_y=\frac{{\mathrm{\πD}}^2}{64}$

将①②③④带入步骤5中的频响函数计算公式H_ij、L_ij、N_ij和P_ij，可得到子结构I-VI自由状态下的频响函数矩阵。

附图说明

图1为本发明专利的机床-主轴-刀柄-铣刀***的子结构划分；

图2为本发明专利的RCSA原理图；

图3为本发明专利的H₇₇和H₇₈锤击试验测量的实验原理图；

图4为本发明专利的2刃铣刀截面；

图5为本发明专利的4刃铣刀截面；

图6为本发明专利的2刃铣刀刀齿端面(截面)坐标系；

图7为本发明专利的4刃铣刀刀齿端面(截面)坐标系；

图8为本发明专利的2刃铣刀过渡段截面；

图9为本发明专利的4刃铣刀过渡段截面；

图10为本发明专利的铣刀过渡段截面的三角形区域；

图11为本发明专利的铣刀过渡段截面的扇形区域；

图12为本发明专利的铣刀刀尖点频响特性测试图；

图13为本发明专利的2刃和4刃铣刀的x向、y向实测刀尖点频响函数对比；

图14为本发明专利的5把铣刀预测刀尖点频响函数和实测刀尖点频响函数的对比。

其中1—机床-主轴，2—刀柄，3—铣刀，4—加速度传感器，5—力锤。

具体实施方式

下面结合附图详细说明本发明专利的优选实施方式。

本发明专利的具体实施方式如下：

本发明专利的目的在于提供基于铣刀精确建模的刀尖点频响函数预测方法，解决铣刀精确建模，准确预测刀尖点频响函数的问题。具体的步骤如下。

步骤1：如图1所示，将铣刀一端安装到刀柄中，将刀柄带有锥度的一端安装在与刀柄尾部锥度相配合的机床-主轴中，从而组成机床-主轴-刀柄-铣刀***。

步骤2：如图1所示，将机床-主轴-刀柄-铣刀***划分为机床-主轴-刀柄-部分刀杆、剩余刀杆、铣刀过渡段、刀齿四个部分，并将刀齿划分为多段，其中子结构I-IV为刀齿部分，子结构V为铣刀过渡段部分，子结构VI为剩余刀杆部分，子结构VII为机床-主轴-刀柄-部分刀杆部分。

步骤3：如图2所示，利用RCSA方法耦合子结构I-VI获得刀尖点频响函数矩阵，两端自由的子结构的频响函数矩阵表示为如下通式：

$R_{ij}=\left[\begin{array}{cc}{H}_{ij}& L_{ij}\\ N_{ij}& P_{ij}\end{array}\right]=\left[\begin{array}{cc}{X}_i/F_j& X_i/M_j\\ {\mathrm{\θ}}_i/F_j& {\mathrm{\θ}}_i/M_j\end{array}\right]$

式中，i，j＝1，2a或i，j＝2b，3a或i，j＝3b，4a或i，j＝4b，5a或i，j＝5b，6a或i，j＝6b，7a；

H_ij、L_ij、N_ij和P_ij依次表示单个两端自由的子结构为在j点激励，在i点获得响应的位移/力、位移/力矩、转角/力以及转角/力矩的频响函数。F_j、M_j依次表示为j点处所受到的外力和外力矩；X_i、Θ_i依次表示为在外力和外力矩的作用下，i点处的平动位移和转角；1，2a，2b，3a，3b，4a，4b，5a，5b，6a，6b，7a分别表示子结构I-VII从铣刀刀尖到机床的两端端点；

其中，自由状态下相邻两个子结构耦合后的频响函数矩阵，可表示为以下通式：

${\mathrm{RC}}_{ij}=\left[\begin{array}{cc}{\mathrm{HC}}_{ij}& {\mathrm{LC}}_{ij}\\ {\mathrm{NC}}_{ij}& {\mathrm{PC}}_{ij}\end{array}\right]$

式中，RC_ij表示自由状态下相邻两个子结构耦合后的频响函数矩阵，i，j＝1，3a，或i，j＝1，4a，或i，j＝1，5a，或i，j＝1，6a，或i，j＝1，7a，，HC_ij、LC_ij、NC_ij和PC_ij依次表示子结构I-II在j点激励，在i点获得响应的位移/力、位移/力矩、转角/力以及转角/力矩的频响函数。

根据自由状态下相邻两个子结构在结合处的兼容条件和平衡方程，最终计算得到自由状态下相邻两个子结构的频响函数矩阵与自由状态下相邻两个子结构频响函数矩阵之间的关系，以子结构I和子结构II为例说明。

$\left\{\begin{array}{c}{\mathrm{RC}}_{11}=R_{11}-R_{12a}{(R_{2a2a}+R_{2b2b})}^{-1}{R}_{2a1}\\ {\mathrm{RC}}_{13a}=R_{12a}{(R_{2a2a}+R_{2b2b})}^{-1}{R}_{2b3a}\\ {\mathrm{RC}}_{3a1}=R_{3a2b}{(R_{2a2a}+R_{2b2b})}^{-1}{R}_{2a1}\\ {\mathrm{RC}}_{3a3a}=R_{3a3a}-R_{3a2b}{(R_{2a2a}+R_{2b2b})}^{-1}{R}_{2b3a}\end{array}\right.$

同理可通过上式顺序耦合子结构I-VII，最终得到刀尖点频响函数矩阵G₁₁：

G₁₁＝RC₁₁-RC_17a(RC_7a7a+R_7b7b)^-1RC_7a1

式中，G₁₁还表示为：

$G_{11}=\left[\begin{array}{cc}{\mathrm{HC}}_{11}& {\mathrm{LC}}_{11}\\ {\mathrm{NC}}_{11}& {\mathrm{PC}}_{11}\end{array}\right]$

式中，HC₁₁、LC₁₁、NC₁₁和PC₁₁依次为刀尖点的位移/力、位移/力矩、转角/力以及转角/力矩频响函数，其中HC₁₁简称为刀尖点频响函数。

步骤4：计算R_7b7b，R_7b7b是子结构VII末端的频响函数矩阵，可表示为

$R_{7b7b}=\left[\begin{array}{cc}{H}_{7b7b}& L_{7b7b}\\ N_{7b7b}& P_{7b7b}\end{array}\right]$

式中，H_7b7b可通过对子结构VII进行锤击试验得到。并采用一阶有限差分法计算R_7b7b中其余参数L_7b7b，N_7b7b，P_7b7b。

$R_{7b7b}=\left[\begin{array}{cc}{H}_{7b7b}& \frac{H_{77}-H_{78}}{S}\\ \frac{H_{77}-H_{78}}{S}& \frac{{N_{7b7b}}^2}{H_{7b7b}}\end{array}\right]$

式中，S表示点7和点8之间的距离，其中点7是刀杆上靠近刀柄末端的点，点8是刀杆上靠近刀齿的点，点7和点8均为激励点。点7位置的选择需要注意两点：(1)便于力锤进行锤击实验；(2)尽可能地接近刀柄末端，因为如果点7距离刀柄太远会影响S的取值。点8和点7之间的距离在保证实验效果的情况下，尽可能大些，从而使H₇₇、H₇₈之间的差异明显，保证R_7b7b计算的准确性。H₇₇表示在点7的原点平动位移频响函数，H₇₈表示在点8的跨点平动位移频响函数，H₇₇、H₇₈可通过锤击试验得到。如图3所示，在点7和点8处分别进行模态锤击实验，使用力锤分别对点7和点8施加激振力，并用安装在铣刀刀杆处与点7正对的加速度传感器记录振动加速度响应，利用数据采集***同步采集力信号和加速度信号，经过多次测量，使用计算机中的ModelVIEW软件对信号进行处理，得到频响函数H₇₇和H₇₈。

步骤5：采用Timoshenko梁模型计算单个子结构I-VI在自由状态下的频响函数矩阵R_ij，其中

$H_{ij}=\underset{r=0}{\overset{\mathrm{\∞}}{\Σ}}\frac{{\mathrm{\Φ}}_r\left(x_i\right){\mathrm{\Φ}}_r\left(x_j\right)}{(1+i\mathrm{\γ}){\mathrm{\ω}}_r^2-{\mathrm{\ω}}^2}$

$N_{ij}=\underset{r=0}{\overset{\mathrm{\∞}}{\Σ}}\frac{{\mathrm{\Φ}}_r^{\′}\left(x_i\right){\mathrm{\Φ}}_r\left(x_j\right)}{(1+i\mathrm{\γ}){\mathrm{\ω}}_r^2-{\mathrm{\ω}}^2}$

$L_{ij}=\underset{r=0}{\overset{\mathrm{\∞}}{\Σ}}\frac{{\mathrm{\Φ}}_r\left(x_i\right){\mathrm{\Φ}}_r^{\′}\left(x_j\right)}{(1+i\mathrm{\γ}){\mathrm{\ω}}_r^2-{\mathrm{\ω}}^2}$

$P_{ij}=\underset{r=0}{\overset{\mathrm{\∞}}{\Σ}}\frac{{\mathrm{\Φ}}_r^{\′}\left(x_i\right){\mathrm{\Φ}}_r^{\′}\left(x_j\right)}{(1+i\mathrm{\γ}){\mathrm{\ω}}_r^2-{\mathrm{\ω}}^2}$

式中， ${\mathrm{\Φ}}_r\left(x\right)=A_r\left(C_{1}sin\right(\frac{{\mathrm{\α}}_r}{L}x)+C_{2}cos(\frac{{\mathrm{\α}}_r}{L}x)+C_3\sinh (\frac{{\mathrm{\β}}_r}{L}x)+C_4\cosh (\frac{{\mathrm{\β}}_r}{L}x\left)\right)$

$C_1=L,C_2=-\frac{D_{11}}{D_{12}}{C}_1,C_3=\frac{{\mathrm{\α}}_r-{\mathrm{\λ}}_r}{{\mathrm{\δ}}_r-{\mathrm{\β}}_r}{C}_1,C_4=-\frac{{\mathrm{\λ}}_r{\mathrm{\α}}_r}{{\mathrm{\β}}_r{\mathrm{\δ}}_r}\frac{D_{11}}{D_{12}}{C}_1$

D₁₁＝(α-λ)(cosα-coshβ)

$D_{12}=(\mathrm{\λ}-\mathrm{\α})sin\mathrm{\α}+\frac{\mathrm{\λ}\mathrm{\α}}{\mathrm{\β}\mathrm{\δ}}(\mathrm{\β}-\mathrm{\δ})\sinh \mathrm{\β}$

$\mathrm{\α}=\sqrt{\mathrm{\Ω}+\mathrm{\ϵ}},\mathrm{\β}=\sqrt{-\mathrm{\Ω}+\mathrm{\ϵ}},\mathrm{\Ω}=\frac{b^2(s^2+R^2)}{2}$

$\mathrm{\ϵ}=b\sqrt{\frac{1}{4}{b}^2{(s^2+R^2)}^2-(b^2{s}^2{R}^2-1)},b^2=\frac{{\mathrm{\ρA\ω}}^2{L}^4}{EI}$

$s^2=\frac{EI}{k^{\′}{\mathrm{AGL}}^2},R^2=\frac{1}{{\mathrm{AL}}^2},\mathrm{\λ}=\mathrm{\α}-\frac{b^2{s}^2}{\mathrm{\α}},k^{\′}=\frac{6(1+\mathrm{\υ})}{7+6\mathrm{\υ}}$

根据振型函数的正交性可得

${\∫}_{x=0}^L{\left\{{\mathrm{\Φ}}_r\left(x\right)\right\}}^T\mathrm{\ρ}A\left\{{\mathrm{\Φ}}_r\left(x\right)\right\}dx=\left\{\begin{array}{c}1,r=s\\ 0,r\≠s\end{array}\right.$

令r，s＝1,2,3…n可得A_r(r＝1,2,3…n)

式中，Φ,Φ'分别表示两端自由的Timoshenko梁的平动振型函数和平动振型函数对x的导数；α_r,β_r均是是无量纲频率参数，A_r是振型函数的模态系数，γ表示阻尼因子，ω表示角频率，r表示模态阶数；ρ、E和G依次表示梁材料的密度、杨氏模量和剪切模量，A、I、L、k’和υ依次表示梁截面的面积、截面惯性矩、长度、剪切系数和梁材料的泊松比。

步骤6：考虑铣刀的精确建模，计算子结构I-VI的截面惯性矩

如图4、图5所示，子结构I-VI的截面形状不同。根据步骤5中H_ij、L_ij、N_ij和P_ij的计算公式，对于同一把铣刀，各子结构I-VI的频响函数计算中除了截面惯性矩不同外，其他参数均为常数，因此，子结构I-VI频响函数计算的差异主要在于子结构I-VI截面惯性矩的不同。考虑到铣刀的精确化建模，子结构I-VI的截面形状不尽相同，因此其截面惯性矩也不同。根据铣刀截面惯性矩的计算可以分为对称型铣刀和非对称型铣刀，即在铣刀截面建立x-y笛卡尔二维坐标系，如果铣刀刀齿截面对x轴的惯性矩和对y轴的惯性矩不等，则该铣刀属于非对称型铣刀，否则属于对称型铣刀。以2刃铣刀和4刃铣刀为例说明两种铣刀各子结构截面惯性矩的计算方法，其中2刃铣刀属于非对称型铣刀，4刃铣刀属于对称型铣刀。

①铣刀刀齿子结构I-IV截面惯性矩的计算

2刃和4刃铣刀截面具有中心对称的特点，先求解1/2的2刃铣刀和1/4的4刃铣刀端面截面惯性矩的计算公式，在铣刀端面处建立2刃铣刀和4刃铣刀刀齿截面的x-y笛卡尔二维坐标系并计算铣刀刀齿的截面惯性矩：

$\left\{\begin{array}{c}{I}_{x,2}={\∫}_0^{\mathrm{\π}}{\∫}_0^{R_{eq2}}{\mathrm{\ρ}}^3{\sin }^2\mathrm{\θ}d\mathrm{\ρ}d\mathrm{\θ}-\frac{{{\mathrm{\πf}}_d}^4}{128}\\ I_{y,2}={\∫}_0^{\mathrm{\π}}{\∫}_0^{R_{eq2}}{\mathrm{\ρ}}^3{\cos }^2\mathrm{\θ}d\mathrm{\ρ}d\mathrm{\θ}-\[\frac{{{\mathrm{\πf}}_d}^4}{128}+\frac{{{\mathrm{\πf}}_d}^2}{8}{\left(\frac{D-f_d}{2}\right)}^2\]\end{array}\right.$

$\left\{\begin{array}{c}{I}_{x,4}={\∫}_0^{\mathrm{\π}/2}{\∫}_0^{R_{eq4}}{\mathrm{\ρ}}^3{\sin }^2\mathrm{\θ}d\mathrm{\ρ}d\mathrm{\θ}-\[\frac{{{\mathrm{\πf}}_d}^4}{128}+\frac{{{\mathrm{\πf}}_d}^2}{8}{\left(\frac{D-f_d}{2}\right)}^2\]\\ I_{y,4}={\∫}_0^{\mathrm{\π}/2}{\∫}_0^{R_{eq4}}{\mathrm{\ρ}}^3{\cos }^2\mathrm{\θ}d\mathrm{\ρ}d\mathrm{\θ}-\frac{{{\mathrm{\πf}}_d}^4}{128}\end{array}\right.$

$\begin{array}{cc}{R}_{eq2}=-acos\mathrm{\θ}+\sqrt{(r^2-a^2)+a^2{\cos }^2\mathrm{\θ}}& 0\≤\mathrm{\θ}\≤\mathrm{\π}\end{array}$

$\begin{array}{cc}{R}_{eq4}=asin\mathrm{\θ}+\sqrt{r^2-a^2{\cos }^2\mathrm{\θ}}& 0\≤\mathrm{\θ}\≤\mathrm{\π}/2\end{array}$

式中，r表示铣刀切削刃圆弧半径，a表示铣刀切削刃圆弧中心到铣刀中心的距离，f_d表示铣刀螺旋槽直径，R_eqi表示i刃铣刀刀齿截面的等效半径，i＝2,4，I_x,j表示j刃铣刀刀齿截面对x轴的惯性矩，I_y,j表示j刃铣刀刀齿截面对y轴的惯性矩，j＝2,4；

根据上式将I_x,2、I_y,2展开，由于I_x,2、I_y,2的展开式较为复杂，因此通过中间量I₂₁、I₂₂、I₂₃和矩阵A、B、C表示：

$\left\{\begin{array}{c}{I}_{x,2}=\frac{AB}{4}-\frac{{\mathrm{\πf}}_d^4}{128}\\ I_{y,2}=\frac{AC}{4}-\[\frac{{\mathrm{\πf}}_d^4}{128}+\frac{{\mathrm{\πf}}_d^2}{8}{\left(\frac{D-f_d}{2}\right)}^2\]\end{array}\right.$

$\left\{\begin{array}{c}{I}_{21}=\frac{2r^4}{a^3}(\frac{t}{8}-\frac{\cos t\sin t}{8}+\frac{\cos {\mathrm{tsin}}^{3}t}{4}){|}_0^{\arcsin \frac{a}{r}}\\ I_{22}=\frac{4r^6}{a^5}(\frac{t}{16}-\frac{\sin 2t}{24}+\frac{\sin 4t}{192}+\frac{\cos t{\sin }^{5}t}{6}){|}_0^{\arcsin \frac{a}{r}}\\ I_{23}=\frac{2r^2}{a}(\frac{t}{2}+\frac{\sin }{4}){|}_0^{\arcsin \frac{a}{r}}\end{array}\right.$

$\left\{\begin{array}{c}A=\[r^4,a^4,2r^2{a}^2,4a^2{r}^2,-4a^4,-4{\mathrm{ar}}^2,-4a^3\]\\ B={\[\frac{\mathrm{\π}}{2},\frac{\mathrm{\π}}{4},\frac{-\mathrm{\π}}{4},\frac{\mathrm{\π}}{8},\frac{\mathrm{\π}}{16},I_{21},I_{21}-I_{22}\]}^T\\ C={\[\frac{\mathrm{\π}}{2},\frac{\mathrm{\π}}{4},\frac{\mathrm{\π}}{4},\frac{3\mathrm{\π}}{8}\frac{\mathrm{\π}}{16},I_{23}-I_{21},I_{23}-3I_{21}+I_{22}\]}^T\end{array}\right.$

同理将I_x,4、I_y,4展开，由于I_x,4、I_y,4的展开式较为复杂，因此通过中间量I₄₁、I₄₂、I₄₃和矩阵D、E、F表示：

$\left\{\begin{array}{c}{I}_{x,4}=\frac{DE}{4}-\[\frac{{\mathrm{\πf}}_d^4}{128}+\frac{{\mathrm{\πf}}_d^2}{8}{\left(\frac{D-f_d}{2}\right)}^2\]\\ I_{y,4}=\frac{DF}{4}-\frac{{\mathrm{\πf}}_d^4}{128}\end{array}\right.$

$\left\{\begin{array}{c}{I}_{41}=\frac{r^2}{a}(\frac{t}{2}+\frac{\sin 2t}{4}){|}_0^{\arcsin \frac{a}{r}}\\ I_{42}=\frac{r^4}{a^3}(\frac{t}{8}-\frac{\cos t\sin t}{8}+\frac{\cos {\mathrm{tsin}}^{3}t}{4}){|}_0^{\arcsin \frac{a}{r}}\\ I_{43}=\frac{2r^6}{a^5}(\frac{t}{16}-\frac{\sin 2t}{24}+\frac{\sin 4t}{192}+\frac{\cos t{\sin }^{5}t}{6}){|}_0^{\arcsin \frac{a}{r}}\end{array}\right.$

$\left\{\begin{array}{c}D=\[r^4,a^4,-2r^2{a}^2,4a^2{r}^2,-4a^4,4{\mathrm{ar}}^2,-4a^3\]\\ E={\[\frac{\mathrm{\π}}{4},\frac{\mathrm{\π}}{8},\frac{-\mathrm{\π}}{8},\frac{3\mathrm{\π}}{16},\frac{\mathrm{\π}}{32},I_{41}-I_{42},-2I_{41}+4I_{42}-I_{43}\]}^T\\ F={\[\frac{\mathrm{\π}}{4},\frac{\mathrm{\π}}{8},\frac{\mathrm{\π}}{8},\frac{\mathrm{\π}}{8},\frac{\mathrm{\π}}{32},I_{42},I_{43}-I_{42}\]}^T\end{array}\right.$

式中通过铣刀截面的几何边界条件确定r，a，对于2刃铣刀通过求解得到：

$\left\{\begin{array}{c}r=(D-f_d)/2\\ a=f_d/2\end{array}\right.$

同理可得4刃铣刀刀齿的参数r，a：

$\left\{\begin{array}{c}r=\frac{D^2-2{\mathrm{Df}}_d+2{f_d}^2}{2D}\\ a=\frac{{\mathrm{Df}}_d-{f_d}^2}{D}\end{array}\right.$

然后根据中心对称可得到2刃铣刀和4刃铣刀刀齿端面剩余部分的截面惯性矩，经过叠加最终得到2刃和4刃铣刀刀齿端面总的截面惯性矩：

$\begin{array}{cc}\left\{\begin{array}{c}{I}_{x,2tot}=2I_{x,2}\\ I_{y,2tot}=2I_{y,2}\end{array}\right.& \left\{\begin{array}{c}{I}_{y,4tot}=2(I_{y,4}+I_{x,4})\\ I_{x,4tot}=I_{y,4tot}\end{array}\right.\end{array}$

②计算距离铣刀刀齿端面不同距离处的截面惯性矩

距离铣刀刀齿端面不同距离处的截面和刀齿端面的形状相同，不同的是其在位置上相对刀齿端面转动了一个滞后角ψ，为便于方便，在距离铣刀刀齿端面不同距离处的截面建立局部坐标系u-v，在刀齿端面建立全局坐标系x-y，并将局部坐标系u-v上的截面惯性矩转换到全局坐标系x-y上，转换公式如下所示：

$\left\{\begin{array}{c}{I}_{xtot}={\cos }^2{\mathrm{\ψI}}_u+{\sin }^2{\mathrm{\ψI}}_v-I_{uv}\sin 2\mathrm{\ψ}\\ I_{ytot}={\cos }^2{\mathrm{\ψI}}_u+{\sin }^2{\mathrm{\ψI}}_v-I_{uv}\sin 2\mathrm{\ψ}\end{array}\right.$

$\mathrm{\ψ}=\frac{2ztan\mathrm{\β}}{D}$

式中，z表示刀齿任意截面距离铣刀刀齿端面的距离；β表示铣刀螺旋角；D表示铣刀刀齿圆柱基体直径。对于2刃和4刃铣刀，I_uv＝0，将①中的I_x,2tot和I_y,2tot、I_x,4tot和I_y,4tot替换I_u和I_v，可分别得到2刃铣刀和4刃铣刀刀齿部分距离铣刀刀齿端面不同距离处的截面惯性矩；

③铣刀过渡段子结构截面惯性矩的计算

铣刀过渡段截面惯性矩的计算采用分块计算，然后叠加的方法。根据铣刀截面的特点，提出了一种通用的划分方法，即在铣刀过渡段截面建立以其形心为坐标原点的笛卡尔二维坐标系，然后以其原点为中心，划分出若干块三角形区域和扇形区域，通过分别计算每块区域的截面惯性矩，利用坐标转换公式，经过叠加，最后得到整体的截面惯性矩。

I_tx表示三角形截面对x轴的惯性矩；I_ty表示三角形截面对y轴的惯性矩；I_sx表示扇形截面对x轴的惯性矩；I_sy表示扇形截面对y轴的惯性矩；

根据坐标转化，经过叠加，得到2刃铣刀过渡段总的截面惯性矩的计算公式：

$\left\{\begin{array}{c}{I}_{x-2tot}=2(I_{tx}+I_{sx})\\ I_{y-2tot}=2(I_{ty}+I_{sy})\end{array}\right.$

同理可得4刃铣刀过渡段总的截面惯性矩的计算公式：

$\left\{\begin{array}{c}{I}_{x-4tot}=2(I_{tx}+I_{sx}+I_{ty}+I_{sy})\\ I_{y-4tot}=I_{x-4tot}\end{array}\right.$

④铣刀刀杆子结构截面惯性矩的计算

铣刀刀杆部分是一段圆柱段，其截面为圆形截面，故其截面惯性矩的计算公式为

$I_x=I_y=\frac{{\mathrm{\πD}}^2}{64}$

将①②③④带入步骤5中的频响函数计算公式H_ij、L_ij、N_ij和P_ij，可得到子结构I-VI自由状态下的频响函数矩阵。

步骤7：将步骤6中计算得到的子结构I-VI的频响函数矩阵代入步骤3中的RCSA耦合公式，得到RC_7a7a、RC_7a1、RC_17a、RC₁₁。然后将RC_7a7a、RC_7a1、RC_17a、RC₁₁和R_7b7b代入步骤3中的刀尖点频响函数矩阵计算公式，最终得到刀尖点频响函数矩阵G₁₁，即可获得刀尖点频响函数HC₁₁。

步骤8：按照图11，对铣刀进行模态锤击试验获得铣刀在x、y向上的实测刀尖点频响函数，测试结果如图13所示，并对比步骤7中预测得到的铣刀刀尖点频响函数，测试结果如图14所示。

由前述理论可知，2刃铣刀属于非对称型铣刀，4刃铣刀属于对称型铣刀，因此2刃铣刀的理论刀尖点频响函数在x、y向上不同，4刃铣刀的理论刀尖点频响函数在x、y向上相同。为了验证理论的正确性，对2刃铣刀和4刃铣刀进行试验，如图13所示。通过实验可知2刃铣刀的实测刀尖点频响函数在x、y向上不同，4刃铣刀的实测刀尖点频响函数在x、y向上基本相同。说明将铣刀分为对称型铣刀和非对称型铣刀分别进行精确建模是有必要的。为验证本发明所提方法的正确性和适用性，结合实施例，采用5把铣刀对比本发明预测的铣刀刀尖点频响函数和实测刀尖点频响函数的差异，以验证本发明对不同材料、直径、长度以及刃数的铣刀的适用性，如图14所示。

采用本发明专利的有益效果在于，本发明提出了一种基于铣刀精确建模的刀尖点频响函数的预测方法。该方法将机床-主轴1-刀柄2-铣刀划分为机床-主轴1-刀柄2-部分刀杆、剩余刀杆、铣刀过渡段以及多段刀齿；考虑铣刀的实际建模，根据实际情况将铣刀分为对称型铣刀和非对称型铣刀，以2刃铣刀和4刃铣刀为例分别详细计算了铣刀各部分的截面惯性矩，并提出一种铣刀过渡段的通用建模方法，提高铣刀刀尖点频响函数的预测精度。

以下结合实施例介绍本发明专利的实用性。

实施例一：铣刀T1，直径为10mm，总长为95mm，材料为高速钢。7、8点之间的距离为20mm。

按照前述步骤得到铣刀T1刀尖点频响函数的预测值和实测值。如图1所示，并得到T1铣刀固有频率的预测值和实测值，如表1所示。

表1T1铣刀固有频率的预测值和实测值对比

实施例二：铣刀T2，2刃铣刀，直径为10mm，总长为100mm，材料为硬质合金。7、8点之间的距离为25mm。

按照前述步骤得到铣刀T2刀尖点频响函数的预测值和实测值。如图1所示，并得到T2铣刀固有频率的预测值和实测值，如表2所示。

表2T2铣刀固有频率的预测值和实测值对比

实施例三：铣刀T3，2刃铣刀，直径为12mm，总长为100mm，材料为硬质合金。7、8点之间的距离为25mm。

按照前述步骤得到铣刀T3刀尖点频响函数的预测值和实测值。如图1所示，并得到T3铣刀固有频率的预测值和实测值，如表3所示。

表3T3铣刀固有频率的预测值和实测值对比

实施例四：铣刀T4，4刃铣刀，直径为10mm，总长为100mm，材料为硬质合金。7、8点之间的距离为25mm。

按照前述步骤得到铣刀T4刀尖点频响函数的预测值和实测值。如图1所示，并得到T4铣刀固有频率的预测值和实测值，如表4所示。

表4T4铣刀固有频率的预测值和实测值对比

实施例五：铣刀T5，4刃铣刀，直径为10mm，总长为75mm，材料为硬质合金。7、8点之间的距离为20mm。

表5T5铣刀固有频率的预测值和实测值对比

由以上五个实施例可以看出，本发明预测和实测固有频率前几阶的相对误差在5％以内，说明本发明预测的结果具有更高的精度，且本发明对不同铣刀参数(材料、直径、长度以及刃数)铣刀都具有适用性。本发明对于不同材料、直径、长度以及刃数的铣刀均试用。从而证明了本发明的实用性。

以上所述仅为本发明专利的实施例，并非因此限制本发明专利的专利范围，凡是利用本发明专利说明书及附图内容所作的等效结构或等效流程变换，或直接或间接运用在其他相关的技术领域，均同理包括在本发明专利的专利保护范围内。

Claims (5)

1.基于铣刀精确建模的刀尖点频响函数预测方法，其特征在于包括以下步骤：

步骤1：将铣刀一端安装到刀柄中，将刀柄带有锥度的一端安装在与刀柄尾部锥度相配合的机床-主轴中，从而组成机床-主轴-刀柄-铣刀***；

步骤2：将机床-主轴-刀柄-铣刀***划分为机床-主轴-刀柄-部分刀杆、剩余刀杆、铣刀过渡段、刀齿四个部分，并将刀齿划分为多段，其中子结构I-IV为刀齿部分，子结构V为铣刀过渡段部分，子结构VI为剩余刀杆部分，子结构VII为机床-主轴-刀柄-部分刀杆部分；

步骤3：利用RCSA方法耦合子结构I-VII获得刀尖点频响函数矩阵；

步骤4：计算R_7b7b，R_7b7b是子结构VII末端的频响函数矩阵；

步骤5：采用Timoshenko梁模型得到单个子结构I-VI在自由状态下的频响函数矩阵R_ij的表达式；

步骤6：考虑铣刀的精确建模，计算子结构I-VI的截面惯性矩，结合步骤5，得到子结构I-VI的频响函数矩阵；

步骤7：将步骤6中计算得到的子结构I-VI的频响函数矩阵代入步骤3中的RCSA耦合公式，得到RC_7a7a、RC_7a1、RC_17a、RC₁₁，然后将RC_7a7a、RC_7a1、RC_17a、RC₁₁和R_7b7b代入步骤3中的刀尖点频响函数矩阵计算公式，最终得到刀尖点频响函数矩阵G₁₁，即可获得刀尖点频响函数HC₁₁；

步骤8：对铣刀进行模态锤击试验获得铣刀在x、y向上的实测刀尖点频响函数，并对比步骤7中预测得到的铣刀刀尖点频响函数。

2.根据权利要求1所述的预测方法，其特征在于：

步骤3具体为：利用RCSA方法耦合子结构I-VI获得刀尖点频响函数矩阵，两端自由的子结构的频响函数矩阵表示为如下通式：

$R_{ij}=\left[\begin{array}{cc}{H}_{ij}& L_{ij}\\ N_{ij}& P_{ij}\end{array}\right]=\left[\begin{array}{cc}{X}_i/F_j& X_i/M_j\\ {\mathrm{\θ}}_i/F_j& {\mathrm{\θ}}_i/M_j\end{array}\right]$

式中，i，j＝1，2a或i，j＝2b，3a或i，j＝3b，4a或i，j＝4b，5a或i，j＝5b，6a或i，j＝6b，7a；

H_ij、L_ij、N_ij和P_ij依次表示为单个两端自由的子结构在j点激励，在i点获得响应的位移/力、位移/力矩、转角/力以及转角/力矩的频响函数，F_j、M_j依次表示为j点处所受到的外力和外力矩；X_i、Θ_i依次表示为在外力和外力矩的作用下，i点处的平动位移和转角；1，2a，2b，3a，3b，4a，4b，5a，5b，6a，6b，7a分别表示子结构I-VII从铣刀刀尖到机床的两端端点；

其中，自由状态下相邻两个子结构耦合后的频响函数矩阵，表示为以下通式：

${\mathrm{RC}}_{ij}=\left[\begin{array}{cc}{\mathrm{HC}}_{ij}& {\mathrm{LC}}_{ij}\\ {\mathrm{NC}}_{ij}& {\mathrm{PC}}_{ij}\end{array}\right]$

式中，RC_ij表示自由状态下相邻两个子结构耦合后的频响函数矩阵，i，j＝1，3a，或i，j＝1，4a，或i，j＝1，5a，或i，j＝1，6a，或i，j＝1，7a，HC_ij、LC_ij、NC_ij和PC_ij依次表示子结构I-II在j点激励，在i点获得响应的位移/力、位移/力矩、转角/力以及转角/力矩的频响函数；

根据自由状态下相邻两个子结构在结合处的兼容条件和平衡方程，最终计算得到自由状态下相邻两个子结构的频响函数矩阵与自由状态下相邻两个子结构频响函数矩阵之间的关系，以子结构I和子结构II进行说明：

$\left\{\begin{array}{c}{\mathrm{RC}}_{11}=R_{11}-R_{12a}{(R_{2a2a}+R_{2b2b})}^{-1}{R}_{2a1}\\ {\mathrm{RC}}_{13a}=R_{12a}{(R_{2a2a}+R_{2b2b})}^{-1}{R}_{2b3a}\\ {\mathrm{RC}}_{3a1}=R_{3a2b}{(R_{2a2a}+R_{2b2b})}^{-1}{R}_{2a1}\\ {\mathrm{RC}}_{3a3a}=R_{3a3a}-R_{3a2b}{(R_{2a2a}+R_{2b2b})}^{-1}{R}_{2b3a}\end{array}\right.$

同理通过上式顺序耦合子结构I-VII，最终得到刀尖点频响函数矩阵G₁₁：

G₁₁＝RC₁₁-RC_17a(RC_7a7a+R_7b7b)^-1RC_7a1

式中，G₁₁还表示为：

$G_{11}=\left[\begin{array}{cc}{\mathrm{HC}}_{11}& {\mathrm{LC}}_{11}\\ {\mathrm{NC}}_{11}& {\mathrm{PC}}_{11}\end{array}\right]$

式中，HC₁₁、LC₁₁、NC₁₁和PC₁₁依次为刀尖点的位移/力、位移/力矩、转角/力以及转角/力矩频响函数，其中HC₁₁简称为刀尖点频响函数。

3.根据权利要求1所述的预测方法，其特征在于：

步骤4具体为：计算R_7b7b，R_7b7b是子结构VII末端的频响函数矩阵，表示为：

$R_{7b7b}=\left[\begin{array}{cc}{H}_{7b7b}& L_{7b7b}\\ N_{7b7b}& P_{7b7b}\end{array}\right]$

式中，H_7b7b通过对子结构VII进行锤击试验得到，并采用一阶有限差分法计算R_7b7b中其余参数L_7b7b，N_7b7b，P_7b7b：

$R_{7b7b}=\left[\begin{array}{cc}{H}_{7b7b}& \frac{H_{77}-H_{78}}{S}\\ \frac{H_{77}-H_{78}}{S}& \frac{{N_{7b7b}}^2}{H_{7b7b}}\end{array}\right]$

式中，S表示刀杆上接近刀柄末端的点7和刀杆上接近刀刃的点8之间的距离，其中点7是刀杆上靠近刀柄末端的点，点8是刀杆上靠近刀齿部分的点，H₇₇表示在点7的原点平动位移频响函数，H₇₈表示在点8的跨点平动位移频响函数，H₇₇、H₇₈通过锤击试验得到，即在点7和点8处分别进行模态锤击实验，使用力锤分别对点7和点8施加激振力，并用安装在铣刀刀杆处与点7正对的加速度传感器记录振动加速度响应，利用数据采集***同步采集力信号和加速度信号，经过多次测量，使用计算机中的ModelVIEW软件对信号进行处理，得到频响函数H₇₇和H₇₈。

4.根据权利要求1所述的预测方法，其特征在于：

步骤5具体为：采用Timoshenko梁模型计算单个子结构I-VI在自由状态下的频响函数矩阵R_ij，其中

$H_{ij}=\underset{r=0}{\overset{\mathrm{\∞}}{\Σ}}\frac{{\mathrm{\Φ}}_r\left(x_i\right){\mathrm{\Φ}}_r\left(x_j\right)}{(1+i\mathrm{\γ}){\mathrm{\ω}}_r^2-{\mathrm{\ω}}^2}$

$N_{ij}=\underset{r=0}{\overset{\mathrm{\∞}}{\Σ}}\frac{{\mathrm{\Φ}}_r^{\′}\left(x_i\right){\mathrm{\Φ}}_r\left(x_j\right)}{(1+i\mathrm{\γ}){\mathrm{\ω}}_r^2-{\mathrm{\ω}}^2}$

$L_{ij}=\underset{r=0}{\overset{\mathrm{\∞}}{\Σ}}\frac{{\mathrm{\Φ}}_r\left(x_i\right){\mathrm{\Φ}}_r^{\′}\left(x_j\right)}{(1+i\mathrm{\γ}){\mathrm{\ω}}_r^2-{\mathrm{\ω}}^2}$

$P_{ij}=\underset{r=0}{\overset{\mathrm{\∞}}{\Σ}}\frac{{\mathrm{\Φ}}_r^{\′}\left(x_i\right){\mathrm{\Φ}}_r^{\′}\left(x_j\right)}{(1+i\mathrm{\γ}){\mathrm{\ω}}_r^2-{\mathrm{\ω}}^2}$

式中， ${\mathrm{\Φ}}_r\left(x\right)=A_r\left(C_{1}sin\right(\frac{{\mathrm{\α}}_r}{L}x)+C_{2}cos(\frac{{\mathrm{\α}}_r}{L}x)+C_3\sinh (\frac{{\mathrm{\β}}_r}{L}x)+C_4\cosh (\frac{{\mathrm{\β}}_r}{L}x\left)\right)$

C₁＝L, $C_2=-\frac{(\mathrm{\α}-\mathrm{\λ})(cos\mathrm{\α}-\cosh \mathrm{\β})}{(\mathrm{\λ}-\mathrm{\α})sin\mathrm{\α}+\frac{\mathrm{\λ}\mathrm{\α}}{\mathrm{\β}\mathrm{\δ}}(\mathrm{\β}-\mathrm{\δ})sinh\mathrm{\β}}L,$ $C_3=\frac{{\mathrm{\α}}_r-{\mathrm{\λ}}_r}{{\mathrm{\δ}}_r-{\mathrm{\β}}_r}L,$

$C_4=-\frac{{\mathrm{\λ}}_r{\mathrm{\α}}_r}{{\mathrm{\β}}_r{\mathrm{\δ}}_r}\frac{(\mathrm{\α}-\mathrm{\λ})(cos\mathrm{\α}-\cosh \mathrm{\β})}{(\mathrm{\λ}-\mathrm{\α})sin\mathrm{\α}+\frac{\mathrm{\λ}\mathrm{\α}}{\mathrm{\β}\mathrm{\δ}}(\mathrm{\β}-\mathrm{\δ})sinh\mathrm{\β}}L$

$\mathrm{\α}=\sqrt{\mathrm{\Ω}+\mathrm{\ϵ}},$ $\mathrm{\β}=\sqrt{-\mathrm{\Ω}+\mathrm{\ϵ}},$ $\mathrm{\Ω}=\frac{b^2(s^2+R^2)}{2}$

$\mathrm{\ϵ}=b\sqrt{\frac{1}{4}{b}^2{(s^2+R^2)}^2-(b^2{s}^2{R}^2-1)},$ $b^2=\frac{{\mathrm{\ρA\ω}}^2{L}^4}{EI}$

$s^2=\frac{EI}{k^{\′}{\mathrm{AGL}}^2},$ $R^2=\frac{1}{{\mathrm{AL}}^2},$ $\mathrm{\λ}=\mathrm{\α}-\frac{b^2{s}^2}{\mathrm{\α}},$ $k^{\′}=\frac{6(1+\mathrm{\υ})}{7+6\mathrm{\υ}}$

根据振型函数的正交性可得

${\∫}_{x=0}^L{\left\{{\mathrm{\Φ}}_r\left(x\right)\right\}}^T\mathrm{\ρ}A\left\{{\mathrm{\Φ}}_r\left(x\right)\right\}dx=\left\{\begin{array}{c}1,r=s\\ 0,r\≠s\end{array}\right.$

令r，s＝1,2,3…n可得A_r(r＝1,2,3…n)

式中，Φ,Φ'分别表示两端自由的Timoshenko梁的平动振型函数和平动振型函数对x的导数；α_r,β_r均是是无量纲频率参数，A_r是振型函数的模态系数，γ表示阻尼因子，ω表示角频率，r表示模态阶数；ρ、E和G依次表示梁材料的密度、杨氏模量和剪切模量，A、I、L、k’和υ依次表示梁截面的面积、截面惯性矩、长度、剪切系数和梁材料的泊松比。

5.根据权利要求1所述的预测方法，其特征在于：

步骤6具体为：考虑铣刀的精确建模，计算子结构I-VI的截面惯性矩，结合步骤5，得到子结构I-VI的频响函数矩阵：

根据铣刀截面惯性矩的计算将铣刀分为对称型铣刀和非对称型铣刀，即在铣刀截面建立x-y笛卡尔二维坐标系，如果铣刀刀齿截面对x轴的惯性矩和对y轴的惯性矩不等，则该铣刀属于非对称型铣刀，否则属于对称型铣刀；以下计算以2刃铣刀和4刃铣刀各子结构截面惯性矩，其中2刃铣刀属于非对称型铣刀，4刃铣刀属于对称型铣刀；

①铣刀刀齿子结构I-IV截面惯性矩的计算

2刃和4刃铣刀截面具有中心对称的特点，先求解1/2的2刃铣刀和1/4的4刃铣刀端面截面惯性矩的计算公式，在铣刀端面处建立2刃铣刀和4刃铣刀刀齿截面的x-y笛卡尔二维坐标系并计算铣刀刀齿的截面惯性矩：

$\left\{\begin{array}{c}{I}_{x,2}={\∫}_0^{\mathrm{\π}}{\∫}_0^{R_{eq2}}{\mathrm{\ρ}}^3{\sin }^2\mathrm{\θ}d\mathrm{\ρ}d\mathrm{\θ}-\frac{{{\mathrm{\πf}}_d}^4}{128}\\ I_{y,2}={\∫}_0^{\mathrm{\π}}{\∫}_0^{R_{eq2}}{\mathrm{\ρ}}^3{\cos }^2\mathrm{\θ}d\mathrm{\ρ}d\mathrm{\θ}-\[\frac{{{\mathrm{\πf}}_d}^4}{128}+\frac{{{\mathrm{\πf}}_d}^2}{8}{\left(\frac{D-f_d}{2}\right)}^2\]\end{array}\right.$

$\left\{\begin{array}{c}{I}_{x,4}={\∫}_0^{\mathrm{\π}/2}{\∫}_0^{R_{eq4}}{\mathrm{\ρ}}^3{\sin }^2\mathrm{\θ}d\mathrm{\ρ}d\mathrm{\θ}-\[\frac{{{\mathrm{\πf}}_d}^4}{128}+\frac{{{\mathrm{\πf}}_d}^2}{8}{\left(\frac{D-f_d}{2}\right)}^2\]\\ I_{y,4}={\∫}_0^{\mathrm{\π}/2}{\∫}_0^{R_{eq4}}{\mathrm{\ρ}}^3{\cos }^2\mathrm{\θ}d\mathrm{\ρ}d\mathrm{\θ}-\frac{{{\mathrm{\πf}}_d}^4}{128}\end{array}\right.$

$\begin{array}{cc}{R}_{eq2}=-acos\mathrm{\θ}+\sqrt{(r^2-a^2)+a^2{\cos }^2\mathrm{\θ}}& 0\≤\mathrm{\θ}\≤\mathrm{\π}\end{array}$

$\begin{array}{cc}{R}_{eq4}=asin\mathrm{\θ}+\sqrt{r^2-a^2{\cos }^2\mathrm{\θ}}& 0\≤\mathrm{\θ}\≤\mathrm{\π}/2\end{array}$

式中，r表示铣刀切削刃圆弧半径，a表示铣刀切削刃圆弧中心到铣刀中心的距离，f_d表示铣刀螺旋槽直径，R_eqi表示i刃铣刀刀齿截面的等效半径，i＝2,4，I_x,j表示j刃铣刀刀齿截面对x轴的惯性矩，I_y,j表示j刃铣刀刀齿截面对y轴的惯性矩，j＝2,4；

根据上式将I_x,2、I_y,2展开，由于I_x,2、I_y,2的展开式较为复杂，因此通过中间量I₂₁、I₂₂、I₂₃和矩阵A、B、C表示：

$\left\{\begin{array}{c}{I}_{x,2}=\frac{AB}{4}-\frac{{\mathrm{\πf}}_d^4}{128}\\ I_{y,2}=\frac{AC}{4}-\[\frac{{\mathrm{\πf}}_d^4}{128}+\frac{{\mathrm{\πf}}_d^2}{8}{\left(\frac{D-f_d}{2}\right)}^2\]\end{array}\right.$

$\left\{\begin{array}{c}{I}_{21}=\frac{2r^4}{a^3}(\frac{t}{8}-\frac{\cos t\sin t}{8}+\frac{\cos {\mathrm{tsin}}^{3}t}{4}){|}_0^{\arcsin \frac{a}{r}}\\ I_{22}=\frac{4r^6}{a^5}(\frac{t}{16}-\frac{\sin 2t}{24}+\frac{\sin 4t}{192}+\frac{\cos t{\sin }^{5}t}{6}){|}_0^{\arcsin \frac{a}{r}}\\ I_{23}=\frac{2r^2}{a}(\frac{t}{2}+\frac{\sin 2t}{4}){|}_0^{\arcsin \frac{a}{r}}\end{array}\right.$

$\left\{\begin{array}{c}A=\[r^4,a^4,2r^2{a}^2,4a^2{r}^2,-4a^4,-4{\mathrm{ar}}^2,-4a^3\]\\ B={\[\frac{\mathrm{\π}}{2},\frac{\mathrm{\π}}{4},\frac{-\mathrm{\π}}{4},\frac{\mathrm{\π}}{8},\frac{\mathrm{\π}}{16},I_{21},I_{21}-I_{22}\]}^T\\ C={\[\frac{\mathrm{\π}}{2},\frac{\mathrm{\π}}{4},\frac{\mathrm{\π}}{4},\frac{3\mathrm{\π}}{8}\frac{\mathrm{\π}}{16},I_{23}-I_{21},I_{23}-3I_{21}+I_{22}\]}^T\end{array}\right.$

同理将I_x,4、I_y,4展开，由于I_x,4、I_y,4的展开式较为复杂，因此通过中间量I₄₁、I₄₂、I₄₃和矩阵D、E、F表示：

$\left\{\begin{array}{c}{I}_{x,4}=\frac{DE}{4}-\[\frac{{\mathrm{\πf}}_d^4}{128}+\frac{{\mathrm{\πf}}_d^2}{8}{\left(\frac{D-f_d}{2}\right)}^2\]\\ I_{y,4}=\frac{DF}{4}-\frac{{\mathrm{\πf}}_d^4}{128}\end{array}\right.$

$\left\{\begin{array}{c}{I}_{41}=\frac{r^2}{a}(\frac{t}{2}+\frac{\sin 2t}{4}){|}_0^{\arcsin \frac{a}{r}}\\ I_{42}=\frac{r^4}{a^3}(\frac{t}{8}-\frac{\cos t\sin t}{8}+\frac{\cos {\mathrm{tsin}}^{3}t}{4}){|}_0^{\arcsin \frac{a}{r}}\\ I_{43}=\frac{2r^6}{a^5}(\frac{t}{16}-\frac{\sin 2t}{24}+\frac{\sin 4t}{192}+\frac{\cos t{\sin }^{5}t}{6}){|}_0^{\arcsin \frac{a}{r}}\end{array}\right.$

$\left\{\begin{array}{c}D=\[r^4,a^4,-2r^2{a}^2,4a^2{r}^2,-4a^4,4{\mathrm{ar}}^2,-4a^3\]\\ E={\[\frac{\mathrm{\π}}{4},\frac{\mathrm{\π}}{8},\frac{-\mathrm{\π}}{8},\frac{3\mathrm{\π}}{16},\frac{\mathrm{\π}}{32},I_{41}-I_{42},-2I_{41}+4I_{42}-I_{43}\]}^T\\ F={\[\frac{\mathrm{\π}}{4},\frac{\mathrm{\π}}{8},\frac{\mathrm{\π}}{8},\frac{\mathrm{\π}}{8},\frac{\mathrm{\π}}{32},I_{42},I_{43}-I_{42}\]}^T\end{array}\right.$

式中通过铣刀截面的几何边界条件确定r，a，对于2刃铣刀通过求解得到：

$\left\{\begin{array}{c}r=(D-f_d)/2\\ a=f_d/2\end{array}\right.$

同理可得4刃铣刀刀齿的参数r，a：

$\left\{\begin{array}{c}r=\frac{D^2-2{\mathrm{Df}}_d+2{f_d}^2}{2D}\\ a=\frac{{\mathrm{Df}}_d-{f_d}^2}{D}\end{array}\right.$

然后根据中心对称可得到2刃铣刀和4刃铣刀刀齿端面剩余部分的截面惯性矩，经过叠加最终得到2刃和4刃铣刀刀齿端面总的截面惯性矩：

$\left\{\begin{array}{c}{I}_{x,2tot}=2I_{x,2}\\ I_{y,2tot}=2I_{y,2}\end{array}\right.$ $\left\{\begin{array}{c}{I}_{y,4tot}=2(I_{y,4}+I_{x,4})\\ I_{x,4tot}=I_{y,4tot}\end{array}\right.$

②计算距离铣刀刀齿端面不同距离处的截面惯性矩

距离铣刀刀齿端面不同距离处的截面和刀齿端面的形状相同，不同的是其在位置上相对刀齿端面转动了一个滞后角ψ，为便于方便，在距离铣刀刀齿端面不同距离处的截面建立局部坐标系u-v，在刀齿端面建立全局坐标系x-y，并将局部坐标系u-v上的截面惯性矩转换到全局坐标系x-y上，转换公式如下所示：

$\left\{\begin{array}{c}{I}_{xtot}={\cos }^2{\mathrm{\ψI}}_u+{\sin }^2{\mathrm{\ψI}}_v-I_{uv}\sin 2\mathrm{\ψ}\\ I_{ytot}={\cos }^2{\mathrm{\ψI}}_u+{\sin }^2{\mathrm{\ψI}}_v-I_{uv}\sin 2\mathrm{\ψ}\end{array}\right.$

$\mathrm{\ψ}=\frac{2ztan\mathrm{\β}}{D}$

式中，z表示刀齿任意截面距离铣刀刀齿端面的距离；β表示铣刀螺旋角；D表示铣刀刀齿圆柱基体直径；对于2刃和4刃铣刀，I_uv＝0，将①中的I_x,2tot和I_y,2tot、I_x,4tot和I_y,4tot替换I_u和I_v，可分别得到2刃铣刀和4刃铣刀刀齿部分距离铣刀刀齿端面不同距离处的截面惯性矩；

③铣刀过渡段子结构截面惯性矩的计算

铣刀过渡段截面惯性矩的计算采用分块计算，然后叠加的方法。根据铣刀截面的特点，提出了一种通用的划分方法，即在铣刀过渡段截面建立以其形心为坐标原点的笛卡尔二维坐标系，然后以其原点为中心，划分出若干块三角形区域和扇形区域，通过分别计算每块区域的截面惯性矩，利用坐标转换公式，经过叠加，最后得到整体的截面惯性矩；

I_tx表示三角形截面对x轴的惯性矩；I_ty表示三角形截面对y轴的惯性矩；I_sx表示扇形截面对x轴的惯性矩；I_sy表示扇形截面对y轴的惯性矩；

根据坐标转化，经过叠加，得到2刃铣刀过渡段总的截面惯性矩的计算公式：

$\left\{\begin{array}{c}{I}_{x-2tot}=2(I_{tx}+I_{sx})\\ I_{y-2tot}=2(I_{ty}+I_{sy})\end{array}\right.$

同理可得4刃铣刀过渡段总的截面惯性矩的计算公式：

$\left\{\begin{array}{c}{I}_{x-4tot}=2(I_{tx}+I_{sx}+I_{ty}+I_{sy})\\ I_{y-4tot}=I_{x-4tot}\end{array}\right.$

④铣刀刀杆子结构截面惯性矩的计算

铣刀刀杆部分是一段圆柱段，其截面为圆形截面，故其截面惯性矩的计算公式为

$I_x=I_y=\frac{{\mathrm{\πD}}^2}{64}$

将①②③④带入步骤5中的频响函数计算公式H_ij、L_ij、N_ij和P_ij，可得到子结构I-VI自由状态下的频响函数矩阵。