并联机构的各零件在加工、使用及生产过程中由于误差的存在，导致连接各构件的关节处存在间隙，当机构工作时，运动关节处的磨损会产生，进而引起系统振动，稳定性降低。基于上述，关节元素间必然存在的关节间隙已成为国内外学者研究的热点问题^[1-11]。

Koshy等^[12]基于Hertz理论对一种平面连杆机构进行了转动关节处含间隙的动力学建模与实验研究。Bai等^[13-14]利用混合非线性接触力模型对一种含间隙的双轴定位机构进行了动力学建模与数值分析。Marques等^[15]提出一种具有径向和轴向间隙的空间转动关节数学模型，并通过实例对模型进行了验证。Luka等^[16]对一种含间隙的滑槽机构进行了动力学建模，并将其作为验证模型与不同的接触模型进行比较，最后与实验结果进行对比。

在基于动力学模型的控制器设计中，关节间隙引起的关节效应等不确定因素会影响机构的控制效果，Li等^[17]对一种2自由度并联机构进行动力学建模，并根据机构动力学模型设计了一种滑模变结构控制器，实验结果验证了其控制器设计的正确性；文献[18]基于平面2自由度并联机器人的位置控制，结合自适应控制设计了一种动态模糊神经网络控制器。上述控制方法需要满足对动力学模型进行参数化处理的要求或需建立在机构关节效应及模型参数不确定项上界已知的情况下。为此，许多学者将自适应控制、神经网络、模糊控制，遗传算法以及滑模控制等相互组合成混合控制器以提升控制效果^[19-20]。

本文基于3-CPaRR机构，对该系统开展含间隙的动力学建模，比较了不同间隙量对其产生的影响，最后设计了该并联机构自适应滑模运动控制器，以提高控制系统的鲁棒性。

1 3-CPaRR机构含间隙关节运动学及接触力模型

3-CPaRR样机模型如图 1所示，各支链间的位置约束方程及各构件参数与其铰接转动副处夹角之间的关系请参见文献[21]。

1.1 间隙关节运动学建模

当考虑驱动构件与定机架铰接处的圆柱副中存在间隙时，驱动关节A_i及其运动模型如图 2所示，由于A_i处含间隙，当系统运动时，其在O_r-X_rY_rZ_r坐标系中沿X_r、Y_r、Z_r方向做随机运动，故引入概率密度函数模型描述A_i关节处的矢量变化^[22]。

$ {e_r} = {e_{r\max }}\int {{f_r}} \left( {{r_r}\left( {{x_r},{y_r}} \right)} \right){\rm{d}}{r_r} $

机构运动时，A_i关节中心O_r在O_r-X_rY_rZ_r中沿直径、轴方向上随机位置矢量值f_r(r_r)、f_a(z_r)可分别表示为：

$ {f_r}\left( {{r_r}} \right) = \left\{ {\begin{array}{*{20}{l}} {k{{\rm{e}}^{ - {{\left( {\frac{{{r_r} - \mu }}{{\sqrt {2\sigma } }}} \right)}^2}}},}&{0 \le r_r^2 = x_r^2 + y_r^2 \le e_{r\max }^2}\\ {0,}&{其他} \end{array}} \right. $

(1)

$ {f_a}\left( {{z_r}} \right) = \left\{ {\begin{array}{*{20}{l}} {k{{\rm{e}}^{ - {{\left( {\frac{{{z_r} - \mu }}{{\sqrt {2\sigma } }}} \right)}^2}}},}&{0 \le {z_r} \le {e_{a\max }}}\\ {0,}&{其他} \end{array}} \right. $

(2)

间隙引起的偏心距、轴向位移变化量分别为：

$ {e_r} = {e_{r\max }}\int {{f_r}\left( {{r_r}\left( {{x_r},{y_r}} \right)} \right){\rm{d}}r} $

(3)

$ {e_a} = {e_{a\max }}\int {{f_a}} \left( {{z_r}} \right){\rm{d}}{z_r} $

(4)

驱动关节A_i运动模型如图 3所示，Q为构件A_iB_i上的点，当A_iB_i伴随机构运动时，其在O_r-X_rY_rZ_r坐标系中的坐标x_Q、y_Q、z_Q为：

$ \left\{ {\begin{array}{*{20}{l}} {{x_Q} = {e_r}\cos \zeta + {l_{k + 1}}\cos \varepsilon }\\ {{y_Q} = {e_r}\sin \zeta + {l_{k + 1}}\sin \varepsilon }\\ {{z_Q} = {z_r} + {e_a}} \end{array}} \right. $

(5)

1.2 间隙关节接触力模型

图 4为关节A_i处的接触模型，基于图 4可得：

$ \left\{ {\begin{array}{*{20}{l}} {{\mathit{\boldsymbol{r}}_{r1}} = {\mathit{\boldsymbol{r}}_{rr}} + {\mathit{\boldsymbol{n}}_r}R}\\ {{\mathit{\boldsymbol{r}}_{r2}} = {\mathit{\boldsymbol{r}}_{ra}} + {\mathit{\boldsymbol{n}}_r}r} \end{array}} \right. $

(6)

关节A_i处的速度矢量可表示为：

$ \left\{ {\begin{array}{*{20}{l}} {{\mathit{\boldsymbol{v}}_{nr}} = \left[ {{{\left( {{{\mathit{\boldsymbol{\dot r}}}_{r1}} - {{\mathit{\boldsymbol{\dot r}}}_{r2}}} \right)}^{\rm{T}}}{\mathit{\boldsymbol{n}}_r}} \right]{\mathit{\boldsymbol{n}}_r}}\\ {{\mathit{\boldsymbol{v}}_{tr}} = \left( {{{\mathit{\boldsymbol{\dot r}}}_{r1}} - {{\mathit{\boldsymbol{\dot r}}}_{r2}}} \right) - {\mathit{\boldsymbol{v}}_{nr}}}\\ {{\mathit{\boldsymbol{v}}_{na}} = \left[ {{{\left( {{{\mathit{\boldsymbol{\dot r}}}_{a1}} - {{\mathit{\boldsymbol{\dot r}}}_{a2}}} \right)}^{\rm{T}}}{\mathit{\boldsymbol{n}}_a}} \right]{\mathit{\boldsymbol{n}}_a}} \end{array}} \right. $

(7)

在3-CPaRR并联机构中，其在刚性且A_i关节处存在间隙条件下时，关节处的接触力与摩擦力模型分别为：

$ {\mathit{\boldsymbol{F}}_n} = K{\delta ^{1.5}}\left[ {1 + \frac{{8\left( {1 - {c_r}} \right)}}{{5{c_r}}}\frac{{\dot \delta }}{{{{\mathit{\boldsymbol{\dot \delta }}}^{( - )}}}}} \right] $

(8)

$ {\mathit{\boldsymbol{F}}_t} = - {\mu _d}{c_d}{F_n}{\mathit{\boldsymbol{v}}_t} $

(9)

式(1)~(9)中相关变量意义及A_i处是否产生碰撞的条件判别等见文献[22]。

2 含关节间隙3-CPaRR并联机构的动力学模型 2.1 含间隙动力学模型

因间隙而使关节A_i处产生的接触力F_if需向系统构件质心进行等效转化，设其为f_i，则：

$ {\mathit{\boldsymbol{F}}_{if}} = {F_{in}}\mathit{\boldsymbol{n}} + {F_{it}}\mathit{\boldsymbol{t}} + {F_{ia}}\mathit{\boldsymbol{n}},i = 1,2,3 $

(10)

$ {\mathit{\boldsymbol{f}}_i} = \sum\limits_{i = 1}^3 {{\mathit{\boldsymbol{R}}_i}} {\mathit{\boldsymbol{F}}_{if}} = {\left[ {{f_{ix}}\quad {f_{iy}}\quad {f_{ix}}} \right]^{\rm{T}}},i = 1,2,3 $

(11)

构件l_1i、l_2i、l_3i和l_4i (i=1, 2, 3)在坐标系O-XYZ中质心坐标为(x_jic, y_jic, z_jic)，速度为v_ji=[ẋ_jic  ẏ_jic  ż_jic]^T，转动惯量为J_ji，动平台相对速度为${{\mathit{\boldsymbol{v}}}_{0}}={{\left[ \begin{array}{*{35}{l}} {{{\dot{q}}}_{1}} & {{{\dot{q}}}_{2}} & {{{\dot{q}}}_{3}} \\ \end{array} \right]}^{\text{T}}}$，则系统动能T、势能V分别表示为：

$ \begin{array}{l} T = \frac{1}{2}\sum\limits_{j = 1}^4 {\sum\limits_{i = 1}^3 {\left[ {{m_{ji}}\left( {\dot x_{jic}^2 + \dot y_{jic}^2 + \dot z_{jic}^2} \right)} \right]} } + \frac{1}{2}\sum\limits_{i = 1}^3 {{J_{1i}}} \dot \theta _i^2 + \\ \;\;\;\;\;\frac{1}{2}\sum\limits_{i = 1}^3 {{J_{2i}}} \dot \theta _{(i + 3)}^2 + \frac{1}{2}\sum\limits_{i = 1}^3 {{J_{3i}}} \dot \theta _i^2 + \frac{1}{2}{m_0}\mathit{\boldsymbol{v}}_0^{\rm{T}}{\mathit{\boldsymbol{v}}_0} \end{array} $

(12)

$ V = \sum\limits_{j = 1}^4 {\sum\limits_{i = 1}^3 {{m_{ji}}} } g{z_{jic}} + {m_0}g\left( {{Z_p} + {q_3}} \right) $

(13)

将式(12)、(13)代入非完整系统的Lagrange方程：

$ \frac{d}{{dt}}\left( {\frac{{\partial T}}{{\partial {{\dot q}_i}}}} \right) - \frac{{\partial T}}{{\partial {q_i}}} + \frac{{\partial V}}{{\partial {q_i}}} = {F_i} + {Q_i},i = 1,2,3 $

(14)

化简后可得：

$ \mathit{\boldsymbol{M}}\left( \mathit{\boldsymbol{q}} \right)\mathit{\boldsymbol{\ddot q}} + \mathit{\boldsymbol{C}}(\mathit{\boldsymbol{q}},\mathit{\boldsymbol{\dot q}})\mathit{\boldsymbol{\dot q}} + \mathit{\boldsymbol{G}}(\mathit{\boldsymbol{q}}) = \mathit{\boldsymbol{F}} + \mathit{\boldsymbol{Q}} $

(15)

式中：

$ \mathit{\boldsymbol{M}}\left( \mathit{\boldsymbol{q}} \right) = \left[ {\begin{array}{*{20}{c}} {{{\mathit{\boldsymbol{\hat M}}}_{11}}}&{{{\mathit{\boldsymbol{\hat M}}}_{12}}}&{{{\mathit{\boldsymbol{\hat M}}}_{13}}}\\ {{{\mathit{\boldsymbol{\hat M}}}_{12}}}&{{{\mathit{\boldsymbol{\hat M}}}_{22}}}&{{{\mathit{\boldsymbol{\hat M}}}_{23}}}\\ {{{\mathit{\boldsymbol{\hat M}}}_{13}}}&{{{\mathit{\boldsymbol{\hat M}}}_{23}}}&{{{\mathit{\boldsymbol{\hat M}}}_{33}}} \end{array}} \right], $

$ \mathit{\boldsymbol{C}}(\mathit{\boldsymbol{q}},\mathit{\boldsymbol{\dot q}}) = \left[ {\begin{array}{*{20}{l}} {{\mathit{\boldsymbol{C}}_{11}}}&{{\mathit{\boldsymbol{C}}_{12}}}&{{\mathit{\boldsymbol{C}}_{13}}}\\ {{\mathit{\boldsymbol{C}}_{21}}}&{{\mathit{\boldsymbol{C}}_{22}}}&{{\mathit{\boldsymbol{C}}_{23}}}\\ {{\mathit{\boldsymbol{C}}_{31}}}&{{\mathit{\boldsymbol{C}}_{32}}}&{{\mathit{\boldsymbol{C}}_{33}}} \end{array}} \right], $

$ \mathit{\boldsymbol{G}}(\mathit{\boldsymbol{q}}) = \left[ {\begin{array}{*{20}{l}} {\frac{{\partial V}}{{\partial {q_1}}}}&{\frac{{\partial V}}{{\partial {q_2}}}}&{\frac{{\partial V}}{{\partial {q_3}}}} \end{array}} \right], $

$ \mathit{\boldsymbol{Q}} = {\left[ {\begin{array}{*{20}{c}} {\sum\limits_{i = 1}^3 {{f_{ix}}} }&{\sum\limits_{i = 1}^3 {{f_{iy}}} }&{\sum\limits_{i = 1}^3 {{f_{iz}}} } \end{array}} \right]^{\rm{T}}}, $

$ {\mathit{\boldsymbol{C}}_{11}} = \frac{1}{2}\frac{{\partial {{\hat M}_{11}}}}{{\partial {q_1}}}\dot q_1^2 + \frac{{\partial {{\hat M}_{11}}}}{{\partial {q_2}}}\dot q_2^2, $

$ {\mathit{\boldsymbol{C}}_{12}} = \left( {\frac{{\partial {{\hat M}_{12}}}}{{\partial {q_2}}} - \frac{1}{2}\frac{{\partial {{\hat M}_{22}}}}{{\partial {q_1}}}} \right)\dot q_2^2 + \left( {\frac{{\partial {{\hat M}_{12}}}}{{\partial {q_3}}} + \frac{{\partial {{\hat M}_{13}}}}{{\partial {q_2}}} - \frac{{\partial {{\hat M}_{23}}}}{{\partial {q_1}}}} \right)\dot q_3^2, $

$ {\mathit{\boldsymbol{C}}_{13}} = \left( {\frac{{\partial {{\hat M}_{13}}}}{{\partial {q_3}}} - \frac{1}{2}\frac{{\partial {{\hat M}_{33}}}}{{\partial {q_1}}}} \right)\dot q_3^2 + \frac{{\partial {{\hat M}_{11}}}}{{\partial {q_3}}}{{\dot q}_1}{{\dot q}_3}, $

$ {\mathit{\boldsymbol{C}}_{21}} = \left( {\frac{{\partial {{\hat M}_{12}}}}{{\partial {q_1}}} - \frac{1}{2}\frac{{\partial {{\hat M}_{11}}}}{{\partial {q_2}}}} \right)\dot q_1^2 + \frac{{\partial {{\hat M}_{22}}}}{{\partial {q_1}}}\dot q_2^2, $

$ {\mathit{\boldsymbol{C}}_{22}} = \frac{1}{2}\frac{{\partial {{\hat M}_{22}}}}{{\partial {q_2}}}\dot q_2^2 + \frac{{\partial {{\hat M}_{22}}}}{{\partial {q_3}}}\dot q_3^2, $

$ {\mathit{\boldsymbol{C}}_{23}} = \left( {\frac{{\partial {{\hat M}_{23}}}}{{\partial {q_3}}} - \frac{1}{2}\frac{{\partial {{\hat M}_{33}}}}{{\partial {q_2}}}} \right)\dot q_3^2 + \left( {\frac{{\partial {{\hat M}_{12}}}}{{\partial {q_3}}} - \frac{{\partial {{\hat M}_{13}}}}{{\partial {q_2}}} + \frac{{\partial {{\hat M}_{23}}}}{{\partial {q_1}}}} \right)\dot q_1^2, $

$ {\mathit{\boldsymbol{C}}_{31}} = \left( {\frac{{\partial {{\hat M}_{{\rm{13}}}}}}{{\partial {q_1}}} - \frac{1}{2}\frac{{\partial {{\hat M}_{{\rm{11}}}}}}{{\partial {q_3}}}} \right)\dot q_{\rm{l}}^2 + \left( {\frac{{\partial {{\hat M}_{{\rm{13}}}}}}{{\partial {q_2}}} - \frac{{\partial {{\hat M}_{{\rm{13}}}}}}{{\partial {q_3}}} + \frac{{\partial {{\hat M}_{{\rm{23}}}}}}{{\partial {q_1}}}} \right)\dot q_2^2, $

$ {\mathit{\boldsymbol{C}}_{32}} = \left( {\frac{{\partial {{\hat M}_{23}}}}{{\partial {q_2}}} - \frac{1}{2}\frac{{\partial {{\hat M}_{22}}}}{{\partial {q_3}}}} \right)\dot q_2^2 + \frac{{\partial {{\hat M}_{33}}}}{{\partial {q_2}}}\dot q_3^2, $

$ {\mathit{\boldsymbol{C}}_{33}} = \frac{1}{2}\frac{{\partial {{\hat M}_{33}}}}{{\partial {q_3}}}\dot q_3^2 + \frac{{\partial {{\hat M}_{33}}}}{{\partial {q_1}}}\dot q_1^2 $

2.2 数值分析

基于上述模型，3-CPaRR机构具体参数如下：材料为碳素钢，密度ρ=7 800 kg/m³，弹性模量分别为E=210 GPa，G=80 GPa，l_1i=160 mm，l_2i=150 mm，l_3i=160 mm，l_4i=55 mm，截面尺寸为b=15 mm，h=10 mm，m₀=0.875 kg，周期T=5 s，步长Δt=0.01 s。J_ix、J_iy、J_iz为构件绕X、Y、Z转动惯量，如表 1所示。

给定轨迹如下：

$ \left\{ {\begin{array}{*{20}{l}} {X = {x_0} - 5\sin (0.4{\rm{ \mathsf{ π} }}t)}\\ {Y = {x_0} - 6\sin (0.4{\rm{ \mathsf{ π} }}t)\quad (0 \le t \le 5)}\\ {Z = {z_0} - 5t} \end{array}} \right. $

(16)

图 5、6为模型的求解结果，由图可得，当构件铰接A_i处为理想运动副时，该处的驱动力F_i为光滑曲线，且不产生波动。当考虑间隙的存在时，F_i开始产生不同幅值的震荡，这说明当机构运行时，间隙的存在使A_i处有碰撞产生，且随着间隙量的增大而愈加剧烈，降低系统稳定性，缩短构件的使用寿命。

2.3 基于ADAMS软件的虚拟样机仿真结果对比

为验证基于含关节间隙转动副模型所建立的支链运动副模型以及考虑关节间隙的并联机构动力学建模的合理性，本文基于ADAMS软件对含关节间隙的3-CPaRR并联机构虚拟样机仿真结果进行原理上的对比分析，以验证该方法的有效性与可参考性。

基于上述对该机构进行仿真验证时，需对间隙关节进行接触定义，在ADAMS软件中有如下2类接触模型：二维接触模型与三维接触模型。其中接触中的法向力模型有回归接触模型和碰撞函数接触模型。其中碰撞函数接触是指ADAMS/Solver运用函数库中Impact函数来计算接触力。对于3-CPaRR并联机构，根据其结构及加工制造时采用的不同材料，定义法向接触力模型参数与库伦摩擦力模型参数，如表 2所示。

在进行了含关节间隙的该3自由度并联机构虚拟样机仿真后(定义关节间隙量c_r=c_a=0.5 mm)，取运动副A₃处的Z方向的广义驱动力F₃作为对照，其广义驱动力变化曲线如图 7所示。

由图 7可知，含关节间隙的3-CPaRR并联机构在运动过程中，含关节间隙的转动副关节元素间的相对运动存在随机性，基于此结论，验证了本文选取概率统计数学模型理论为基点来描述转动副关节元素在其径向与轴向运动规律数学模型的合理性，同时，图 7也可验证前述基于正态分布概率统计模型来计算分析含关节间隙3-CPaRR并联机构的运动学与动力学一定的正确性与合理性。

此外，可见提出的关关节间隙数学模型的分析结果与ADAMS仿真结果存在一定的偏差，经分析其原因有如下2点：1)含关节间隙运动副模型中的法向接触力模型与切向接触力模型与ADAMS软件中的接触力模型不完全一致，这是不可避免的原因；2)建模计算过程中，参考文献[23-24]，定义了含关节间隙并联机构的各转动关节的广义坐标与封闭方程求得角θ_i近似相等；在本文含关节间隙并联机构的动力学方程建模过程中，因考虑机构动力学模型的计算复杂问题，为提高计算效率，进行了部分简化计算。

3 3-CPaRR并联机构控制系统设计

本文以滑模控制理论为基础设计其运动控制器。如机构系统有模型参数误差以及间隙等关节效应，这些不确定项的存在均会导致未知扰动W，因此，式(15)可表示为：

$ \left( {{\mathit{\boldsymbol{M}}_0} + \Delta \mathit{\boldsymbol{M}}} \right)\mathit{\boldsymbol{\ddot q}} + \left( {{\mathit{\boldsymbol{C}}_0} + \Delta \mathit{\boldsymbol{C}}} \right)\mathit{\boldsymbol{\dot q}} + \left( {{\mathit{\boldsymbol{G}}_0} + \Delta \mathit{\boldsymbol{G}}} \right) = \mathit{\boldsymbol{F}} + \mathit{\boldsymbol{W}} $

(17)

式中：M₀、C₀、G₀、W₀分别为M、C、G、W的名义值；ΔM=M-M₀，ΔC=C-C₀，ΔG=G-G₀，ΔW=W-W₀表示为建模误差。该动力学模型式(17)具有如下一些性质：1) M∈R^3×3为正定矩阵，且有界；2)Ṁ-2C为反对称矩阵，即存在x^T(Ṁ-2C)x=0其中x∈Rⁿ。

由理想多体机构系统动力学方程中的参数不确定项和负载扰动集中起来作为总的不确定项，记为d∈R^3×3，表达式为：

$ \mathit{\boldsymbol{d}} = \Delta \mathit{\boldsymbol{M\ddot q}} + \Delta \mathit{\boldsymbol{C\dot q}} + \Delta \mathit{\boldsymbol{G}} + \Delta \mathit{\boldsymbol{f}} + \mathit{\boldsymbol{W}} $

(18)

式中：∣d_i∣≤b_i，b_i>0。

常规滑模控制的滑模面定义为：定义q_d、${{{\mathit{\boldsymbol{\dot{q}}}}}_{d}}$∈R^3×3分别为动平台位姿、速度的期望轨迹，跟踪误差表达式为$\mathit{\boldsymbol{\tilde q}} = \mathit{\boldsymbol{q}} - {\mathit{\boldsymbol{q}}_d}$，$\mathit{\boldsymbol{\dot {\tilde q} = \dot q}} - {{\mathit{\boldsymbol{\dot q}}}_d}$，Λ为正定矩阵。当机构模型不确定上界b已知时，可采用常规滑模控制器控制：

$ \begin{array}{*{20}{c}} {\mathit{\boldsymbol{F}} = {\mathit{\boldsymbol{u}}_{eq}} - \mathit{\boldsymbol{ \boldsymbol{\varGamma} }}{\rm{sgn}}(\mathit{\boldsymbol{s}})}\\ {{\mathit{\boldsymbol{u}}_{eq}} = {\mathit{\boldsymbol{M}}_0}\left( {{{\mathit{\boldsymbol{\ddot q}}}_d} + \mathit{\boldsymbol{ \boldsymbol{\varLambda} \dot {\tilde q}}}} \right) + {\mathit{\boldsymbol{C}}_0}\left( {{{\mathit{\boldsymbol{\dot q}}}_d} + \mathit{\boldsymbol{ \boldsymbol{\varLambda} \tilde q}}} \right) + {\mathit{\boldsymbol{G}}_0}\mathit{\boldsymbol{q}} - {\mathit{\boldsymbol{W}}_0}} \end{array} $

(19)

式中，控制律可分为2个部分：

1) u_eq是基于机构系统建立的名义动力学模型，作为控制器前馈控制项，用于在线补偿各个构件的惯性力及运动关节之间的冲击力与摩擦力。

2) -ρsgn(s)是针对d的补偿项，sgn(s)是符号函数，ρ=diag(ρ_i)是正定切换增益矩阵，当ρ_i>b_i时，系统可趋于稳定。

3.1 自适应滑模控制

常规滑模控制律(19)存有以下问题：该并联机构运动过程中关节间隙的影响是复杂多变的，因不确定项的上界b难以准确预估。针对此问题，设计该并联机构的自适应滑模控制器，利用自适应算法进行实时补偿，可省去对不确定项b的不确定估计，提高控制器的准确性与鲁棒性。考虑并联机构的模型参数及外界扰动等不确定项，设定α、β和ε为并联机构物理参数常数。定义并联机构的系统参数为a，a=[α β ε]^T，â为a的估计，$\mathit{\boldsymbol{\hat a}} = {\left[ {\begin{array}{*{20}{c}} {\hat \alpha }&{\hat \beta }&{\hat \varepsilon } \end{array}} \right]^{\rm{T}}}$。取ã=â-a，由于a为常数向量，则$\mathit{\boldsymbol{\dot {\tilde a}}} = \mathit{\boldsymbol{\dot {\hat a}}}$。

由上述可知α、β和ε为未知常数，运动轨迹误差${\mathit{\boldsymbol{\tilde q}}}$(t)=q(t)-q_d(t)，定义该并联机构的滑模控制器参数：${{\mathit{\boldsymbol{\dot q}}}_r} = {{\mathit{\boldsymbol{\dot q}}}_d} - \mathit{\boldsymbol{ \boldsymbol{\varLambda} \tilde q}}, {{\mathit{\boldsymbol{\ddot q}}}_r} = {{\mathit{\boldsymbol{\ddot q}}}_d} - \mathit{\boldsymbol{ \boldsymbol{\varLambda} }}\dot {\tilde q}$，式中Λ=diag(λ₁  λ₂  λ₃)，λ_i>0，(i=1, 2, 3)。

滑模函数为：

$ \mathit{\boldsymbol{s}} = \mathit{\boldsymbol{\dot {\tilde q}}} - \mathit{\boldsymbol{ \boldsymbol{\varLambda} \tilde q}} $

(20)

设计控制器为：

$ \mathit{\boldsymbol{F}} = \mathit{\boldsymbol{\hat M}}\left( \mathit{\boldsymbol{q}} \right){{\mathit{\boldsymbol{\ddot q}}}_r} + \mathit{\boldsymbol{\hat C}}(\mathit{\boldsymbol{q}},\mathit{\boldsymbol{\dot q}}){{\mathit{\boldsymbol{\dot q}}}_r} + \mathit{\boldsymbol{\hat G}}(\mathit{\boldsymbol{q}}) - {\mathit{\boldsymbol{K}}_d}\mathit{\boldsymbol{s}} $

(21)

其中：K_d=diag(k_d1  k_d2  k_d3)，k_d>0，i=1, 2, 3。

由于M为正定矩阵，设计稳定性判定的Lyapunov函数为：

$ V = \frac{1}{2}{\mathit{\boldsymbol{s}}^{\rm{T}}}\mathit{\boldsymbol{Ms}} + \frac{1}{2}{\mathit{\boldsymbol{a}}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varGamma} a}} $

(22)

其中：Γ=diag(γ₁  γ₂  γ₃)，γ_i>0，i=1, 2, 3。

因此，式(22)表示为：

$ \begin{array}{l} \mathit{\boldsymbol{\dot V}} = \frac{1}{2}{\mathit{\boldsymbol{s}}^{\rm{T}}}\mathit{\boldsymbol{M\dot s}} + \frac{1}{2}{\mathit{\boldsymbol{s}}^{\rm{T}}}\mathit{\boldsymbol{\dot Ms}} + {{\mathit{\boldsymbol{\tilde a}}}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varGamma} \dot {\tilde a}}} = \\ \;\;\;\;\;\;{\mathit{\boldsymbol{s}}^{\rm{T}}}\left( {\mathit{\boldsymbol{M\ddot q}} - \mathit{\boldsymbol{M}}{{\mathit{\boldsymbol{\ddot q}}}_r}} \right) + \frac{1}{2}{\mathit{\boldsymbol{s}}^{\rm{T}}}\mathit{\boldsymbol{\dot Ms}} + {{\mathit{\boldsymbol{\tilde a}}}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varGamma} \dot {\tilde a}}} = \\ \;\;\;\;\;\;{\mathit{\boldsymbol{s}}^{\rm{T}}}\left( {\mathit{\boldsymbol{F}} - \mathit{\boldsymbol{C\dot q}} - \mathit{\boldsymbol{G}} - \mathit{\boldsymbol{M}}{{\mathit{\boldsymbol{\ddot q}}}_r}} \right) + \frac{1}{2}{\mathit{\boldsymbol{s}}^{\rm{T}}}\mathit{\boldsymbol{\dot Ms}} + {{\mathit{\boldsymbol{\tilde a}}}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varGamma} \dot {\tilde a}}} = \\ \;\;\;\;\;\;{\mathit{\boldsymbol{s}}^{\rm{T}}}\left( {\mathit{\boldsymbol{F}} - \mathit{\boldsymbol{C}}\left( {\mathit{\boldsymbol{s}} + {{\mathit{\boldsymbol{\dot q}}}_r}} \right) - \mathit{\boldsymbol{G}} - \mathit{\boldsymbol{M}}{{\mathit{\boldsymbol{\ddot q}}}_r}} \right) + \frac{1}{2}{\mathit{\boldsymbol{s}}^{\rm{T}}}\mathit{\boldsymbol{\dot Ms}} + {{\mathit{\boldsymbol{\tilde a}}}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varGamma} \dot {\tilde a}}} \end{array} $

(23)

将控制律式(19)代入式(22)，可得：

$ \begin{array}{l} \dot V = {\mathit{\boldsymbol{s}}^{\rm{T}}}\left( {\mathit{\boldsymbol{\hat M}}{{\mathit{\boldsymbol{\ddot q}}}_r} + \mathit{\boldsymbol{\hat C}}{{\mathit{\boldsymbol{\dot q}}}_r} + \mathit{\boldsymbol{\hat G}} - {\mathit{\boldsymbol{K}}_d}\mathit{\boldsymbol{s}} - \mathit{\boldsymbol{C}}\left( {\mathit{\boldsymbol{s}} + {{\mathit{\boldsymbol{\dot q}}}_r}} \right) - \mathit{\boldsymbol{G}} - } \right.\\ \;\;\;\;\;\;\left. {\mathit{\boldsymbol{M}}{{\mathit{\boldsymbol{\ddot q}}}_r}} \right) + \frac{1}{2}{\mathit{\boldsymbol{s}}^{\rm{T}}}\mathit{\boldsymbol{\dot Ms}} + {{\mathit{\boldsymbol{\tilde a}}}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varGamma} \dot {\tilde a}}} = \\ \;\;\;\;\;\;{\mathit{\boldsymbol{s}}^{\rm{T}}}\left( {\mathit{\boldsymbol{\tilde M}}{{\mathit{\boldsymbol{\ddot q}}}_r} + \mathit{\boldsymbol{\tilde C}}{{\mathit{\boldsymbol{\dot q}}}_r} + \mathit{\boldsymbol{G}} - {\mathit{\boldsymbol{K}}_d}\mathit{\boldsymbol{s}} - \mathit{\boldsymbol{Cs}}} \right) + \frac{1}{2}{\mathit{\boldsymbol{s}}^{\rm{T}}}\mathit{\boldsymbol{\dot Ms}} + \\ \;\;\;\;\;\;{{\mathit{\boldsymbol{\tilde a}}}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varGamma} \dot {\tilde a}}} \end{array} $

(24)

根据该并联机构动力学方程的线性化特性，有：

$ \mathit{\boldsymbol{\tilde M}}{{\mathit{\boldsymbol{\ddot q}}}_r} + \mathit{\boldsymbol{\tilde C}}{{\mathit{\boldsymbol{\dot q}}}_r} + \mathit{\boldsymbol{\tilde G}} = \mathit{\boldsymbol{Y}}\left( {\mathit{\boldsymbol{q}},\mathit{\boldsymbol{\dot q}},{{\mathit{\boldsymbol{\dot q}}}_r},{{\mathit{\boldsymbol{\ddot q}}}_r}} \right)\mathit{\boldsymbol{\tilde a}} $

(25)

$ \begin{array}{l} \dot V = {\mathit{\boldsymbol{s}}^{\rm{T}}}\left( {\mathit{\boldsymbol{Y\tilde a}} - {\mathit{\boldsymbol{K}}_d}\mathit{\boldsymbol{s}} - \mathit{\boldsymbol{Cs}}} \right) + \frac{1}{2}{\mathit{\boldsymbol{s}}^{\rm{T}}}\mathit{\boldsymbol{\dot Ms}} + {{\mathit{\boldsymbol{\tilde a}}}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varGamma} \dot {\tilde a}}} = \\ \;\;\;\;\;\;{\mathit{\boldsymbol{s}}^{\rm{T}}}\left( {\mathit{\boldsymbol{Y\tilde a}} - {\mathit{\boldsymbol{K}}_d}\mathit{\boldsymbol{s}}} \right) + \frac{1}{2}{\mathit{\boldsymbol{s}}^{\rm{T}}}\left( {\mathit{\boldsymbol{\dot M}} - 2\mathit{\boldsymbol{C}}} \right)\mathit{\boldsymbol{s}} + {{\mathit{\boldsymbol{\tilde a}}}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varGamma} }}\dot {\tilde a} = \\ \;\;\;\;\;\;{\mathit{\boldsymbol{s}}^{\rm{T}}}\left( {\mathit{\boldsymbol{Y\tilde a}} - {\mathit{\boldsymbol{K}}_d}\mathit{\boldsymbol{s}}} \right) + {{\mathit{\boldsymbol{\tilde a}}}^{\rm{T}}}\mathit{\boldsymbol{ \boldsymbol{\varGamma} \dot {\tilde a}}} = \\ \;\;\;\;\;\;{{\mathit{\boldsymbol{\tilde a}}}^{\rm{T}}}\left( {{\mathit{\boldsymbol{Y}}^{\rm{T}}}\mathit{\boldsymbol{s}} + \mathit{\boldsymbol{ \boldsymbol{\varGamma} \dot {\tilde a}}}} \right) - {\mathit{\boldsymbol{s}}^{\rm{T}}}{\mathit{\boldsymbol{K}}_d}\mathit{\boldsymbol{s}} \end{array} $

(26)

设计自适应律为：

$ \mathit{\boldsymbol{\dot {\hat a}}} = - {\mathit{\boldsymbol{ \boldsymbol{\varGamma} }}^{ - 1}}{\mathit{\boldsymbol{Y}}^{\rm{T}}}\mathit{\boldsymbol{s}} $

(27)

则${\dot V}$=-s^T K_ds≤0。

取${\dot V}$≡0则s≡0，根据LaSalle不变集定理，当t→∞时，s→0，即$\mathit{\boldsymbol{\dot q}}\left( t \right) \to 0$，$\mathit{\boldsymbol{\dot {\tilde q}}}\left( t \right) \to 0$。由于当且仅当s=0时，${\dot V}$=0，此时V不变，因此无法保证ã→0。

为实现3-CPaRR并联机构的高精度运动控制效果，本节将含关节间隙的转动副引起的并联机构冲击负载定义为外界扰动，基于上文的自适应滑模控制器，对含间隙的该并联机构进行数值仿真，该3自由度并联机构自适应滑模变结构控制的系统框图如图 8所示。

3.2 仿真分析

基于上述，在该机构的自适应滑模控制仿真过程中，其参数取：Λ=diag(5  5  5)，K_d=diag(10  10  10)，Γ=diag(1  1  1)。为验证误差的收敛性，假设该机构存在初始误差。机构输出端在初始时刻的状态为：q₁(0)=q₂(0)=q₃(0)=1.2 mm，${{\dot q}_1}$(0)=${{\dot q}_2}$(0)=${{\dot q}_3}$(0)=0。仿真时间为t=20 s，Δt=0.000 1 s。结果如图 9所示。

为评定控制器的控制效果，引入2项评定指标：平均控制误差AVE和最大控制误差MAX，其中将AVE表示为:

$ {\rm{AVE}} = \frac{1}{{{t_m}}}\int_0^{{t_m}} {\left| {e\left( t \right)} \right|{\rm{d}}t} $

(28)

式中t_m表示运行时间。

图 9对比了常规滑模控制和自适应滑模控制下3自由度并联机构动平台正交方向上运动控制误差曲线，尽管机构系统存在初始误差以及时变负载扰动影响，但系统依然可以保证在一定时间内趋近平稳，从而验证了上述的2种控制方法都可以使该并联机构的运动精度误差收敛至稳定值。自适应滑模控制下，动平台轨迹误差能在0.9 s内收敛到稳定状态，收敛的速度明显比前者控制更快。

表 3和表 4分别是2种控制方法的AVE指标和MAX指标的对比。在运动精度达到稳定状态后，采用常规滑模控制q₁、q₂、q₃的平均跟踪误差(10^-5 mm)分别为2.72、3.87、5.51，最大控制误差值(10^-5 mm)分别为6.00、4.78、5.83，与之相比，采用自适应滑模控制时的误差值(10^-5 mm)分别为0.53、1.34、0.67，最大控制误差值(10^-5 mm)分别为1.35、2.82、1.72，说明了自适应滑模控制不仅收敛速度快，而且控制精度明显比常规滑模控制更高，对模型不确定性及外界扰动的抑制能力更强。

图 9是2种控制方法的机构输出运动轨迹q₁、q₂、q₃的跟踪误差，从中可得该机构输出的运动轨迹虽然有抖振现象但并没有出现奇异性，因此证明自适应滑模控制器较常规滑模控制器的优越性。

如图 9所示，在t=1 s时刻左右，该机构实际输出曲线开始稳定跟踪期望轨迹，这说明该控制器的设计是合理正确的。

4 结论

1) 基于正态分布、Flores接触力、Coulomb摩擦模型及非完整Lagrange方程实现了该机构间隙动力学建模与求解，并基于ADAMS软件进行了含间隙模型仿真，验证理论建模可行性。

2) 针对常规滑模变结构控制的抖振与无法准确抑制时变扰动问题，设计了该并联机构自适应滑模运动控制器，控制器省去了对不确定项上界的估计，对系统状态实时补偿，对比仿真计算结果显示提高了系统的鲁棒性。

本文从理论上对3-CPaRR并联机构进行了含间隙动力学建模及自适用滑模控制器的设计，在后续的研究工作中将对该机构的样机展开实验研究，同时考虑润滑、摩擦等因素对该系统相关特性产生的影响。