为提高机床定位精度，降低成本，需要对五轴机床各轴平衡度、垂直度和定位误差等进行精度优化设计，机床精度分配优化设计是决定五轴数控机床制造成本以及质量的重要因素之一. 传统的精度分配设计方法是通过标准与设计人员的经验来决定，近年来，有许多研究着眼于将智能算法应用于精度指标设计领域，于涛等^[1]将制造成本及质量损失设定为目标函数，基于粒子群算法和遗传算法分别求解装配公差，从而实现精度优化分配；熊珍琦等^[2]利用模拟退火算法，把工序公差与工序加工成本分别作为目标函数，获得工序公差优化结果，在很大程度上优化了产品加工精度；张岩等^[3]基于产品装配成本构建基于装配流程的公差优化模型，通过改进的粒子群优化（PSO）法求解制造成本和装配公差的多目标优化问题，并实现精度优化分配；匡兵等^[4]综合考虑成本−公差函数的公差优化分配方法，提出一种综合尺寸公差和形位公差，利用公差混合优化分配法和传统方法进行非线性优化分配求解的精度分配方法.

然而，现有研究通常将精度指标当作某个常量纳入模型进行研究，但在实际生产过程中，误差是符合正态分布的随机量. 作为常量的误差量与作为随机量的误差量经由尺寸链传递后，由于误差间相互耦合作用，最终计算得到的机床空间定位误差并不完全准确^[5]，该误差会影响精度指标分配结果.

针对这一问题，本研究将二阶矩误差模型引入精度分配模型，假设误差服从正态分布^[6]，以矩阵形式表达各误差项，从而消除由于忽略误差分布情况导致的分配误差.

1.   五轴数控机床误差模型建立

转台摆头型数控机床如图1所示. 基于多体运动学理论简化机床，由于各部分之间为串联结构，所以可将机床抽象为具有两条分支的树系统拓扑构型. 第1分支由床身、纵向滑台、转台构成，第2分支由床身、横向滑台、升降台和摆头构成.

图  1  C100P五轴机床结构

Figure  1.  Structure of C100P five axis machine tool

下载: 全尺寸图片 幻灯片

假设工件坐标系中工件坐标为${p_{\text{w}}}$，刀具坐标系中刀尖坐标为${p_{\rm{t}}}$，刀具与工件应重叠；${{\boldsymbol{T}}_{i({\rm{s}})}}$为刀具侧第i项平移矢量，${{\boldsymbol{T}}_{i({\rm{p}})}}$为刀具侧第i项旋转矢量，${{\boldsymbol{T}}_{j({\rm{s}})}}$为工件侧第j项平移矢量，${{\boldsymbol{T}}_{j({\rm{p}})}}$为工件侧第j项旋转矢量. 机床拓扑构型如图2所示.

图  2  机床拓扑构型

Figure  2.  Machine tool topology

下载: 全尺寸图片 幻灯片

假设工件坐标系共有n项空间变换，刀具坐标系共有m项空间变换，则工件坐标系与刀具坐标系之间的空间变换关系式为

刀尖点相对于工件对刀点的相对坐标为

根据多体系统坐标变换规则，空间变换矢量T可表示为

式中：$\mu \left( i \right)$为第i项坐标变换的均值；$\sigma \left( i \right)$为第i项坐标变换的方差；${\delta _i}$为i轴为旋转中心的角度变换；${\rho _i}$为沿i轴的平移变换.

2.   二阶矩误差模型定义

二阶矩误差模型本质上是各误差以某种分布形式表达，体现某一加工工序或尺寸精度指标的实际概率分布情况，从而实现对机床空间误差进行分配的一种误差建模方法. 根据泰勒展开式保留项的个数可分为一次二阶矩法^[7]、二次二阶矩法. 其中，一次二阶矩法能够表达离散误差变量的均值和方差特征，二次二阶矩法虽然加入了二次项，能够表达连续变量的变化率特征，但是由于误差测量得到的结果具有离散非连续的特点，二次二阶矩法不适用于机床精度分配. 因此，本研究选择一次二阶矩法定义二阶矩误差.

以Z表示误差测量值的集合，X表示误差数量，假设误差值符合正态随机分布，则其均值和方差可表示为

则该组误差可写成二阶矩形式为$Z = [\mu ,\sigma ]$. 特别地，当某误差为定值时，$\sigma = {\text{0}}$，误差二阶矩表示为$Z = [\mu ,{\text{0}}] = \mu$.

3.   二阶矩运算规则

应用二阶矩误差的齐次坐标变换法中，由于各矩阵中元素均以二阶矩形式表达，在运算过程中引入二阶矩阵运算，要想进行齐次坐标变换，需要定义基于二阶矩误差的矩阵运算规则.

二阶矩矩阵运算按照运算符号可分类为乘法与加法两类，根据运算符两端元素的类别可进一步分类为常量项与常量项、常量项与误差项、误差项与常量项、误差项与误差项等4类^[8]. 4种运算的数学表达式见表1.

表  1  二阶矩运算规则

Table  1.  Second-order moment operation rules

运算类型	二阶矩运算表达式
常数项与误差项加法	${c } + \Delta x = \left( {a,{b^2} } \right) + \left( {c,0} \right)$
常数项与误差项乘法	${c } {\text{•} } \Delta x = \left( {a,{b^2} } \right) {\text{•} } \left( {c,0} \right)$
误差项与误差项加法	$\Delta {x_1} + \Delta {x_2} = \left( {{a_{\text{1}}},{b_{\text{1}}}^2} \right) + \left( {{a_{\text{2}}},{b_{\text{2}}}^2} \right)$
误差项与误差项乘法	$\Delta {x_{\text{1} } } {\text{•} } \Delta {x_{\text{2} } } = \left( { {a_{\text{1} } },{b_{\text{1} } }^2} \right) {\text{•} } \left( { {a_{\text{2} } },{b_{\text{2} } }^2} \right)$

下载: 导出CSV 

| 显示表格

表中：$c$为常数项，$\Delta x$为满足正态分布的误差项，$a$为误差均值，$b$为误差标准差. 由表中运算表达式可知，常数项与误差项的加法和乘法运算本质上都是误差项与误差项运算在方差为0情况下的特例.

3.1   误差二阶矩与误差二阶矩加法运算

假设存在两个二阶矩阵${{x}_{\text{1}}}$和${{x}_2}$，误差项${{x}_{\text{1}}}$为${{{x}}_{\text{1}}} = \left( {{a_{\text{1}}},{b_{\text{1}}}^2} \right)$，表示${{x}_{\text{1}}}$满足均值为${a_{\text{1}}}$、标准差为${b_{\text{1}}}$的正态分布，误差项${{x}_2}$为${x_{\text{2}}} = \left( {{a_{\text{2}}},{b_{\text{2}}}^2} \right)$.

根据正态分布性质，两个独立正态分布相加得到的分布均值为两者均值之和，方差为两者方差之和，可表示为

3.2   误差二阶矩与误差二阶矩乘法运算

两个误差项为独立正态分布，即满足$\dfrac{{\phi ({x_{\text{1}}}) - {a_{\text{1}}}}}{{{b_{\text{1}}}}}\sim$$N\left( {0,1} \right),\dfrac{{\phi ({x_{\text{2}}}) - {a_{\text{2}}}}}{{{b_{\text{2}}}}}\sim N\left( {0,1} \right)$.

通过拟合正态分布可得，$N\left(\text{0}，\text{1}\right){\text{•}} N\left(\text{0}，\text{1}\right)$的均值为0，方差为1，表达式为

式（9）左侧为二次多项式形式，最高次项为$\phi ({x_{\text{1}}}) • \phi ({x_{\text{2}}})$且仅有一项，可得$\phi ({x_{\text{1}}}) • \phi ({x_{\text{2}}})$满足

根据式（7）、式（8）可得$\phi ({x_{\text{1}}}) {\text{•}} \phi ({x_{\text{2}}})\sim \Big( {a_1}{a_2}$，$\sqrt {{{\left( {{a_1}{b_2}} \right)}^2} + {{\left( {{a_2}{b_1}} \right)}^2} + {{\left( {{b_1}{b_2}} \right)}^2}} \Big)$

即$\phi ({x_{\text{1}}}) {\text{•}} \phi ({x_{\text{2}}})$的运算规则可表示为

如图3所示，以均值为0.1、方差为0.01的常量二阶矩阵${x_{\text{1}}} = \left( {{\text{0}}{\text{.1}},0.{\text{01}}} \right)$和均值为0、方差为0.04的二阶矩阵${x_{\text{2}}} = \left( {{\text{0}},0.{\text{04}}} \right)$为例，从${x_{\text{1}}}$和${x_{\text{2}}}$中各随机任取20组数值. 将两组数据相乘得到一组误差二阶矩与误差二阶矩的积，表示为

图  3  二阶矩乘法运算验证

Figure  3.  Verification of second-order moment multiplication

下载: 全尺寸图片 幻灯片

部分运算结果见表2.

表  2  误差二阶矩乘法运算部分结果

Table  2.  Partial results of error second-order moment multiplication

i	1	2	3	4	5	6	7
$\phi ({x_{\text{1}}}\left( i \right)){\text{/mm}}$	−0.0206	0.0509	0.0049	0.0588	−0.0062	−0.0202	−0.0083
$\phi ({x_{\text{2}}}\left( i \right)){\text{/mm}}$	0.0902	0.0976	0.0965	0.1049	0.1099	0.0986	0.1117
$\phi (y\left( i \right)){\text{/mm}}$	−0.0019	0.0050	0.0005	0.0062	−0.0007	−0.0020	−0.0009

下载: 导出CSV 

| 显示表格

$y\left( i \right)$的频率直方分布图和正态分布拟合曲线如图4所示. 以长方形的高度表示该数值区间出现的频率. 图4中$\phi ({C_{\text{1}}}\left( i \right))$和$\phi ({x_{\text{1}}}\left( i \right))$之和满足均值为0、方差为0.0044的正态分布，符合式（11）、式（12）的预测，证明了二阶矩阵常量项与误差项之和的测试集满足式（11）和式（12）.

图  4  $y\left( i \right)$频率直方分布及正态分布拟合曲线

Figure  4.  Fitting curve of frequency square distribution and normal distribution for$y\left( i \right)$

下载: 全尺寸图片 幻灯片

4.   机床精度指标−装配成本模型

以图1所示的C100P转台摆头式五轴数控机床为例，机床17项精度指标见表3.

表  3  C100P五轴数控机床精度指标

Table  3.  Accuracy index of C100P five axis machine tool

序号	变量表达式	变量含义	误差类型
1	$\Delta {\ell _X} = [{u_{\Delta {\ell _X}}},{o_{\Delta {\ell _X}}}]$	X直线轴轴线运动直线度	沿垂直X轴平动
2	$\Delta {\ell _Y} = [{u_{\Delta {\ell _Y}}},{o_{\Delta {\ell _Y}}}]$	Y直线轴轴线运动直线度	沿垂直Y轴平动
3	$\Delta {\ell _Z} = [{u_{\Delta {\ell _Z}}},{o_{\Delta {\ell _Z}}}]$	Z直线轴轴线运动直线度	沿垂直Z轴平动
4	$\Delta {\theta _{YZ}} = [{u_{\Delta {\theta _{YZ}}}},{o_{\Delta {\theta _{YZ}}}}]$	YZ直线轴垂直度误差	绕X轴转动
5	$\Delta {\theta _{XZ}} = [{u_{\Delta {\theta _{XZ}}}},{o_{\Delta {\theta _{XZ}}}}]$	XZ直线轴垂直度误差	绕Y轴转动
6	$\Delta {\theta _{XY}} = [{u_{\Delta {\theta _{XY}}}},{o_{\Delta {\theta _{XY}}}}]$	XY直线轴垂直度误差	绕Z轴转动
7	$\Delta {\theta _{AL}} = [{u_{\Delta {\theta _{AL}}}},{o_{\Delta {\theta _{AL}}}}]$	A轴与直线轴垂直度误差	绕Y，Z轴转动
8	$\Delta {\theta _{BL}} = [{u_{\Delta {\theta _{BL}}}},{o_{\Delta {\theta _{BL}}}}]$	B轴与直线轴垂直度误差	绕X，Z轴转动
9	$\Delta {\ell _{MZ}} = [{u_{\Delta {\ell _{MZ}}}},{o_{\Delta {\ell _{MZ}}}}]$	主轴轴线与Z轴平行度	沿X，Y轴平动
10	${T_{\text{p}}} = [{u_{{T_{\text{p}}}}},{o_{{T_{\text{p}}}}}]$	台面平面度	沿Z轴平动
11	$\Delta {\theta _{PX}} = [{u_{\Delta {\theta _{PX}}}},{o_{\Delta {\theta _{PX}}}}]$	台面与X轴运动平行度	绕Y，Z轴转动
12	$\Delta {\theta _{PY}} = [{u_{\Delta {\theta _{PY}}}},{o_{\Delta {\theta _{PY}}}}]$	台面与Y轴运动平行度	绕X，Z轴转动
13	${\varepsilon _X} = [{u_{{\varepsilon _X}}},{o_{{\varepsilon _X}}}]$	X直线轴定位精度	沿X轴平动
14	${\varepsilon _Y} = [{u_{{\varepsilon _Y}}},{o_{{\varepsilon _Y}}}]$	Y直线轴定位精度	沿Y轴平动
15	${\varepsilon _Z} = [{u_{{\varepsilon _Z}}},{o_{{\varepsilon _Z}}}]$	Z直线轴定位精度	沿Z轴平动
16	$\Delta {\theta _A} = [{u_{\Delta {\theta _A}}},{o_{\Delta {\theta _A}}}]$	A旋转轴定位精度	绕X轴转动
17	$\Delta {\theta _B} = [{u_{\Delta {\theta _B}}},{o_{\Delta {\theta _B}}}]$	B旋转轴定位精度	绕Y轴转动

下载: 导出CSV 

| 显示表格

定义机床参数：刀尖到A轴转动中心距离为$L = 400$；工件对刀点到转台转动中心沿X轴向距离为${D_X} = 0$；工件对刀点到转台中心沿Y轴向距离为${D_Y} = {\text{1}}0$；工件对刀点到转台中心沿Z轴向距离${D_Z} = 0$，直线轴每米定位精度定义为${\varepsilon _{{i}}}'(i = x,y,z)$，直线轴P在i方向上的轴角偏差表示为${\theta _{P{{i}}}}(P=X,Y,Z; i = X, Y, Z)$，3个直线轴运动量分别表示为X、Y、Z，转动轴运动量表示为A、B，旋转轴A、B沿X、Y、Z轴向偏移分别表示为${o_{ji}}(i = A,B; j = X,Y,Z)$，以上参数根据实际测量值作为定值纳入误差模型. 在五轴数控机床中涉及各轴直线度、垂直度等误差的控制，因此机床装配成本模型可表示为

式中：$C$为装配成本；${T_{i}}$为机床误差模型中第i项精度指标；${a_i}$为影响因子；$n$为精度参数数量；e为自然对数.

根据装配误差的类型、装配方式和装配条件等，不同精度指标的装配成本模型不尽相同^[9]. 为简便计算，本研究使用车削成本参数${a_0} = 0.96, {a_1} = 7.{\text{08}},{a_2} = - 2{\text{0}}.84$.

5.   实例分析

5.1   约束条件定义

在机床装配过程中，需要考虑两方面的约束，即装配能力约束和装配功能要求，在图1所示五轴数控机床中，各精度指标约束条件见表4.

表  4  各精度指标约束区间

Table  4.  Constraint interval of each accuracy index

序号	变量名	下极限偏差	上极限偏差
1	$\Delta {\ell _X}$	0.020 mm	0.030 mm
2	$\Delta {\ell _Y}$	0.020 mm	0.030 mm
3	$\Delta {\ell _Z}$	0.014 mm	0.015 mm
4	$\Delta {\theta _{YZ}}$	0.005/600 mm	0.025/600 mm
5	$\Delta {\theta _{XZ}}$	0.005/600 mm	0.025/600 mm
6	$\Delta {\theta _{XY}}$	0.005/600 mm	0.025/600 mm
7	$\Delta {\theta _{AL}}$	0.006 mm	0.010 mm
8	$\Delta {\theta _{BL}}$	0.006 mm	0.010 mm
9	$\Delta {\ell _{MZ}}$	0.004 mm	0.020 mm
10	${T_{\text{p}}}$	0.020 mm	0.020 mm
11	$\Delta {\theta _{PX}}$	0.025 mm	0.050 mm
12	$\Delta {\theta _{PY}}$	0.020 mm	0.030 mm
13	${\varepsilon _X}$	0.0064/1000 mm	0.015/1000 mm
14	${\varepsilon _Y}$	0.0118/1000 mm	0.015/1000 mm
15	${\varepsilon _Z}$	0.0073/1000 mm	0.010/1000 mm
16	$\Delta {\theta _A}$	7"	10"
17	$\Delta {\theta _B}$	4.6"	10"

下载: 导出CSV 

| 显示表格

根据式（14）构建误差传递模型分配机床精度指标，分配模型表示为

其中，目标函数设定为机床装配成本与机床空间定位精度，e为装配合格率，装配合格率指标${{e}_0}$取99.8%.

5.2   优化分配结果分析

应用快速非支配排序遗传算法进行多目标参数优化^[10]. 首先创建初始父代种群，每个决策变量在上下限间随机取值，取得设计种群数量的决策变量组. 然后将各组决策变量代入建立的二阶矩机床精度模型，取得装配成本与其对应的机床空间定位误差，进行非支配快速排序与拥挤度排序. 再对初始父代种群进行交叉变异，生成大小与父代相同的子代种群，合并两代种群后再次进行非支配快速排序和拥挤度排序，选取其中排名靠前的决策变量组，生成新的父代，循环直到结果收敛或达到循环次数上限gen，上限次数设定为200. 最终得到Pareto最优解曲线，如图5所示.

图  5  优化前后Pareto前沿对比

Figure  5.  Comparison of Pareto front before and after optimization

下载: 全尺寸图片 幻灯片

选取曲线上精度与成本均衡的一组解作为优化方案，具体分配数值见表5.

表  5  精度分配结果

Table  5.  Precision allocation results mm

变量名	优化前	优化后		变量名	优化前	优化后
$\Delta {\ell _X}$	0.021	0.022		$\Delta {\theta _{PX}}$	0.025	0.031
$\Delta {\ell _Y}$	0.022	0.024		$\Delta {\theta _{PY}}$	0.020	0.020
$\Delta {\ell _Z}$	0.014	0.015		${\varepsilon _X}$	0.010/1000	0.015/1000
$\Delta {\theta _{YZ}}$	0.020/600	0.003/600		${\varepsilon _Y}$	0.010/1000	0.013/1000
$\Delta {\theta _{XZ}}$	0.020/600	0.004/600		${\varepsilon _Z}$	0.010/1000	0.009/1000
$\Delta {\theta _{XY}}$	0.030/600	0.023/600		$\Delta {\theta _A}$	10"	7.8"
$\Delta {\theta _{AL}}$	0.006	0.007		$\Delta {\theta _B}$	10"	7.9"
$\Delta {\theta _{BL}}$	0.006	0.007		${p_{{t}}}$	0.055	0.055
$\Delta {\ell _{MZ}}$	0.010	0.008		$C$	161.97元	146.71元
${T_{\text{p}}}$	0.020	0.020

下载: 导出CSV 

| 显示表格

分析表5中精度优化分配结果，Pereto曲线上机床空间定位精度为0.055 mm的点较为均衡，因此选取该点精度分配结果进行对比分析，将曲线上定位精度为0.055 mm的点对应的各项精度值代入式（14），得到装配成本为146.71元. 对照组精度指标不考虑误差分布，精度指标以定值形式代入优化模型，在相同约束条件下选取机床空间定位精度相近的解，代入式（14）得到机床空间定位精度为0.055 mm时，装配成本为161.97元. 比较可知，本研究的优化模型在保证精度的前提下，降低了机床装配成本.

6.   结　语

利用二阶矩误差模型对五轴数控机床精度指标进行优化分配设计研究，得到如下结论.

1）结合五轴数控机床结构特点，定义五轴数控机床关键精度指标的二阶矩误差模型，分析二阶矩误差的各类运算规则，根据误差二阶矩频率直方分布情况，对误差二阶矩的运算规则进行了验证.

2）结合文中的精度指标优化分配模型，以C100P五轴数控机床为例，根据各项精度指标的装配成本模型和约束条件，求解五轴数控机床精度分配模型，得到各项精度指标的最优解集.

3）对优化前后五轴数控机床的装配成本与空间定位精度进行对比，证明了基于二阶矩误差模型的精度分配优化算法在保证机床精度的情况下，降低了机床装配成本.