机床的加工精度主要受导轨误差的影响：导轨误差映射到运动副上，运动副误差再传递至末端刀具，最终造成工件误差。研究误差传递过程有助于在机床设计阶段确定导轨刮研曲线，提高机床的加工精度。导轨误差与运动副误差的映射模型为过约束条件下的非线性模型，导轨与滑块的刚度为非线性，且滑块固定件受力存在过约束特性，模型较为复杂。对于磨削、铣削等刀具以体的形式参与的加工方式，无法通过刀尖点轨迹直接生成刀具误差与工件误差的映射模型。针对机床加工误差的传递过程，相关学者开展了大量研究。

孙健利^[1]与孙伟等^[2]应用赫兹接触理论，分析了滚动直线导轨副的力学性能，并建立了其沿垂直方向的刚度模型。Shaw等^[3]分析了外载荷及预紧力对直线导轨接触刚度的影响。仇进生等^[4]推导了滚动直线导轨副在转动力矩作用下的静刚度模型。倪国林等^[5]与张巍等^[6]对整个滚动直线导轨副的力学性能进行了分析，实现了任意外载荷下滚动直线导轨副的五维静刚度建模。Hung^[7]、Jastrzębski等^[8]与Pawełko等^[9]将滚动体等效为单向移动的弹簧单元，采用弹簧单元法对滚动直线导轨副结合部进行了有限元建模。Wu等^[10]与鲁明^[11]利用假想材料法模拟滚柱与滚道面结合部的接触特性，并通过设置虚拟材料的弹性模量、泊松比等参数来等效滚动体，从而建立了滚动直线导轨副的有限元模型。Park等^[12]为探究导轨直线度误差与工作台运动误差间的映射关系，运用全解析法构建了油膜作用力与工作台运动误差间的传递模型，并定义了误差传递函数。Ni等^[13]、孙伟等^[14]与孙光明等^[15]采用传递函数来描述导轨直线度误差与导轨承载力之间的关系并分析了其特性，由此量化了滚动直线导轨系统的误差均化效应。Rahman等^[16]与刘又午等^[17]以多体系统理论与机床拓扑结构为基础，构建了多轴数控机床的几何误差模型。陈国达等^[18]与Yu等^[19]通过计算工件加工表面刀具接触点的分布得到了加工表面轮廓，并以此为基础开展了加工精度预测。

然而，关于考虑滑块固定件结构刚度与导轨副非线性刚度的导轨误差与运动副误差间的映射模型以及铣削、磨削工况下刀具误差与工件误差间的映射模型尚未得到深入研究，且已有研究大多针对双导轨四滑块系统，无法推广到多导轨多滑块系统。为解决上述问题，本文通过建立单滑块系统等效刚度模型、多滑块系统误差映射模型及工件加工误差预测模型，综合分析在考虑滑块固定件结构刚度与导轨副非线性刚度条件下多滑块系统中导轨误差与运动副误差间的映射关系、运动副误差与末端刀具误差间的映射关系以及刀具误差与工件误差间的映射关系，最终得到导轨误差与工件误差间的映射模型，旨在为精密机床的设计提供指导。

当滚动导轨副产生误差时，滑块会偏离理想的运动轨迹，导致内力发生变化。为重新达到平衡状态，滑块发生偏摆，进而产生位姿误差。在进行单滑块系统静力分析时，作如下假设：

1）误差仅存在于导轨安装基面，导轨的装配误差及导轨基面误差均等效为导轨滚道面误差。

2）滚柱与滑块裙部均视为线性弹性体，导轨视为刚体。

3）滑块裙部在各滚动体所处位置的变形相同。

滚动导轨副由导轨、滑块与滚动体构成。滚动体在预紧力与外部载荷的作用下产生变形，其与导轨滚道面、滑块滚道面间的接触符合赫兹接触理论。基于赫兹接触理论，圆柱滚动体的变形量δ与所受载荷Q满足以下关系^[20]：

式中：v为滚动体材料的泊松比，E为滚动体材料的弹性模量，l为滚动体长度。

化简式(1)，可得：

其中：

滚动导轨副中单个滑块的尺寸及滚动体的排列情况如图1所示（坐标原点O_s位于滑块顶部平面中心）。图中：j为滚动体所在列，j=1, 2, …, m；i为滚动体所在行，i=1, 2, …, 4；L为滑块长度；D为2行滚动体间水平距离的一半；H为滚动体与滑块顶部平面间的垂直距离。当滑块存在相对位姿误差 e_s=[δ_x δ_y δ_z ε_x ε_y ε_z ]^T时，第i行第j列滚动体对应的滚道面位置处会产生沿水平方向（x向）和垂直方向（y向）的误差分量，分别为e_{x, ji} 和e_{y, ji}。

令第j列滚动体对应的滚道面在z向上与坐标原点O_s的距离为L_j，则有：

根据式(3)，误差分量e_{x, ji} 、e_{y, ji} 可表示为：

以第j列滚动体为例，对单滑块导轨副的静力平衡状态进行受力分析，如图2所示。假设滚动体的预压缩量为δ₀，初始接触角为α，第i行滚道面处滑块裙部的弹性系数为E_hi。

由图2可知，导轨与滚动体间、滚动体与滑块间各存在一对相互作用力，致使导轨、滚动体与滑块均产生变形，最终滑块变形与滚动体变形互相协调，使得导轨副达到静力平衡状态。令第i行第j列滚动体的变形量为δ_ji，第i行滚动体处滑块裙部滚道面的变形量为δ_i，第i行第j列滚动体的弹性力为F_{e, ji}。对于第1行第j列滚动体（初始压缩量为δ₀），当滚道面在x向上存在误差e_{x, j1}时，滚动体压缩量减少e_{x, j1}cos α；当滚道面在y向上存在误差e_{y, j1}时，滚动体压缩量增加e_{y, j1}sin α；当滑块裙部滚道面变形量为δ₁时，滚动体压缩量减少δ₁，则可得以下静力平衡方程：

同理，分析其余3行滚动体的受力情况，可得相应的静力平衡方程：

根据方程组(6)至(9)，可求解得到滑块裙部滚道面的变形量δ₁、δ₂、δ₃、δ₄，进而计算得到各滚动体的变形量δ_j1、δ_j2、δ_j3、δ_j4及弹性力F_{e, j1}、F_{e, j2}、F_{e, j3}、F_{e, j4}。

将滚动体对滑块的作用力分解，得到滚动体对滑块的水平作用力与垂直作用力F_{ex, ji} 、F_{ey, ji} ：

将第i行第j列滚动体沿x、y向的分力F_{ex, ji} 、F_{ey, ji} 等效作用在滑块顶部平面中心，得到x向等效力F_{x, ji} 、y向等效力F_{y, ji} 、x向等效力矩M_{x, ji} 、y向等效力矩M_{y, ji} 、z向等效力矩M_{z, ji}，分别表示为：

滚动体对滑块的作用力在滑块顶部平面中心的等效合力与滑块所受的外载荷平衡，则可得滑块所受外载荷 F 沿x、y、z向的分力F_x 、F_y 、F_z 和分力矩M_x 、M_y 、M_z，分别表示为：

由此可得到滑块相对位姿误差 e_s=[δ_x δ_y δ_z ε_x ε_y ε_z ]^T与滑块所受外载荷 F =[F_x F_y F_z M_x M_y M_z ]^T之间的映射关系。

第i行第j列滚动体的弹性势能可表示为：

滑块裙部在第i行滚动体处的弹性势能可表示为：

故单滑块导轨副的总势能W可表示为：

总势能W为关于变量δ_x 、δ_y 、δ_z 、ε_x 、ε_y 、ε_z 的函数，记为f (δ_x, δ_y, δ_z, ε_x, ε_y, ε_z )，将其分解后可得：

式中：W_δx 、W_δy 、W_δz 、W_εx 、W_εy 、W_εz 分别为由δ_x 、δ_y 、δ_z 、ε_x 、ε_y 、ε_z 引起的势能。

令由δ_x 引起的势能W_δx 表示为：

式中：K_δx 为单滑块导轨副对应δ_x 的等效刚度，K_δx 同样为关于变量δ_x 、δ_y 、δ_z 、ε_x 、ε_y 、ε_z 的函数。

同理，可得到单滑块导轨副对应δ_y 、δ_z 、ε_x 、ε_y 、ε_z 的等效刚度K_δy 、K_δz 、K_εx 、K_εy 、K_εz。令 K_s=diag(K_δx, K_δy, K_δz, K_εx, K_εy, K_εz )，则单滑块导轨副的总势能W可表示为：

由此可得到单滑块导轨副对应相对位姿误差 e_s的等效六维刚度 K_s。

为验证上述等效刚度模型的准确性，对单滑块导轨副进行仿真验证，其仿真模型如图3所示。将滚动体设置为非线性弹簧单元，利用ANSYS软件进行仿真，固定导轨基面，分别给滑块施加不同方向的力或力矩，获得其位移和转角，计算单滑块导轨副的等效刚度并与理论计算值进行对比，结果如图4所示。

由图4可以看出，单滑块导轨副的x向、y正向、y负向刚度及z向扭转刚度的理论计算值与仿真值基本一致，而x向、y向扭转刚度的理论计算值与仿真值存在差异但相差较小，这是由理论建模时假设滑块裙部在各滚动体所处位置的变形相同所引起的。

导轨误差传递到运动副的过程并不是误差直接叠加的过程，在传递过程中误差有所衰减，即存在误差均化效应。多滑块系统的结构如图5所示，其中滑块数量为n。导轨内部力的存在导致滚动体、滑块、导轨和滑块固定件产生弹性变形，抵消了部分导轨副误差，从而产生误差均化效应。如图6所示，在滑块固定件上取远离滑块的一点为原点O，根据机床坐标系方向建立坐标系A（O-xyz），在x、y、z轴上分别取点M、N、P。当多滑块系统受到外载荷作用后，滑块固定件产生形变。根据圣维南定理，O点附近的应力基本为零，即变形后原点O'与M'、N'、P'仍可构成正交坐标系，故以O'M'、O'N'、O'P'为坐标轴建立坐标系B（O'-x'y'z'），以分析导轨误差与运动副误差间的映射关系。

如图7所示，在滑块k处施加x向单位力F_{t, xk}，则滑块固定件各滑块位置的变形量为6n维向量，可表示为：

其中：

式中： e_{Ft, s} （s=1, 2, …, n）为单位力作用下滑块固定件在滑块s处的六维变形量。

同理，在滑块k处施加y向单位力、z向单位力、x向单位力矩、y向单位力矩、z向单位力矩，可得到对应的滑块固定件各滑块位置的变形量 e_{Ft, yk} 、 e_{Ft, zk} 、 e_{Mt, xk} 、 e_{Mt, yk} 、 e_{Mt, zk}。

当滑块k处受到作用力 F_k =[F_xk F_yk F_zk M_xk M_yk M_zk ]^T时，由于线性系统具有叠加性，滑块固定件各滑块位置的变形量 e_Fk 可表示为：

其中：

E_Fk =[ e_{Ft, xk} e_{Ft, yk} e_{Ft, zk} e_{Mt, xk} e_{Mt, yk} e_{Mt, zk} ]

式中： E_Fk 为作用力 F_k 下滑块固定件的柔度矩阵，为6n行6列矩阵。

由此可知，当滑块固定件受到作用力 F_b=[ F₁F₂ … F_n ]时，滑块固定件各滑块位置的总变形量 e_b为：

式中： E_b为滑块固定件的柔度矩阵。

根据式(28)，计算得到滑块固定件柔度矩阵 E_b：

则滑块固定件的刚度矩阵 K_b可表示为：

滑块固定件柔度矩阵 E_b为其各滑块位置在单位力作用下的变形量，故固定点的选择会影响刚度矩阵 K_b。滑块固定件受力 F_b为一组平衡力，当滑块固定件在 F_b的作用下达到静力平衡后，由圣维南定理可知，在远离载荷作用的区域内，应力基本上只与所受载荷的合力和合力矩有关，即应力为0。在此区域内任意选择2个面作为固定面，由于2种情况下受力一致，最终滑块固定件的变形势能相等。由此可知，只要在远离载荷作用的区域选取固定面，滑块固定件势能就与固定面的选取无关。

同理，可计算得到导轨刚度 K_g。

在外力作用下，导轨与滑块的变形如图8所示。

由图8可知，导轨的初始误差为：

式中： u₀为导轨初始误差，为6n维向量；d_xk 、d_yk 、d_zk，e_xk 、e_yk 、e_zk 分别为导轨在滑块k位置处的位置误差和转角误差。

当导轨副存在误差时，为了达到平衡状态，导轨产生变形 e_g，滑块相对导轨偏移 e_s。令平衡后的导轨误差为 u_g，滑块误差为 u_b，滑块固定件变形量为 e_b，其中 e_g、 e_s、 u_g、 u_b、 e_b均为6n维向量，且满足以下关系：

多滑块导轨副等效刚度矩阵 K_s、滑块固定件刚度矩阵 K_b及导轨刚度矩阵 K_g均为6n阶方阵，当多滑块系统达到静力平衡状态后，运动副的六维误差 v =[δ_x δ_y δ_zε_xε_yε_z ]^T。

则坐标系A相对于坐标系B的坐标变换矩阵 T_AB可表示为：

设滑块1在坐标系A下的位置矩阵 H_A1（用齐次坐标表示）为：

则滑块1在坐标系B下的位置矩阵 H_B1可表示为：

忽略滑块固定件变形时，滑块1的位移∆H₁可表示为：

同理，可得其余滑块的各向位移。

当滑块误差为 u_b，运动副误差为 v 时，根据忽略滑块固定件变形时各滑块相对理想位置的误差 Av，可得滑块固定件的变形量 e_b：

其中：

根据滑块固定件的刚度矩阵 K_b，计算滑块固定件受力时所做的功，即可得滑块固定件的弹性势能W_b：

根据最小势能原理^[21]， v 的取值应令W_b最小，即任取 v′= v+Δ v 时，ΔW_b=W′_b-W_b应恒大于等于0，而ΔW_b可表示为：

式(39)的后半部分可看作滑块固定件各滑块位置产生变形量 A Δ v 时所需的能量，恒大于等于0，故令式(39)的前半部分恒等于0即可满足要求，由此可得：

此时，滑块固定件的变形量可表示为：

其中：

式中： E 为单位矩阵。

将式(41)代入式(38)，可得滑块固定件的弹性势能W_b：

其中：

同理，多滑块导轨副的弹性势能W_s可表示为：

已知导轨变形量为 e_g，根据导轨刚度 K_g计算得到导轨受力为 K_ge_g，该力所做的功即为导轨的弹性势能，可表示为：

由此可得，多滑块系统的总势能W_t为：

根据最小势能原理，滑块误差 u_b与导轨变形量 e_g的取值应令多滑块系统的总势能W_t最小，即任取 u_b′= u_b+Δ u_b， e_g′= e_g+Δ e_g时，势能变化量ΔW_t=W_t′-W_t应恒大于等于0，其中ΔW_t可表示为：

令：

其中：

由于ΔW₂可看作刚度为 K_s、 K_g、 K_B的滑块固定件分别在各滑块位置产生变形量（Δ e_g-Δ u_b）、Δ e_g、Δ u_b时的弹性势能，故ΔW₂恒大于等于0，因此令ΔW₁恒等于0即可满足滑块误差 u_b与导轨变形量 e_g的取值要求，使得多滑块系统的总势能W_t最小，即：

由式(48)可解得：

将式(49)和式(50)代入式(40)，即可计算得到 v，从而得到在导轨副等效刚度为线性情况下的导轨误差与滑块固定件误差间的映射模型。当导轨副等效刚度为非线性时，根据式(32)，可得：

式(51)为关于滑块相对导轨的偏移量 e_s的多元方程，通过解方程可求得 e_s，进而得到考虑导轨副等效刚度为非线性时的导轨误差与运动副误差间的映射模型。

为验证上文构建的误差映射模型的准确性，以八滑块系统为例进行有限元仿真。将运动副误差的理论计算结果与仿真结果（包含忽略和不忽略滑块固定件结构刚度两种情况）进行对比，如图9所示。

分析图9数据发现，基于本文方法计算得到的运动副的x、y向位置误差和x、y、z向转角误差与原始仿真值（不忽略滑块固定件结构刚度）的相对误差最高分别为3.51%、2.46%、3.64%、11.90%、2.74%，而忽略滑块固定件结构刚度的仿真值与原始仿真值的相对误差最高分别为0.75%、17.90%、2.59%、57.30%、3.08%。由此可见，相比于传统的忽略滑块固定件弹性变形的方法，本文所提出的建模方法更贴近真实结果。

为方便计算工件的加工范围，需在工件坐标系下进行分析，因此需将末端刀具坐标系转化为工件坐标系。常见机床的拓扑结构如图10所示，其中：工件高度为p，刀具高度为q，刀具半径为R。

根据多体系统理论，刀具坐标系到工件坐标系的齐次坐标变换矩阵 T_wt可表示为：

式中： T_1w、 T_1t分别为机床坐标系到工件坐标系、刀具坐标系的齐次坐标变换矩阵。

则刀具坐标系下的齐次坐标 C_t与工件坐标系下的齐次坐标 C_w的变换公式可表示为：

刀具坐标系下刀具底面圆心位置向量 C_{t, c}与刀轴矢量 C_{t, a}可分别表示为：

联立式(52)和式(53)，可得工件坐标系下刀具底面圆心位置向量 C_{w, c}与刀轴矢量 C_{w, a}：

刀具上的离散点满足点与刀轴矢量的距离不超过R的条件，则有：

式中： D_{w, c}、 D_{w, a}分别为以笛卡儿坐标表示的工件坐标系下的刀具底面圆心位置向量与刀轴矢量，由式(55)计算获得。

此外，刀具上的点在刀轴矢量上的投影位于[0，q]内，则有：

以dr为间隔将工件离散为n_r 个圆柱面，以半径r的圆柱面为例，将柱坐标转换为笛卡儿坐标后代入式(56)与式(57)，可得工件加工区域内离散点的表达式：

假设工件加工时间为T，将时间离散为n_t 段；将圆柱面按角度离散为n_θ 段后，再将每段按高度离散为n_z 个点，并将各点赋予初值1，得到n_θ ×n_z 阶初始加工矩阵，如图11所示。将各点坐标与刀具位姿误差代入式(58)进行计算。若点在被加工区域内，则对应的矩阵元素赋值0，得到此时的加工矩阵。其中：1代表该位置有实体，0代表该位置被加工去除。经过n_t 次重新赋值后，可得到圆柱面的最终加工矩阵，随后将矩阵内值为1的元素转化回柱坐标系，即可得到圆柱面经加工去除后的剩余实体点。将n_r 个圆柱面的剩余实体点求并，即可得到最终工件实体点的分布情况。

利用MATLAB软件进行拟合计算，可得到根据给定导轨误差所加工的工件表面轮廓，具体流程如图12所示，并以此为基础计算工件的尺寸误差与圆柱度误差。

以某立式坐标磨床为例，其结构如图13所示。

考虑到尺寸误差可通过后处理补偿，故本文重点研究圆度误差与圆柱度误差。在实际加工过程中，导轨误差除受到导轨基面误差的影响外，还受到由静力作用导致的静力变形与由加工过程生热导致的热变形的影响^[22]。利用激光自准直仪测量导轨基础误差，并通过有限元仿真得到导轨的静力变形量与热变形量，将三者叠加，得到导轨的综合误差。随后，利用本文所提出的方法对工件表面轮廓进行预测，结果如图14所示。采用最小二乘法计算工件的圆柱度误差，为0.360 μm。

随后，采用最小二乘法分析工件上8、33、58、83、108 mm高度处圆形截面的圆度误差，结果如图15所示。

除导轨误差外，转台的跳动误差也是影响工件圆度与圆柱度的重要因素，而主轴转速远高于工作台转速，可看作其在跳动误差影响下的包络体参与加工，故主轴的跳动误差对工件的圆度误差与圆柱度误差无影响。采用电涡流传感器测量转台的径向跳动误差，并将转台径向跳动误差与导轨误差共同代入计算，得到综合误差影响下的工件表面轮廓，如图16所示。采用最小二乘法计算工件的圆柱度误差，为0.525 μm。

随后，采用最小二乘法分析工件上8、33、58、83、108 mm高度处圆形截面的圆度误差，结果如图17所示。

对比分析图14、图15与图16、图17，发现导轨误差对工件各截面处的圆度误差影响较小，导轨误差主要影响工件不同截面处的尺寸误差，进而影响工件的圆柱度误差，而转台跳动误差显著影响工件的圆度误差，进而增大了工件的圆柱度误差。

1）通过对滚动导轨副进行受力分析，得到考虑滑块裙部变形与滚动体非线性刚度的导轨副等效刚度模型。根据最小势能原理，对滑块固定件结构刚度进行分析，得到考虑滑块固定件结构刚度与导轨副非线性刚度时导轨误差与运动副误差间的映射模型。与传统忽略滑块固定件结构刚度的映射模型相比，本文方法所计算的运动副误差更精确。

2）将工件离散化后，通过布尔运算得到刀具误差与工件误差间的映射关系，最终建立导轨误差与工件误差间的映射模型，可在设计阶段为导轨刮研曲线的确定提供指导，从而提高机床设计效率。

3）通过实例分析了导轨误差以及导轨误差与转台误差综合影响下的工件误差，对比发现导轨误差与转台误差对工件圆度误差与圆柱度误差的影响不同。实例结果表明，本文所提出的方法具有普适性，适用于多导轨多滑块系统。