基于误差分布特性的凹球面铣削轮廓误差预测

2022-12-05葛人杰王仁伟栾晓娜

计算机集成制造系统 2022年11期

葛人杰，张松+，王仁伟，栾晓娜

(1.山东大学机械工程学院高效洁净机械制造教育部重点实验室，山东济南 250061；2.山东大学机械工程国家级实验教学示范中心，山东济南 250061)

1 问题描述

凹形复杂曲面多用于光学镜片[1]，凹形模具[2]等重要零件的表面，是一种应用广泛，加工精度要求高的曲面。而凹球面作为凹形复杂曲面的一种简化曲面，研究其表面的轮廓误差分布规律对提高凹形曲面加工精度具有十分重要的意义。加工曲面的轮廓误差预测是控制曲面加工精度的重要内容。加工曲面的轮廓误差反映了实际加工形成的曲面与理论设计曲面的偏离程度，是衡量曲面加工精度的重要指标之一。而曲面的加工精度对曲面的力学性能和表面功能起到了决定性的作用[3]。三坐标测量仪(Coordinate Measuring Machine，CMM)由于其精确性常用于对加工后的表面进行测量，通过比较实测坐标与理论坐标的差值即可得到零件表面的轮廓误差。测点密度对铣削轮廓误差的测量精度具有直接影响，但是随着测点密度的增大，测量效率和测量成本也相应提高。空间克里金插值法通过考虑数据空间分布相关性，从而根据有限采样点的分布规律对未采样点进行最优无偏估计。因此，研究克里金插值预测轮廓误差对提高铣削轮廓误差的测量效率具有重要意义。

空间克里金插值是空间数据统计分析的一种方法，近些年来空间数据统计在工件表面误差分析方面取得了一定成果。地理学第一定律指出，所有事物或现象在空间上都是有联系的，但相距近的事物或现象之间的联系一般较相距远的事物或现象间的联系要紧密[4]。SURIANO等[5]指出在一定的空间间隔内测量加工表面的数据，其数据确实表现出很大的自相关性。这是因为在空间相近的加工区域之间，它们的加工工艺参数(切削速度、进给量、切削深度)是相近的，且刀具和工件材料的性质也是相近的[6]。YANG等[7]提出一种利用空间统计量估计形状误差的方法，通过对5种常用的加工表面形状误差的估计，验证了该方法的有效性。费兰等[8]运用泛克里金插值法，以少量数据点对发动机缸体表面的平面度误差进行估计，提出的方法比传统的最小二乘估计法精度提高了5%～10%。DU等[9]利用协同克里金插值方法，将刀具振动通过交叉变异函数考虑到对平面和曲面的形状误差预测中，验证了考虑加工条件的协同克里金插值法预测精度优于普通克里金插值法、反距离加权法和三角剖分插值法。YANG等[10]认为坐标测量仪测量的波纹度具有确定的平滑的趋势，而粗糙度部分是随机的，因此将形状误差分为确定和随机的部分，并利用误差随机部分的空间数据是不相关的这一特点，将形状误差区分开来，用非线性最小二乘法对确定性误差进行最优拟合。YAN等[11]基于空间自相关性的Moran’sI指数建立了可以将随机性误差分量从总的误差中分离出来的迭代算法，将该算法应用于球头铣刀加工的表面进行几何误差分析，提取出了确定性误差的拟合曲面。陈岳坪等[12-14]利用空间自相关分析对加工误差分解方面做了大量研究，他们提出加工误差的随机部分是空间不自相关的，从而分解出系统误差和随机误差，实现零件加工过程的系统误差补偿。然而，铣削曲面轮廓误差空间分布的各向异性却鲜有人关注。图1为铣削曲面常用的等高铣削方式。如图1所示，等高铣削中，沿切削进给方向上刀具的侧偏角一直不变；而沿横向进给方向上刀具的侧偏角在缓慢改变。这必然导致曲面铣削轮廓误差在这两个方向上分布的差异性。

因此，本文通过构建凹球面轮廓误差在不同方向上的变异函数，揭示了轮廓误差空间分布的各向异性。基于地统计学的基本假设建立了空间克里金插值模型，实现了对不同空间位置的轮廓误差进行预测。分析了插值范围对普通克里金、简单克里金和泛克里金插值方法预测精度的影响，得到了本问题中针对各种插值方法的最佳插值范围。最后，通过凹球面加工实验，验证了提出的克里金插值方法的有效性，并与线性插值预测方法进行比较，验证了提出的克里金插值预测方法的优异性。

2 凹球面铣削轮廓误差的空间分布特性

2.1 凹球面铣削实验设计及轮廓误差测量

为了分析铣削轮廓误差在不同方向上的分布特点，设计了如图2所示的凹球面加工实验，设计半径为31.8 mm，图中β为曲面倾角，β的加工范围在0°～75°之间。主轴转速恒定为3 000 r/min，每齿进给量0.1 mm，则进给速度为600 mm/min，采用等高轮廓顺铣的方式加工曲面，相邻刀具轨迹的倾角间隔为1°，沿球面的法向切深为0.6 mm。实验在三轴立式加工中心(MXR-460V，OKUMA-BYJC)上进行，为降低能耗和保护环境，实验采用不加切削液的干式切削方法。实验的工件为预硬塑料模具钢AISI P20；切削刀具为整体式球头立铣刀(JH970100-Tribon，Seco)，刀具半径5 mm，刀齿数为2，螺旋角为30°，刀具悬伸60 mm。

加工实验结束后，在三坐标测量仪(RAPID Plus，THOME)上对凹球面进行测量，图3为测量示意图，图中α为圆周角。以加工起点为0°,沿进给方向上每隔22.5°的圆周角对凹球面截面进行一次扫描读取测量点坐标，测量时每0.5 mm测量一个点。

图4a展示了实测点与理论设计曲面的关系，图4b展示了不同测量路径上各测量点的轮廓误差结果。表1列出了测量点的轮廓误差数据，其中1 074个数据用于轮廓误差的空间特性分析及预测模型的建立，剩下的100个带*标记的数据为随机挑选的轮廓误差数据，用于第二部分的轮廓误差预测的待估点。

表1 凹球面轮廓误差

为满足克里金插值的要求,需对凹球面轮廓误差正态化处理，处理后的统计数据如表2所示。正态化处理公式为：ln[e(Pi)+20]。由表2可知,Jarque-Bera统计量为0.079 6<5.991 5，说明在显著性为0.05的条件下，处理后的数据满足正态分布要求。

表2 正态化处理后凹球面的轮廓误差统计表

2.2 轮廓误差的空间统计描述方法

设三坐标测量仪采集铣削曲面轮廓的数据点总数为n，采样点的轮廓误差数据为e(Pi)(i=1，2，3…n)，Pi为采样点的空间位置坐标。在平稳加工过程中、相邻加工位置处的加工工艺参数(切削速度、进给速度、切削深度)是相近的，工件的物理化学性质、机床的状态、刀具的磨损情况等都是相近的。除非工件有疏松孔洞等缺陷，刀具突然的崩刃等加工条件突然改变的情况，否则在较小的空间距离h上的工件表面轮廓误差e(Pi)的增量[e(Pi)-e(Pi+h)]的数学期望为零[8]，即

E[e(Pi)-e(Pi+h)]=0，∀Pi,∀h。

(1)

进而，可用增量的方差函数作为空间关系的一种测量尺度，该方差函数只依赖增量之间的距离，不依赖其确切的空间位置[15]，因此有

D[e(Pi)-e(Pi+h)]=

E{[e(Pi)-e(Pi+h)]2}=2γ(h)。

(2)

式中：D为方差符号；E为数学期望符号；γ(h)为空间距离h步长上的半方差函数，也称变异函数，其只与步长h有关，与其所在的空间位置无关。由此可得增量的协方差函数与变异函数的关系如下。

Cov[e(Pi),e(Pi+h)]=σ2-γ(h)。

(3)

式中：Cov为协方差函数符号；σ2为h为0时的协方差，也称为先验方差。

式(1)与式(2)总称为本征假设，是地统计学中对随机函数的基本假设[16]。

变异函数可以反应数据空间相关性的变化情况，是空间分析的有效工具[17]。变异函数是根据变异函数实验值拟合得到的，而变异函数实验值是通过对实际测量结果计算获得的。为了获取变异函数实验值，首先设定一个分离距离Δh，设采样空间中存在N对数据点(P1,P1+jΔh)，(P2,P2+jΔh)，…(Pi,Pi+jΔh)，…(PN,PN+jΔh)，它们之间的距离d1，d2，…di，…dN在(j-1)Δh～jΔh之间，则变异函数实验值可由下式求得[18]：

(4)

2.3 凹球面轮廓误差空间分布特性

在凹球面等高铣削加工过程中，刀具沿圆周角方向做切削进给运动，沿着曲面倾角方向做横向进给运动。而同一条刀具路径上切削条件相同，不同刀具路径上刀具侧偏角不同。因此，曲面轮廓误差沿着曲面倾角和圆周角方向的分布规律是不同的，下面对曲面倾角和圆周角方向上轮廓误差的变异函数分别拟合建模，考虑这两个方向的各向异性。

为使变异函数满足本征假设，常用球状模型、指数模型和高斯模型等对实验获得的变异函数实验值进行拟合。其表达式如下：

球状模型：

指数模型：

高斯模型：

图5为曲面倾角方向轮廓误差的变异函数实验值以及分别采用球状模型、指数模型和高斯模型对变异函数进行拟合的结果，拟合的参数如表3所示。分离距离设置为1度。采用拟合的决定系数R2和均方根误差(RMSE)对球状模型、指数模型和高斯模型的拟合结果进行评价比较，表达式如下：

(5)

(6)

式中：r(hj)是球状模型、指数模型和高斯模型在hj处拟合的变异函数值；r*(hj)是实验获得的变异函数值；¯r*(hj)是实验获得的变异函数值的均值；M为实验获得的变异函数总数，此处为37。

由表3可知，球状模型的拟合决定系数R2=0.922 0，大于指数模型和高斯模型；且球状模型的均方根误差为1.897×10-3，小于指数模型和高斯模型。这说明球状模型的拟合精度高，拟合误差小。故本文采用球状模型描述沿曲面倾角方向上轮廓误差分布特点。

表3 理论拟合模型的参数

图6为圆周角方向轮廓误差的变异函数拟合模型，步长距离为22.5°，由图可知其变异函数随步长的变化不大，故采用常数模型拟合即γ(h)=0.019 34，该模型称为纯块金效应模型[16]，表明轮廓误差在圆周方向上是随机分布的，不具有空间相关性。一方面是由于采样间距大于轮廓误差在圆周方向上的空间变异尺度；另一方面是由于圆周方向为切削进给方向，在一次走刀过程中，加工参数几乎完全一样，因此沿圆周方向上的轮廓误差为一个定值加上随机部分，随机部分是考虑到刀具振动、工件的硬质点等可能产生的随机情况，正是随机部分导致了纯块金效应。

比较图5和图6可知，凹球面铣削轮廓误差的空间分布呈现各向异性。曲面倾角方向上轮廓误差的变异函数值随着测量点间距离的增加而增大，即空间相关性随着距离的增大而减小。当距离大于变程时，变异函数值稳定，此时轮廓误差不再具有空间相关性。而圆周角方向的轮廓误差分布呈现纯块金效应，表明沿该方向上的轮廓误差分布不具有空间相关性。

3 基于克里金插值法的轮廓误差预测

3.1 克里金插值方法

在许多情况下，区域化变量在研究区域内是非平稳的，其数学期望不是一个常数，即E[e(P)]=m(P)，m(P)称为漂移。曲面的轮廓误差也是符合这个条件的，因为在不同的曲面位置上刀具和工件的接触情况等并不完全一致，从而导致不同位置轮廓误差的数学期望随着空间位置的变化而变化，即产生漂移。在漂移存在的条件下，一般采用泛克里金插值方法对空间位置点进行插值估计[19]。漂移一般采用多项式表示：

(7)

式中：al为未知系数，l为多项式的阶数。

泛克里金插值预测的公式为：

(8)

式中：e*(P0)为待估点轮廓误差的预测值；e(Pi)为已知位置Pi(i=1,2,…n)处的已知轮廓误差，可简记为ei，Pi是通过设定插值范围来选择的，插值范围是指用于插值计算的数据点Pi离待估点P0的距离大小；λi为权重系数，可通过对预测点的无偏最优估计并结合拉格朗日乘数法计算得到，结果如下:

(9)

以矩阵的形式表示为:

(10)

式中:φ为拉格朗日乘子;矩阵表达式中各部分计算方法如下:

由式(3)中协方差函数与变异函数的关系，可将方程中的协方差函数用变异函数表示，进而求得权重系数λi。当漂移为常量漂移时，插值方法称为普通克里金插值；当普通克里金插值中常量漂移的a0已知，且预测公式为式(13)时，插值方法称为简单克里金插值[20]。

(11)

3.2 预测方法的实现

本文提出的凹球面轮廓误差预测流程如图7所示。通常方法是根据测量的误差数据用线性插值等方法直接预测未知点的轮廓误差，而本文所提方法通过分析轮廓误差的空间分布特性对轮廓误差进行预测。其步骤如下：

(1)对数据进行正态化处理，使其满足克里金插值的要求。

(2)计算实验的变异函数实验值并拟合，分析轮廓误差的空间分布特性。

(3)依次选择不同的克里金插值方法，并逐步增加插值范围，计算克里金方程式(10)中的权重系数，得到待估点的预测值。

(4)比较不同的插值方法在不同插值范围下的预测精度，选择最佳的插值方法及插值范围进行预测。

由于所提方法考虑了轮廓误差在采样区间的空间分布特点，使得预测的轮廓误差更符合空间变化规律，比通常的线性插值预测精度更高。

需要注意的是，插值范围是指用于插值估计的数据点离待估点的距离。当插值范围较小时，由于用于插值的数据点较少，对待测点的预测精度较低；随着插值范围的增加，预测精度会上升，但当用于插值的数据点距离待估点较远时，数据之间的相关性降低，可能会使预测精度降低，因此确定一个合适的插值范围对提高插值预测精度很关键[21]。

3.3 不同克里金插值模型预测结果对比分析

为验证克里金插值预测的有效性，并对比分析不同插值方法的预测精度，在所测量的数据点中随机选择100个作为待估点，并分为两个方向分别进行预测插值。曲面倾角方向采用普通克里金、简单克里金、一次漂移的泛克里金、二次漂移的泛克里金插值方式与线性插值方式对待估点的轮廓误差进行预测；圆周角方向采用同样的克里金方式与平均值预测方式(即同一曲面倾角其他圆周角处的均值作为待估点的预测值)对待估点的轮廓误差进行预测。考虑到距离越近，空间相关性越大，曲面倾角方向的插值范围取为2°，3°，4°，5°，6°，7°，8°，9°，10°，15°，30°，45°和全部数据点；圆周角方向的插值范围为45°，90°，135°，180°和全部数据点。并用均方误差(MSE)对预测的误差进行评价，表达式如下：

(12)

式中：epre和emea分别为待估点的预测值和测量值;m为待估点的样本总量；m=100。

图8为曲面倾角方向4种克里金插值方法在不同插值范围下预测的均方误差变化情况。由图8可知，普通克里金插值法在插值范围为2°时，预测的均方误差最小，之后随着插值范围的增加而缓慢增加，最后趋于平稳；简单克里金插值法预测的均方误差从2°～3°有一个明显的降低，之后随着插值范围的增加而缓慢减少，当插值范围大于6°时，预测效果均优于线性插值；一次漂移和二次漂移的泛克里金插值方法分别在插值范围为3°和4°时预测的均方误差最低，之后都随着插值范围的增加而缓慢增加，最后趋于平稳。由图8可知，克里金插值方法的预测效果整体优于线性插值。表4列出了各个克里金插值方法预测的均方误差最小的情况，由表可知，采用插值范围为2°的普通克里金插值法预测的均方误差最低为5.222，比线性插值预测的均方误差6.294降低了17%。

表4 各插值方法最优插值结果

图9为圆周角方向各个插值方法在不同插值范围下预测的均方误差对比情况。由图9可知，克里金插值预测效果明显劣于平均值预测方法，只有在插值范围为180°时，采用普通克里金和简单克里金插值方法得到的预测效果优于平均值预测方法，但预测效果提升的程度并不明显。对比图8，曲面倾角方向插值的最小预测均方误差为5.22，圆周角方向插值的最小预测均方误差为14.029，可以看出沿圆周方向的插值方法预测误差明显高于沿曲面倾角方向的预测误差。这是因为沿圆周角方向上的变异函数呈现纯块金效应，这表明沿圆周角方向的轮廓误差是完全不相关的，因此沿着圆周角方向进行插值预测的误差较高。这也解释了沿圆周角方向的克里金插值预测效果明显劣于平均值预测方法，因为克里金插值是基于数据间的空间相关性而建立的，当数据没有空间相关性时，其预测误差必然很高。