本发明属于多轴数控加工的技术领域，涉及到三轴加工切屑几何与切削力集成建模的问题。
背景技术：
在制造业飞速发展的今天，数控加工技术不断应用于各种制造活动中。对于制造活动而言，零件的成形精度和制造质量是评价加工质量合格与否的标准，而数控仿真技术在这其中扮演了重要的角色。数控加工仿真技术的方式有几何仿真和物理仿真两种。几何仿真是机械加工的前提，它是对刀具运动的仿真，单从几何角度考虑，验证数控程序是否正确。而物理仿真是对加工本质的揭示，考虑切削参数和切削力等物理因素，对切削过程中的力、热、振动等进行分析与计算。我们可以通过数控仿真对切削参数进行优化，从而使切削过程更为高效。
在现有的研究中，更多的是单独从几何仿真或物理仿真的角度出发，并未将二者进行有效地集成。在几何仿真方面，模型精度和仿真速度的提高，以及模型的通用性，一直是研究的热点问题。而在物理仿真方面，现在的研究多以平均切削力的计算为研究内容，而不是瞬时切削力，从而无法很好地模拟动态的切削过程。
为了解决上述问题，本文将提出一种基于切屑几何模型和切削力模型的集成模型，达到同时表达加工中工件几何形状及未变形切屑形状，并预测加工中瞬时切削力的目的。
技术实现要素：
本发明需要解决的技术问题是：给出一种分层方法求解加工中未变形切屑形状和工件实时模型，并精确计算未变形切屑厚度，最后采用一种应用广泛的切削力预测模型进行加工中瞬时切削力的预测，解决三轴加工中几何仿真和物理仿真相结合的问题。
本发明的技术方案是：一种三轴加工切屑几何与切削力集成建模方法，其特征在于，包括以下步骤：
步骤一：刀具和工件建模
子步骤1)通用刀具建模，完成通用刀具截面形状的参数化表达；其中对于刀轴方向上某一高度l所对应的切削圆半径r(l)可由下式表达：
其中，l1＝lR-R2·cosα和l2＝lR-R2·sinβ；α表示下锥面与水平向的夹角、l1表示下锥面的高度、R1表示圆弧段圆心到刀轴间的距离、R2表示圆弧段的半径、lR表示圆弧段圆心到刀底的距离、l2表示圆弧段顶端到底端的距离、β表示上锥面与竖直方向的夹角。
子步骤2)工件建模，用分层思想对工件进行建模。即假定有一系列垂直于刀轴，且以ΔL为偏移距离的平行平面，用这些平面去切工件，将两个平面之间的工件称之为一层。用这些二维的层来表示三维工件。
步骤二：求解局部坐标系下切削圆表达式，包括以下子步骤：
子步骤1)：建立局部坐标系
局部坐标系到工件坐标系的坐标变换矩阵为：
其中，γ表示局部坐标系x轴与工件坐标系X之间的夹角、X1，Y1，Z1分别表示局部坐标系三个轴与工件坐标系三个轴各自的偏移距离。
子步骤2)在局部坐标系下，取工件的某层称之为层Π，根据层Π被刀具下锥面切削、层Π被刀具圆环面切削和层Π被刀具上锥面切削这三种情况，分别求解给定某层Π的切削圆半径；
第一种：当0≤zΠ≤l1且α＞0时，刀尖点位于z，层Π上的切削圆半径为：
第二种：当l1＜zΠ≤l2且α＞0时，刀尖点位于z，层Π上的切削圆半径为：
第三种：当l2＜zΠ且α＞0时，刀尖点位于z，层Π上的切削圆半径为：
其中，z表示刀尖点位置、zΠ表示层Π在局部坐标系下的z轴坐标。z2表示第二个刀位点在局部坐标系下的z坐标。
子步骤3)：根据子步骤2)中所得到的切削圆半径，求出刀尖点在局部坐标系下的切削圆表达式；
其中，θ∈[0,2π]，z∈[0,z2]，x2、z2分别表示在局部坐标系下第二个刀位点的坐标。
步骤三：根据步骤二中求得的切削圆表达式，结合现有包络理论，求出切削圆包络E：
步骤四：结合步骤三中切削圆包络的表达式，精确计算瞬时切削力和瞬时材料去除率，包括以下子步骤：
子步骤1)：计算在局部坐标系z轴每个刀步间的增量为：
其中，NR表示在一个刀步中刀具旋转的圈数、ns表示刀轴转速、fT表示进给率。
因此，第i个刀步所对应的的z值为：
zi＝Δz·i (9)
子步骤2)：在局部坐标系下求解切削角
φj(z,t)＝φ(t)+j·φp-ψ(z) (10)
其中，j＝1,2，…，NE，NE为刀具切削刃数、ψ(z)为滞后角函数、φp＝2π/NE为齿间角、φ(t)＝(2π·ns·t)/60为在加工中的t时，参考点的切削角。
子步骤3)：在局部坐标系下求解未变形切屑几何厚度
其中，φ啮合角、δ(φ)角度φ时层Π上的切屑厚度、ω为水平方向与切削表面法向方向的夹角。其可以表示为：
子步骤4)：采用线性切削力预测模型，计算切削盘上的切向力ΔFt,j，径向力ΔFr,j和轴向力ΔFa,j可由下式计算得到：
ΔFt,j(z)＝[ktc·hj(φj(z))+kte]·ΔL
ΔFr,j(z)＝[krc·hj(φj(z))+kre]·ΔL (13)
ΔFa,j(z)＝[kac·hj(φj(z))+kae]·ΔL
其中，hj(φj(z))为未变形切屑厚度，ΔL为工件层厚度，ktc,krc,kac分别为材料在切向、径向和轴向的的剪切系数，kte,kre,kae为对应的边缘系数；由此能够推导出层Π在x，y，z方向上的切削力：
ΔFx,j(z)＝-ΔFt,j·cos(φj(z))-(ΔFr,j·cos(ω(z))-ΔFa,j·sin(ω(z)))·sin(φj(z))
ΔFy,j(z)＝ΔFt,j·sin(φj(z))-(ΔFr,j·cos(ω(z))-ΔFa,j·sin(ω(z)))·cos(φj(z)) (14)
ΔFz,j(z)＝ΔFr,j·sin(ω(z))+ΔFa,j·cos(ω(z))
发明效果
本发明的技术效果在于：达到同时表达加工中工件几何形状及未变形切屑形状，并预测加工中瞬时切削力的目的。
附图说明
图1为通用刀具截面的参数化表达。
图2为一个刀步内工件材料的去除过程。
图3为层Π上工件初始二维轮廓。
图4为工件坐标系(X,Y,Z)与局部坐标系(x，y，z)的关系。
图5为局部坐标系下切削圆轮廓的表示。
图6为切削力模型原理。
图7为计算为变形切屑厚度图解。
图8左为加工的工件模型
图9左为工件处于位置1时切屑形状，右为预测的和测量的瞬时切削力。
图10左为工件处于位置2时切屑形状，右为预测的和测量的瞬时切削力。
图11左为工件处于位置3时切屑形状，右为预测的和测量的瞬时切削力。
图12为一个刀具旋转周期内的预测的合成力与实验值的误差。
具体实施方式
下面结合具体实施实例，对本发明技术方案进一步说明。
本发明首先给出一种分层方法求解加工中未变形切屑形状和工件实时模型，并精确计算未变形切屑厚度，然后采用一种应用广泛的切削力预测模型进行加工中瞬时切削力的预测。最后进行系统计算和加工实验，并做数据分析和比对。
步骤1：通用刀具的参数化表达，将加工中各类常见刀具的截面线统一划分为直线段和圆弧段，采用一种7参数表示法，给出对应参数的定义，然后根据不同刀具的截面形状设置不同的参数，从而完成通用刀具截面形状的统一表达。
步骤2：求解刀具加工过程中形成的包络和未变形切屑的几何形状，给出局部坐标系下切削圆的数学表达式，通过坐标变换，对应得到工件坐标系下的切削圆表示式。然后将工件模型进行分层，根据包络理论，求解出某一工件层对应的顶部平面上，一个刀步内一组切削圆的包络轮廓，从而获得该刀步切削后该平面上的工件二维轮廓以及该刀步内未变形切屑的二维轮廓。
步骤3：通过分层的建模方法将工件三维实体模型转化为二维轮廓，并计算出刀具沿刀具路径运动时形成的切入切出角，计算未变形切屑厚度。最后结合现有的切削力预测模型，计算获得加工中的瞬时切削力。具体相关步骤为：
步骤1：
刀具截面轮廓的表达及参数设定如图1所示。根据实际加工中所用的刀具，定义其切削刃数目NE和螺旋角ψ，将刀具的切削表面看做由无穷多个切削圆沿刀轴方向堆积而成，如图1所示。对于刀轴方向上某固定高度l，其所对应的切削圆半径r(l)可由下式表达：
其中，l1＝lR-R2·cosα和l2＝lR-R2·sinβ。
将α＝0,β＝0,R1＝5,R2＝0,l1＝0,l2＝10,lR＝10,带入式(1)完成刀具建模。
步骤2：
工件坐标系下切削圆表达式的求解。图2表示刀具以线性移动的方式从一个刀位点沿特定的进给方向移动至下一刀位点的过程，将它定义为一个刀步。为了更为直观的解释工件材料的去除过程，将工件模型沿Z方向进行分层，假定水平平面Π上初始二维形状如图3所示，由几条直线段和圆弧段组成，即{B1B2,B2B3,B3B4,B4B5,B5B6,B6B1}。将刀具的线性运动离散成若干个刀步，刀步的设定与刀具的旋转周数NR有关。NR的值根据具体加工而设定，以图2为例，其NR值为1。当加工开始时，切削圆1在平面Π上，平面Π在其内部的工件材料被去除，当刀具沿其刀轴旋转并沿进给方向移动时，切削圆1离开平面Π，在切削圆1下面的切削圆依次到达平面Π，并将其中的材料切除。在一个刀步中，所有切削圆依次连续的到达平面Π，可在该平面Π上得到一个由这些切削圆组成的二维包络轮廓E，在这个轮廓内的材料为该工步所去除的所有材料。将包络轮廓E与平面Π的最初形状求交可以得到该平面上未变形切屑几何的二维形状{B1C1,C1C2(E),C2B6,B6B1}及该工步结束后的工件二维形状{C1B2,B2B3,B3B4,B4B5,B5C2,C2C1(E)}。
在三轴加工中，刀具从刀位点CL1移动到点CL2，两点在工件坐标系(X,Y,Z)下的坐标分别为和在CL1处构建一个局部坐标系(x，y，z)。z轴与工件坐标系中的Z轴平行，将刀具移动的向量投影至XOY面，得到向量以点CL1为局部坐标系的坐标原点，x轴平行于投影向量x轴与z轴求向量积，并定义该方向为y轴。工件坐标系与局部坐标系的关系如图4所示。
投影与X轴之间的夹角为γ，则由局部坐标系到工件坐标系的坐标变换矩阵为：
那么，局部坐标系下的CL1及CL2点的坐标为：
在局部坐标系下，某一平面的切削圆轮廓可以很容易的通过计算得到，然后将其通过坐标变换矩阵M转到工件坐标系下进行模型更新。
为了得到局部坐标系下的切削圆轮廓的一般表达式，将水平平面Π的Z轴坐标ZΠ转到工件坐标系下：zΠ＝ZΠ-Z1。当刀具从CL1运动至CL2，刀尖的坐标是关于其z坐标的函数。层Π最顶层平面的切削圆圆心o和刀尖有同样的x坐标和y坐标，可以表示为：
在这个刀具运动的过程中，层Π上的切削圆半径同样是关于z坐标的函数，如图5所示，具体可分为以下三种情况：
第一种：当0≤zΠ≤l1且α＞0时2，刀尖点位于z，层Π上的切削圆半径为：
第二种：当l1＜zΠ≤l2且α＞0时，刀尖点位于z，层Π上的切削圆半径为：
第三种：当l2＜zΠ且α＞0时，刀尖点位于z，层Π上的切削圆半径为：
因此，当刀尖点在z，在局部坐标系下的相应的层Π上的切削圆表达式为：
其中，θ取值范围为[0,2π]，z的取值范围为[0,z2]。
当一组切削圆通过层Π的顶部平面时，一部分工件材料被去除掉。被去除掉的材料轮廓和加工后的工件轮廓均可通过该组切削圆的包络轮廓来得到。该组切削圆的包络轮廓需满足下面的条件：
通过上式可以得到参数θ和z的关系：
对于第三种情况，等式两边分别对z求导可得:
那么，包络轮廓E的表达式为：
将X1＝80,Y1＝25,Z1＝-7,X2＝52,Y2＝25,Z2＝-7带入两点在工件坐标系(X,Y,Z)下的坐标分别为和通过式(3)中旋转矩阵M求得上两点分别在局部坐标系下的坐标和为了得到局部坐标系下的切削圆轮廓的一般表达式，将水平平面Π的Z轴坐标ZΠ转到工件坐标系下：zΠ＝ZΠ-Z1。然后根据式(4)计算刀尖点在局部坐标系下的坐标oL(z)，根据zΠ的取值范围结合式(5)、(6)、(7)计算对应的切削圆半径，由式(8)获得层Π上的切削圆表达式最后根据包络理论并结合式(12)可求得层Π上的包络轮廓
步骤3：
1、为了精确计算瞬时切削力和瞬时材料去除率，将相邻两个刀位点间刀具的线性运动离散为若干个刀步。为了获得较高精度，将刀具从CL1运动至CL2，刀具在每一个刀步中旋转周数NR设为1。因此，在局部坐标系中每个工步间z的增量为：
因此，第i个刀步所对应的的z值为：
zi＝Δz·i (28)
2、在一个刀步中，刀具切过一组工件层，每层所对应的刀具切削层有相同的轴向高度z和厚度ΔL。因此，在每层上刀具与工件作用的切入角和切出角可以通过计算得到，并可以通过切入切出角判断在工步中参与切削的切削刃。刀步中刀具与工件的作用可以用图6表示。
某一时刻切削盘所在的位置不同，所以参与切削的切削刃及未变形切屑厚度、瞬时切削力也不同。因此，需要建立关于时间和工件层对应的z坐标的切削角函数。由于切削刃是螺旋状的，故相邻的两个切削盘的切削刃线段间有一个相位差，定义为滞后角ψ。在实际中，铣刀的切削刃有恒定的导程和螺旋角，滞后角ψ可以由此计算得出。在这里，定义滞后角函数为ψ(z)，包括所有可能参与切削的切削刃。
在任一时刻，切削盘上第j个切削刃的切削角φj可由滞后角函数ψ(z)和z坐标计算得到。对于一个有NE个切削刃的铣刀而言，齿间角为φp＝2π/NE。假设最底端的第一个切削刃在y轴上，在切削前其切削角为0。这个点可作为工步中计算任意切削刃上的点的切削角的参考点。在加工中的t时，参考点的切削角为φ(t)＝(2π·ns·t)/60。那么，对应z坐标为z的切削盘上第j个切削刃在时刻t所对应的切削角为：
φj(z,t)＝φ(t)+j·φp-ψ(z) (29)
3、在一个工序内，未变形切屑的形状对于计算该工序内的瞬时未变形切屑厚度十分重要。在该工序内，切削圆连续的去除层Π上的工件材料并形成未变形切屑。如图1-6所示，从一个切削圆中心o引一条射线，交切屑边界于A和B点。该射线与y轴的夹角φ称为浸入角，点A和点B间的距离定义为δ。当浸入角φ变化时，δ随之变化。图中的φs、φe为别为刀具——工件相作用的切入角和切出角。假定每个切削刃均匀地去除材料，且每个切削刃旋转一周产生一个切屑。所以每个工步的平均切屑厚度可以通过距离δ(φ)和切削刃数目NE计算得到。因为平均切屑厚度为水平方向，而不是实际加工中切削表面的法向方向，所以将平均切屑厚度沿法向方向投影，从而得到真实切削中的切屑厚度，进一步提高计算精度。第j个切削刃的未变形切屑几何厚度可以表示为：
其中，j＝1,2，…，NE，ω为水平方向与切削表面法向方向的夹角。其可以表示为：
4、采用线性切削力预测模型，切削盘上的切向力ΔFt,j，径向力ΔFr,j和轴向力ΔFa,j可由下式计算得到：
ΔFt,j(z)＝[ktc·hj(φj(z))+kte]·ΔL
ΔFr,j(z)＝[krc·hj(φj(z))+kre]·ΔL (32)
ΔFa,j(z)＝[kac·hj(φj(z))+kae]·ΔL
其中，hj(φj(z))为未变形切屑厚度，ΔL为工件层厚度，ktc,krc,kac分别为材料在切向、径向和轴向的的剪切系数，kte,kre,kae为对应的边缘系数。这些系数由特定的刀具和工件材料所决定。
可以推导出x，y，z方向上的切削力：
ΔFx,j(z)＝-ΔFt,j·cos(φj(z))-(ΔFr,j·cos(ω(z))-ΔFa,j·sin(ω(z)))·sin(φj(z))
ΔFy,j(z)＝ΔFt,j·sin(φj(z))-(ΔFr,j·cos(ω(z))-ΔFa,j·sin(ω(z)))·cos(φj(z)) (33)
ΔFz,j(z)＝ΔFr,j·sin(ω(z))+ΔFa,j·cos(ω(z))
将ns＝3000,ft＝30带入式(13)中计算每个工步间z的增量，将NE＝4带入φp＝2π/NE，可求得齿间角。通过式(15)求得对应z坐标为z的切削盘上第j个切削刃在时刻t所对应的切削角φj(z,t)，根据式(16)求得第j个切削刃的未变形切屑几何厚度hj(φ)，根据式(17)求水平方向与切削表面法向方向的夹角ω(z)，根据(18)可以推导出x，y，z方向上的切削力。
通过实验标定切削力系数，结果如表1所示。
表1切削力系数标定结果