一种基于二分法的机床加工稳定性边界快速求解方法与流程

本发明涉及机床加工参数优化和铣削稳定性分析领域，尤其涉及稳定性边界的求解，更具体地，涉及一种基于二分法的机床加工稳定性边界快速求解方法。
背景技术：
：复杂自由曲面类零件(如叶轮、叶片、螺旋桨等)在加工过程中由于薄壁、悬臂等结构特征，自身刚度弱，在加工过程中如果加工参数选择不当，常常会导致加工过程失稳，发生颤振等现象，进而引起加工缺陷，造成设备故障等问题，并且会加速刀具的磨损。其中，加工失稳是指在加工过程中由于刀具与工件之间的相互作用发生强烈振动，其主要由共振和颤振引起，相对于加工共振，加工颤振难以避免。因此，为了减少加工颤振的影响，需要选择合理的加工参数，并进行加工稳定性分析，求解加工稳定性边界(又称叶瓣图)。目前，求解加工稳定性有多种方法，如Altintas、Budak等人提出的零阶系数平均法(ZOA)，Bayly等人提出的时间有限元法(TFEA)，Insperger、Stepan等人提出的半离散法(SDM)以及丁烨等人提出的全离散法(FDM)和数值积分法(NIM)(YeDing，LiMinZhu，XiaoJianZhangandHanDing，NumericalIntegrationMethodforPredictionofMillingStability，JournalofManufacturingScienceandEngineering，133(3)，031005，Jun08，2011)等。同等条件下，数值积分法的计算效率优于其他方法，且精度较高，是目前优秀的求解加工稳定性的方法之一。但是考虑稳定域的求解需要离散参数域平面，数值积分法的计算速度仍然不够快，尤其在参数域离散点较多的时候，计算耗时仍然较长，难以满足工程实际需求。技术实现要素：针对现有技术的以上缺陷或改进需求，本发明提供了一种基于二分法的机床加工稳定性边界快速求解方法，其采用二维二分法优化参数域内参与计算的节点数目，并利用数值积分法求解参数域离散网格节点处动力学方程传递矩阵的特征根，能大幅缩减参数域内节点的数目，减少计算时间，达到快速求解稳定性边界的目的。为实现上述目的，本发明提出了一种基于二分法的机床加工稳定性边界快速求解方法，包括如下步骤：(1)预设需求解稳定性边界的参数平面的边界条件、判稳函数f(x)和迭代次数，根据所述边界条件获得所需的参数平面，对所述参数平面进行初步划分，将其划分为P×Q个网格；其中：所述判稳函数f(x)具体为给定网格节点处依据Floquet定理给出的判稳条件，其函数值采用如下公式计算：f(x)＝ρ(Φ)-1；式中：ρ(Φ)＝max{|λi|}为Φ的谱半径，Φ为给定网格节点处由数值积分法求得的Floquet传递矩阵，λi为传递矩阵Φ的特征值；(2)利用二维二分法对所述每个网格进行再次划分，将其划分为更小的子网格；(3)在每个子网格的顶点处利用数值积分法求解所述判稳函数f(x)的函数值，对于每个子网格，若四个顶点中f(x)的函数值有异号，则该子网格为包含网格，否则，则该子网格为非包含网格；(4)将所述非包含网格执行步骤(2)～(3)一到两次，如果仍然为非包含网格，则结束；否则，获得新的包含网格，并转向步骤(5)；(5)将步骤(3)和(4)中获得的所有包含网格重复执行步骤(2)～(4)，逐步逼近f(x)曲线，直至达到所述预设的迭代次数；(6)将最后获得的所有包含网格进行线性插值，得到近似的判稳函数f(x)的零点，根据所述零点绘制散点图，获得所需的稳定性边界。作为进一步优选的，所述参数平面为二维平面，该二维平面以切削速度为横轴，以切削深度为纵轴。作为进一步优选的，所述P和Q均小于5，所述二维二分法具体指将每个平面有界区域等分为四个与原区域相似的更小的有界子区域。作为进一步优选的，所述利用数值积分法求解判稳函数f(x)的函数值具体指：利用数值积分法求得在加工中给定网格节点处的Floquet传递矩阵，将传递矩阵的谱半径减去1得到函数值。作为进一步优选的，所述四个顶点中f(x)的函数值有异号具体指每个网格四个顶点处的函数值符号不同时为正号或不同时为负号。作为进一步优选的，所述将最后获得的所有包含网格进行线性插值，得到近似的判稳函数f(x)的零点具体包括如下步骤：设包含网格的两异号顶点坐标及其函数值分别为(x1,y1,f1)，(x2,y2,f2)，由于函数值f1与f2异号，则此两顶点之间存在函数值为0的零点，设该零点的坐标及函数值为(x0,y0,0)，然后利用下式计算零点的坐标值，获得判稳函数f(x)的零点：x0＝x1-f1·(x2-x1)/(f2-f1)＝(x1·f2-x2·f1)/(f2-f1)；y0＝y1-f1·(y2-y1)/(f2-f1)＝(y1·f2-y2·f1)/(f2-f1)。作为进一步优选的，所述步骤(3)还包括子步骤(3.1)：将各个子网格顶点的坐标及其对应的函数值进行储存，并建立查找表。总体而言，通过本发明所构思的以上技术方案与现有技术相比，主要具备以下的技术优点：1.本发明采用二维二分法对整个参数平面进行划分，通过迭代获取包含稳定性边界的部分，以优化参数域内参与计算的节点数目，从而减少计算时间，计算过程中能够大幅缩减参数域内节点的数目，从而快速获取铣削中的稳定性边界，为加工中选择合适的加工工艺参数提供依据，实现加工生产中高效、精密加工。2.本发明中采用的二维二分法具有鲁棒性强，速度快，精度适中的特点，可应用于多种加工环境下的稳定性边界的快速求解，通过对加工稳定性进行分析，可以有效保证加工质量，从而达到对零件进行工艺参数优化以及高速高精加工的目标。附图说明图1(a)是二维二分法中使用的所有网格节点示意图；图1(b)是采用二维二分法得到的函数图像示意图；图1(c)是直接绘制得到的函数图像示意图；图2是本发明方法的具体实施流程图；图3是单自由度铣削实例的二维二分法与数值积分法对比图；图4是二自由度铣削实例的二维二分法与数值积分法对比图。具体实施方式为了使本发明的目的、技术方案及优点更加清楚明白，以下结合附图及实施例，对本发明进行进一步详细说明。应当理解，此处所描述的具体实施例仅仅用以解释本发明，并不用于限定本发明。此外，下面所描述的本发明各个实施方式中所涉及到的技术特征只要彼此之间未构成冲突就可以相互组合。本发明的基本原理是采用二维二分法优化参数域内参与计算的节点数目，并利用数值积分法(NIM)求解参数域离散网格节点处动力学方程传递矩阵的特征根，能够大幅缩减参数域内节点的数目，减少计算时间，达到快速求解稳定性边界的目的，为加工中选择合适的加工工艺参数提供理论依据，实现加工生产中高效、精密的加工。本发明的一种基于二分法的机床加工稳定性边界快速求解方法，具体包括如下步骤：(1)预设需求解稳定性边界的参数平面的边界条件和判稳函数f(x)，根据所述边界条件获得所需的参数平面，对所述参数平面进行初步划分，将其划分为P×Q个网格。步骤(1)中参数平面是由横轴为切削速度，纵轴为切削深度的二维平面，函数f(x)是依据Floquet定理给出的判稳条件，f(x)为给定网格节点处Floquet传递矩阵谱半径(特征根绝对值的最大值)减去1，即设矩阵Φ为给定网格节点处由数值积分法求得的Floquet传递矩阵，λi为其特征值，ρ(Φ)＝max{|λi|}为Φ的谱半径，则f(x)＝ρ(Φ)-1，根据Floquet定理，当f(x)等于0时对应着系统的稳定性边界。其中，对参数平面进行初步划分主要是考虑到稳定性边界的复杂性(可能存在突变、尖峰和孤岛等)，若直接使用二维二分法可能会造成误判，造成部分曲线段丢失(例如：很难搜索确定小孤岛的存在，除非经过多次迭代)，所以，本发明预先将参数平面等分为P×Q个网格(P为横向等分个数，Q为纵向等分个数)，其中P和Q均小于5，相当于对整个参数平面进行细分，然后再分别对每个网格进行分析。(2)利用二维二分法对所述每个网格进行再次划分，将其划分为更小的子网格。上述二维二分法是将一维二分法拓展到二维平面的方法。在一维二分法中，每个区间被等分为两个子区间，类似的，在二维二分法中，每个平面有界区域被等分为四个与原区域相似的更小的有界子区域，将此方法应用于稳定域计算时，也是不断对参数平面的有界区域进行划分，随迭代次数增加，划分的网格逐渐变小，其面积反比于迭代次数的平方。因此，通过二维二分法的不断划分，使得包含网格更为精细，在最后进行线性插值时能获得精度较高的稳定性边界。(3)在每个子网格的顶点(节点)处利用数值积分法求解所述判稳函数f(x)，得到其函数值，对于每个子网格，如果四个顶点中f(x)的函数值有异号(即每个网格四个顶点处的函数值符号不同时为正号或不同时为负号，至少存在两个顶点异号)，则认为该子网格为包含网格，即此子网格内至少包含f(x)曲线的一部分；否则，则认为该子网格为非包含网格。在每个子网格的顶点处利用数值积分法求解判稳函数f(x)的函数值具体为：首先利用数值积分法求解给定节点处的传递矩阵，在给定某个子网格的顶点(包含切削转速Ω＝60·T·N和切削深度信息ap)后，可获取该节点处的Floquet传递矩阵，其表达式为：Φ=(I1-C1-apτ2D1)-1(-apτ2D1+E);]]>其中Φ为传递矩阵，矩阵I1为单位矩阵，矩阵C1，D1和E可以通过以下三个方程获得：其中：τ＝(T-tf)/n，表示刀具切入段时间区间[t0+tf,t0+T]被等分成n份。A=-M-1C/2M-1CM-1C/4-K-CM-1/2;]]>Bi=00Kc(ti)0,(i=1,...,n+1);]]>其中M，C，K分别表示系统质量、阻尼和刚度；Kc(t)为切削力系数周期矩阵函数，即Kc(t)＝Kc(t+T)；在给定加工条件和时间区间离散值等条件后，ap和τ可根据网格节点坐标获得：ap等于网格节点纵轴坐标，即切削深度；τ可由网格节点横轴坐标计算得到；计算获得Φ后，根据f(x)＝ρ(Φ)-1即可得出f(x)的函数值。具体的，评估每个子网格顶点(网格节点)处的函数值，通过分析每个网格四个顶点的函数值的正负性可初步判断该网格是不是包含网格。在一维二分法中，对于连续函数，如果区间两端点的函数值异号，则认为此区间有函数零点存在。类似的，在二维二分法中，如果一个网格的四个顶点中存在其f(x)的函数值异号的情况，则认为此网格包含有f(x)函数曲线的一部分，定义为是包含网格；否则，则认为此网格不包含曲线段，定义为非包含网格。另外，由于存在相邻网格均为包含网格的情况，其公共顶点(网格节点)会被重复计算。为解决上述问题，本发明提出了减少要计算的网格节点数目的方法：对获得的各个网格节点的坐标及其函数值进行储存，并建立查找表，通过查表程序进行查表以进一步减少计算时间，若网格节点的函数值已被计算，则可从查找表中查询得到，否则对该节点进行计算并储存该点坐标和函数值。(4)对于非包含网格，则重复执行步骤(2)～(3)一到两次，即以非包含网格为对象，再次利用二维二分法对非包含网格进行再次划分，分为更小的子网格，然后在每个子网格的顶点处利用数值积分法求解判稳函数f(x)的函数值，判断是否为非包含网络，如果获得的子网格仍然为非包含网格，则结束；否则，获得新的包含网格，并转向步骤(5)。由于判断方法搜索次数少时会造成误判，若将包含网格判断为非包含网格，则会影响逼近的f(x)函数曲线的完整性，因此，本发明对非包含网格再次使用二维二分法迭代一到两次，若子网格仍全是非包含网格，则认为此网格确为非包含网格，并结束此子流程；否则，可获得新的包含网格，以保证逼近的f(x)函数曲线完整性。(5)将步骤(3)和(4)中获得的所有包含网格重复执行步骤(2)～(4)，逐步逼近f(x)曲线，直至达到所述预设的迭代次数。即再次利用二维二分法将步骤(3)和(4)中获得的所有包含网格划分为更小的子网格，然后以这些更小的子网络为对象，再次在每个子网格的顶点处利用数值积分法求解判稳函数f(x)的函数值，判断是否为非包含网络，对于非包含网格重复划分和计算，获得新包含网格，以逐步逼近f(x)曲线，直至达到预设的迭代次数(即重复步骤(2)～(4)的次数)或者所需的计算精度，也就是说对获得的包含网格反复使用二维二分法，不断细化网格，不断迭代直至满足迭代次数，获得以更高精度逼近f(x)函数曲线的包含网格。(6)将最后获得的所有包含网格(即经过预设的迭代次数后最终获得的所有包含网格)进行线性插值，得到近似的判稳函数f(x)的零点(即稳定性边界点)。即对最后一次迭代获取的包含网格进行分析，每个包含网格必定至少有两个顶点异号，在相互异号的两个顶点之间进行线性插值，可获得近似的判稳函数f(x)的零点，对全部包含网格进行线性插值后，便可绘出f(x)函数零点的散点图，即组成稳定性边界。获取稳定型边界后，在加工中可以根据稳定性图来选取合适的加工参数，既能避免颤振产生又能提高加工效率，以优化加工工艺参数。线性插值具体为：设包含网格的两异号顶点坐标及其函数值分别为(x1,y1,f1)，(x2,y2,f2)，由于函数值f1与f2异号，则此两顶点之间存在函数值为0的点，设其坐标及函数值为(x0,y0,0)，由线性插值：x0＝x1-f1·(x2-x1)/(f2-f1)＝(x1·f2-x2·f1)/(f2-f1)；y0＝y1-f1·(y2-y1)/(f2-f1)＝(y1·f2-y2·f1)/(f2-f1)；根据所述零点坐标绘制散点图，获得所需的稳定性边界。以下为本发明的具体实施例。实施例1该实施例为一个单自由度系统，按照本发明的方法进行铣削稳定性求解，其步骤为：(1)设定参数平面，参数平面为：x轴为铣削速度Ω(rpm)，上下边界为5000rpm～25000rpm，y轴为铣削深度ap(m)，上下边界为0～0.008m；函数f(x)为判稳矩阵的特征根的绝对值的最大值减去1；设定径向切深a/D＝0.50；设定迭代次数i＝5。其他相关参数：铣刀齿数N＝2；切向切削力系数Kt＝6×108(N/m2)；法向切削力系数Kn＝2×108(N/m2)；系统固有频率fn＝922(Hz)；相对阻尼ζ＝0.011；模型质量mt＝0.03993(kg)；数值积分法中初始离散数目n＝20；然后，初步划分整个参数平面，划分为4×1个网格，即在x轴(铣削速度)上等分为4段，每隔5000rpm划分一段，在y轴(铣削深度)上保持不变，此时的网格认为均为包含网格；(2)对所有包含网格分别使用二维二分法，每个网格被划分为4个更小的子网格；(3)在每个子网格的四个顶点(网格节点)处求解判稳函数f(x)，得到四个节点的函数值；对所有子网格进行分析，若一个网格的四个节点处的函数值有异号，则认为此网格是包含网格，储存此网格四个节点的坐标及函数值，建立查表程序以进行下次计算；(4)若在步骤(3)中有非包含网格，则对所有非包含网格进行分析，并存储此子流程中获得的新的包含网格的坐标及函数值，为平衡计算效率和计算精度，在本实例中，只对第一次迭代得到的非包含网格进行此子流程的运算。(5)对于步骤(3)和(4)中获得的所有包含网格，重复执行步骤(2)～(4)，逐步逼近f(x)的曲线，即逐步逼近铣削稳定性边界，直到达到设定的迭代次数i＝5；(6)对最后获得的所有包含网格进行线性插值，得到近似的稳定性边界点，并存储边界点坐标，对边界点进行处理并绘图，获得稳定性边界。获取稳定型边界后，在加工中可以根据稳定性图来选取合适的加工参数，既能避免颤振产生又能提高加工效率，以优化加工工艺参数。如果在稳定性散点图上方选取加工参数会引起颤振产生，影响加工后的表面质量，而在稳定性散点图下方选取加工参数可以避免颤振产生，但是在未知稳定性图的情况下，选取的加工参数往往偏小，导致加工效率低，因此根据本发明的方法在获取稳定性图后，可以选择合适的加工参数，以高速高效加工，同时能保证加工质量。将利用本发明获得的稳定性边界结果与与利用现有的数值积分法(NIM)获得结果进行比较，采用本发明的方法获得的稳定性边界图(lobes图)与数值积分法的对比图如图3所示，其中，图3(a)为本发明的结果图，图3(b)为数值积分法的结果图，与数值积分法相比，本发明的计算效率更高(时间大大缩短)，并且计算得到的结果基本一致，虽然使用二维二分法得到的图形会丢失一些细节(如尖峰，小孤岛)，但是对实际加工中加工工艺参数的选取没有影响。实例2该实施例为一个二自由度系统，按照本发明的方法进行铣削稳定性求解，其步骤为：(1)设定参数平面，参数平面为：x轴为铣削速度Ω(rpm)，上下边界为5000rpm～25000rpm；y轴为铣削深度ap(m)，上下边界为0～0.008m；函数f(x)为判稳矩阵的特征根的绝对值的最大值减去1；设定径向切深a/D＝0.50；设定迭代次数i＝4。其他参数：铣刀齿数N＝2；切向切削力系数Kt＝6×108(N/m2)；法向切削力系数Kn＝2×108(N/m2)；系统固有频率ωx＝ωy＝922(Hz)；相对阻尼ζx＝ζy＝0.011；模型质量mx＝my＝0.03993(kg)；数值积分法中初始离散数目n＝15；然后初步划分整个参数平面，划分为4×1个网格，即在x轴(铣削速度)上等分为4段，每隔5000rpm划分一段，在y轴(铣削深度)上保持不变，此时的网格认为均为包含网格；(2)对所有包含网格分别使用二维二分法，每个网格被划分为4个更小的子网格；(3)在每个子网格的四个顶点(网格节点)处求解判稳函数f(x)，得到四个节点的函数值；对所有子网格进行分析，若一个网格的四个节点处的函数值有异号，则认为此网格是包含网格，储存此网格四个节点的坐标及函数值，建立查表程序以进行下次计算；(4)若在步骤(3)中有非包含网格，则对所有非包含网格进行分析，并存储此子流程中获得的新的包含网格的坐标及其函数值，为平衡计算效率和计算精度，在本实例中，只对第一次迭代得到的非包含网格进行此子流程的运算；(5)对于步骤(3)和(4)中获得的所有包含网格，重复步骤(2)～(4)，逐步逼近f(x)的曲线，即逐步逼近铣削稳定性边界，直到达到设定的迭代次数i＝4。(6)对最后获得的所有包含网格进行线性插值，得到近似的稳定性边界点，并存储边界点坐标，对边界点进行处理并绘图。将利用本发明获得的稳定性边界结果与利用现有的数值积分法(NIM)获得结果进行比较，采用本发明的方法获得的稳定性边界图(lobes图)与数值积分法的对比图如图4所示，其中，图4(a)为本发明的结果图，图4(b)为数值积分法的结果图，与数值积分法相比，本发明的计算效率更高(时间大大缩短)，并且计算得到的结果基本一致，虽然使用二维二分法得到的图形会丢失一些细节(如尖峰，小孤岛)，但是对实际加工中加工工艺参数的选取没有影响。本领域的技术人员容易理解，以上所述仅为本发明的较佳实施例而已，并不用以限制本发明，凡在本发明的精神和原则之内所作的任何修改、等同替换和改进等，均应包含在本发明的保护范围之内。当前第1页1 2 3