随着国内人口老龄化的加速，中风发病率呈逐年增加的趋势. 据统计可知，国内中风患者目前已达1 300多万，预计2030年可能会超过3 000万人^[1]. 中风患者的高效康复训练是当前的关注点之一，相对于传统理疗和药物治疗，康复机器人可以节约大量人力、物力，手腕康复训练装置是其中之一^[2-3].

近年来，国内外学者对腕部康复机器人/装置开展了较多的研究工作，如气动人工肌肉手腕康复装置^[4]、绳索驱动腕部康复机构^[5]、混联式腕关节康复机构^[6]以及基于3RRR球面并联机构(spherical parallel mechanism，SPM)的手腕康复装置^[7]. 其中，SPM手腕康复装置具有结构紧凑、贴合手腕运动特征的优点.

SPM具有俯仰、偏航和滚动的运动特点^[8]，Gosselin等^[9-11]对该机构的运动学特性开展了一定的研究工作. 通常将该机构逆运动学的研究转换为3个一元二次方程的求解问题，由于方程组最多存在8个解，给驱动系统控制器的设计带来了一定的挑战^[12]. 目前，3RRR球面并联机构逆运动学唯一解的求解方案有2个：方案1是建立3RRR球面并联机构逆运动学方程，选取一元二次方程正根对应的解作为机构逆运动学解^[13]；方案2是在8组解中选取与初始姿态相同工作模式的一组解作为唯一的逆运动学解^[14-15]. 上述2种方案都存在不足，方案1存在部分动平台姿态无逆运动学解的问题，方案2不能得到具体的逆运动学解析解.

针对共轴3RRR球面并联机构（CSPM）逆运动学唯一解求解方法的不足，提出基于欧拉角的逆运动学分步求解方法. 该方法将基于欧拉角描述的动平台姿态分解为绕Z轴和绕X、Y轴旋转的子姿态A和子姿态B，并分别求解出子姿态A的逆解和子姿态B的最优解，将2个子姿态最优解线性相加，可以获得CSPM唯一逆解. 将该方法与单位四元数相结合，并应用于球面手腕康复装置的轨迹规划中，以提高装置运动控制的实时性，这对CSPM控制器设计具有实际的工程意义.

人体腕部具有3个自由度，以人体腕部建立坐标系，腕部分别可实现绕Z轴转动的屈曲(F)/伸展(E)运动、绕Y轴转动的内收(A)/外展(O)和绕X轴转动的旋内(M)/旋外(L)运动，如图1所示. 由GB/T 15499—1995及其他资料^[16-17]可知，正常情况下实现所述6个运动时手腕转动角度 $\gamma $的范围，如表1所示.

共轴3RRR球面并联机构如图2所示，由3个近端连杆（弯曲弧度为 $ {\alpha }_{1} $）、3个远端连杆（弯曲弧度为 $ {\alpha }_{2} $）、动平台和底座组成，其中近端连杆为驱动杆. 近端连杆与底座、近端连杆与远端连杆、远端连杆与动平台均采用转动副连接，使得CSPM动平台在空间中只能作旋转运动，与手腕运动相贴合. 在初始条件下，CSPM动平台和底座呈水平状态，3根近端连杆呈120°均匀分布. CSPM的特点是所有关节轴线相交于O点，后文将其称为旋转中心.

以图2的O点为坐标原点建立固定坐标系，Z轴垂直于底座，方向向上. 在共轴3RRR球面并联机构中，3根驱动轴存在共轴关系，驱动杆轴向量、近端连杆与远端连杆相连的轴向量、动平台姿态向量分别采用U_i、W_i、V_i（i = 1，2，3）表示. V_i与Z轴的夹角为 $\; \beta $. 定义 $ Z $轴与 $ {\boldsymbol{W}}_{i} $所组成的平面为测量平面P_i（i = 1，2，3），Y轴位于P₁平面内，利用右手定理可以确定X轴方向.

设定驱动轴关节转角 $ \mathit{{\boldsymbol{\theta}} }=\left[{\theta }_{1},{\theta }_{2},{\theta }_{3}\right] $是测量平面 $ {P}_{i} $转动的角度，逆时针方向为正，顺时针为负，如图3所示. 在工作空间内，使用XYZ欧拉角 $ \left({\varphi }_{1},{\varphi }_{2},{\varphi }_{3}\right) $表示动平台姿态，初始姿态下，动平台呈水平状态，设定驱动杆转角为 $ \mathit{{\boldsymbol{\theta}} }=\left[{\theta }_{1},{\theta }_{2},{\theta }_{3}\right] $.

1) 初始姿态. 在固定坐标系下，3个驱动轴单位向量相同，

(1) $ {{\boldsymbol{U}}_{\text{1}}}{\text{ = }}{{\boldsymbol{U}}_2} = {{\boldsymbol{U}}_{\text{3}}} = {\left[ {0,\;0,\;- 1} \right]^{\text{T}}}. $

在初始时刻， $ {\boldsymbol{W}}_{i}^{*}(i=\mathrm{1,2},3) $为

(2) $ {\boldsymbol{W}}_i^{\mathbf{*}} = \left[ {\begin{array}{*{20}{c}} { - \sin\; {\eta _i}\sin\; {\alpha _1}} \\ {\cos \;{\eta _i}\sin\; {\alpha _1}} \\ { - \cos \;{\alpha _1}} \end{array}} \right]. $

式中： $ {\alpha }_{1} $为机构参数； $\eta _{i}={2\left(i-1\right)\pi }/{3},i=\mathrm{1,2},3$.

在初始姿态下， $ {\boldsymbol{V}}_{1} $水平投影与平面 $ {P}_{1}\mathrm{的} $夹角为 $ \delta $，如图3所示. 根据几何关系可知，在初始时刻，

(3) $ {\boldsymbol{V}}_i^{\mathbf{*}} = \left[ {\begin{array}{*{20}{c}} {\sin\; (\delta - {\eta _i})\sin \;\beta } \\ {\cos\; (\delta - {\eta _i})\sin \;\beta } \\ { - \cos \;\beta } \end{array}} \right]. $

式中： $ \;\beta $为机构参数.

2) 非初始姿态. 当输入角度为 $ {\boldsymbol{\theta }}_{i} $时， $ {\boldsymbol{W}}_{i} $为

(4) $ {{\boldsymbol{W}}_i} = \left[ {\begin{array}{*{20}{c}} {\left( { - \sin\; {\eta _i}\cos\; {\theta _i} - \sin\; {\theta _i}\cos\; {\eta _i}} \right)\sin\; {\alpha _1}} \\ {\left( { - \sin\; {\eta _i}\sin\; {\theta _i}+\cos\; {\theta _i}\cos\; {\eta _i}} \right)\sin\; {\alpha _1}} \\ { - \cos\; {\alpha _1}} \end{array}} \right]. $

动平台姿态 $ {\boldsymbol{V}}_{i} $为

(5) $ {{\boldsymbol{V}}_i} = {{\boldsymbol{R}}_Z}({\varphi _3}){{\boldsymbol{R}}_Y}({\varphi _2}){{\boldsymbol{R}}_X}({\varphi _1})\left[ {\begin{array}{*{20}{c}} {\sin\; (\delta - {\eta _i})\sin\; \beta } \\ {\cos\; (\delta - {\eta _i})\sin\; \beta } \\ { - \cos\; \beta } \end{array}} \right]. $

式中：

(6) $\left. \begin{split} & {{\boldsymbol{R}}_Z}({\varphi _3}) = \left[ {\begin{array}{*{20}{c}} {\cos\; {\varphi _3}}&{ - \sin\; {\varphi _3}}&0 \\ {\sin\; {\varphi _3}}&{\cos\; {\varphi _3}}&0 \\ 0&0&1 \end{array}} \right] , \\ & {{\boldsymbol{R}}_Y}({\varphi _2}) = \left[ {\begin{array}{*{20}{c}} {\cos\; {\varphi _2}}&0&{\sin\; {\varphi _2}} \\ 0&1&0 \\ { - \sin\; {\varphi _2}}&0&{\cos\; {\varphi _2}} \end{array}} \right] , \\ &{{\boldsymbol{R}}_X}({\varphi _1}) = \left[ {\begin{array}{*{20}{c}} 1&0&0 \\ 0&{\cos\; {\varphi _1}}&{ - \sin\; {\varphi _1}} \\ 0&{\sin\; {\varphi _1}}&{\cos\; {\varphi _1}} \end{array}} \right] , \end{split} \right\}$

其中， $ {\varphi}_{1}、{\varphi }_{2}、{\varphi }_{3} $为欧拉角参数.

使用欧拉角 $ \left({\mathrm{\varphi }}_{1},{\varphi }_{2},{\varphi }_{3}\right) $表示动平台姿态，CSPM最多存在8组逆运动学解，但只有一组解满足与初始姿态工作模式相同的要求，称为手腕康复装置的唯一逆运动学解.

联立式(4)、(5)，根据几何约束列出方程：

(7) $ {{\boldsymbol{W}}}_{i}\cdot{{\boldsymbol{V}}}_{i}=\mathrm{cos}\;{\alpha }_{2}. $

将式(4)、(5)代入式(7)，可以简化为

(8) $ {A_i}T_i^2+{B_i}{T_i}+{C_i} = 0; {\text{ }}i = 1,2,3. $

(9) $ {T_i} = \tan \;\left( {{{{\theta _i}}}/{2}} \right). $

(10) $\left. \begin{gathered} {A_i} = \sin {\eta _i}\sin {\alpha _1}{V_{ix}} - \cos {\eta _i}\sin {\alpha _1}{V_{iy}} - \cos {\alpha _1}{V_{iz}} - \cos {\alpha _2} , \\ {B_i} = - 2\cos {\eta _i}\sin {\alpha _1}{V_{ix}} - 2\sin {\eta _i}\sin {\alpha _1}{V_{iy}}, \\ {C_i} = \cos {\eta _i}\sin {\alpha _1}{V_{iy}} - \sin {\eta _i}\sin {\alpha _1}{V_{ix}} - \cos {\alpha _1}{V_{iz}} - \cos {\alpha _2}. \\ \end{gathered} \right\}$

由于CSPM具有3根输入轴重合的特征，当动平台绕 $ Z $轴旋转 $ \varphi $时，3个驱动杆也会旋转 $ \varphi $. 由此可知，在CSPM逆运动学的求解过程中，驱动杆旋向差异主要取决于动平台姿态角 $ {\varphi }_{1}\mathrm{、}{\varphi }_{2} $. 基于上述特性，提出基于分步计算CSPM逆运动学求解的新方法. 将动平台姿态分解为绕Z轴和绕X、Y轴旋转的2个部分，分别计算出2个部分姿态近端连杆所需的最优输入角度 $ {{\boldsymbol{\theta}}}_{i}^{{'}} $及 $ {{\boldsymbol{\theta}}}_{i}^{*} $；根据输入轴重合的特点 ${{\boldsymbol{\theta}}}_{i}^{*}=\left[\varphi_{3},\varphi_{3},\varphi _{3}\right]$，将 $ {{\boldsymbol{\theta}}}_{i}^{{'}}、{{\boldsymbol{\theta}}}_{i}^{*} $进行加法运算后，得到完整的逆运动学解.

1)根据CSPM的特点，将动平台姿态 $ \left({\varphi }_{1},{\varphi }_{2},{\varphi }_{3}\right) $分解为绕 $ X、Y $轴旋转的 $ \left({\varphi }_{1},{\varphi }_{2},0\right) $和绕 $ Z $轴旋转的 $ \left(0,0,{\varphi }_{3}\right) $，动平台姿态向量使用 $ {\boldsymbol{V}}_{i} $表示.

(11) $ {{\boldsymbol{V}}}_{i}=\left({{\boldsymbol{Q}}}_{1}\cdot{{\boldsymbol{Q}}}_{2}\right)\left[\begin{array}{c}\mathrm{sin}\;(\delta -{\eta }_{i})\mathrm{sin}\;\beta \\ \mathrm{cos}\;(\delta -{\eta }_{i})\mathrm{sin}\;\beta \\ -\mathrm{cos}\;\beta \end{array}\right]. $

(12) $ {{\boldsymbol{Q}}_{\text{1}}} = {{\boldsymbol{R}}_z}\left( {{\varphi _1}} \right){{\boldsymbol{R}}_y}\left( {{\varphi _2}} \right){{\boldsymbol{R}}_x}\left( 0 \right). $

(13) $ {{\boldsymbol{Q}}_2} = {{\boldsymbol{R}}_z}\left( 0 \right){{\boldsymbol{R}}_y}\left( 0 \right){{\boldsymbol{R}}_x}\left( {{\varphi _3}} \right). $

2)利用式(5)计算 $ \left({\varphi }_{1},{\varphi }_{2},0\right) $姿态下的动平台向量 $ {\boldsymbol{V}}_{i}^{{''}} $，

(14) $ {{\boldsymbol{V}}}_{i}^{''}={{\boldsymbol{R}}}_z(0){{\boldsymbol{R}}}_y({\varphi }_{2}){{\boldsymbol{R}}}_x({\varphi }_{1})\left[\begin{array}{c}\mathrm{sin}\;(\delta -{\eta }_{i})\mathrm{sin}\;\beta \\ \mathrm{cos}\;(\delta -{\eta }_{i})\mathrm{sin}\;\beta \\ -\mathrm{cos}\;\beta \end{array}\right]. $

3)根据式(10)计算：

(15) $ \begin{gathered} \left. \begin{array}{l} A_i^* = \sin {\eta _i}\sin {\alpha _1}V_{ix}^{''} - \cos {\eta _i}\sin {\alpha _1}V_{iy}^{''} - \cos {\alpha _1}V_{iz}^{''} - \cos {\alpha _2}, \\ B_i^* = - 2\cos {\eta _i}\sin {\alpha _1}V_{ix}^{''} - 2\sin {\eta _i}\sin {\alpha _1}V_{iy}^{''} , \\ C_i^* = \cos {\eta _i}\sin {\alpha _1}V_{iy}^{''} - \sin {\eta _i}\sin {\alpha _1}V_{ix}^{''} - \cos {\alpha _1}V_{iz}^{''} - \cos {\alpha _2} . \\ \end{array} \right\} \end{gathered} $

式中： $ V_{ix}^{{''}}、V_{iy}^{{''}}、V_{iz}^{{''}} $分别为 $ {\boldsymbol{V}}_{i}^{{''}} $的X、Y、Z轴分量.

4)结合式(9)，可得

(16) $ {{\boldsymbol{T}}_i} = \left[ {\begin{array}{*{20}{c}} {{T_{i1}}} \\ {{T_{i2}}} \end{array}} \right],\;i = 1,2,3. $

(17) $ \left. \begin{gathered} {T_{i1}} = \frac{{ - B_i^*+\sqrt {B{{_i^*}^2} - {\text{4}}A_i^*C_i^*} }}{{{\text{2}}A_i^*}} , \\ {T_{i2}} = \frac{{ - B_i^* - \sqrt {B{{_i^*}^2} - {\text{4}}A_i^*C_i^*} }}{{{\text{2}}A_i^*}}. \\ \end{gathered} \right\}$

5)利用式(9)、(10)，计算 $ {\boldsymbol{T}}_{i} $对应的关节转角 $ {\boldsymbol{\theta }}_{\mathit{i}}^{\mathit{{''}}} $.

(18) $ {\boldsymbol{\theta }}_i^{''} = \left[ {\begin{array}{*{20}{c}} {\theta _{i1}^{''}} \\ {\theta _{i2}^{''}} \end{array}} \right];\;i = 1,2,3. $

(19) $ \begin{gathered} \left. \begin{array}{l} \theta _{i1}^{''} = \arctan \;{T_{i1}} , \\ \theta _{i2}^{''} = \arctan\; {T_{i2}}. \\ \end{array} \right\} \end{gathered} $

6)由于每个驱动杆旋转角 $ {{\boldsymbol{\theta}} }_{i}^{{''}} $存在2个可能值， CSPM存在8种可能的逆解. 对于手腕康复装置而言，为了避免CSPM连杆发生碰撞及结构处于奇异位形，从每个驱动杆旋转角的2个可能值中选取绝对值较小的一个值，作为逆运动学的唯一解 ：

(20) $ {{\boldsymbol{\theta }}^ * } = [\theta _1^*,\;\theta _2^*,\;\theta _3^*], $

(21) $ \theta _i^* = \min \;({\text{abs}}\;(\theta _{i1}^{''}{\text{,}}\;\theta _{i2}^{''}{\text{))}}{\text{.}} $

式中：min表示取最小值，abs表示取绝对值.

7)当动平台绕 $Z$轴旋转 $ \left(\mathrm{0,0},{\varphi }_{3}\right) $时，近端连杆所需的输入角度为

(22) $ {{\boldsymbol{\theta }}^{'}} = [{\varphi _3},\;{\varphi _3},\;{\varphi _3}] . $

将 $ {{\boldsymbol{\theta}} }^{*} $与 $ {{\boldsymbol{\theta}} }^{{'}} $进行加法运算后，可得唯一的逆运动学解，即当CSPM动平台旋转角度 $ \left({\varphi }_{1},{\varphi }_{2},{\varphi }_{3}\right) $时，近端连杆输入角度为

(23) $ {\boldsymbol{\theta }} = [{\theta _1},\;{\theta _2},\;{\theta _3}]. $

(24) $ {\theta _i} = \theta _i^*+{\varphi _3},\;i = 1,2,3. $

对上述过程归纳整理后，得到如图4所示的唯一逆运动学解流程图.

以图2所示的共轴3RRR球面并联机构为例，机构参数为 $ {\alpha }_{1}={50}^{\circ},{\alpha }_{2}={57.6}^{\circ},\;\beta ={75}^{\circ},\;\delta ={60}^{\circ} $. 机构的装配方式为远端连杆均在近端连杆的左侧，如图3所示. 设动平台的位姿角为(30°, 15°, 20°)，利用式(5)求解得到

(25) $ {{\boldsymbol{V}}_1} = \left[ {\begin{array}{*{20}{c}} {0.576\;2} \\ {0.792\;5} \\ { - 0.199\;8} \end{array}} \right] , {{\boldsymbol{V}}_2} = \left[ {\begin{array}{*{20}{c}} { - 0.942\;4} \\ {0.239\;8} \\ {0.233\;3} \end{array}} \right] , {{\boldsymbol{V}}_3} = \left[ {\begin{array}{*{20}{c}} {0.069\;9} \\ { - 0.727\;1} \\ { - 0.683\;0} \end{array}} \right]. $

根据式(14)~(19)，求解得到

(26) $\begin{split} & {{\boldsymbol{\theta}} }_{1}^{{''}}\text=\left[\begin{array}{c} -{113.143\;9}^{\circ }\\ {1.108\;0}^{\circ }\end{array}\right]\text{，}{{\boldsymbol{\theta}} }_{2}^{{''}}\text=\left[\begin{array}{c}-{87.265\;4}^{\circ }\\ -{41.289\;8}^{\circ }\end{array}\right]\text{，}\\ &{{\boldsymbol{\theta}} }_{3}^{{''}}\text=\left[\begin{array}{c}-{154.548\;8}^{\circ }\\ {5.527\;2}^{\circ }\end{array}\right]. \end{split} $

根据式(20)，得到 $ {\boldsymbol{\theta }}^{*}=\left[{1.11}^{\circ},{-41.29}^{\circ},{5.53}^{\circ}\right] $. 将 $ {\alpha }_{3}={20}^{\circ} $代入步骤7），计算得到逆运动学解为 $ \boldsymbol{\theta }= \left[{21.11}^{\circ},{-21.29}^{\circ},{25.53}^{\circ}\right] $.

利用文献[14]的正运动学方法对逆运动学新解法进行验证，将得到的逆运动学解 ${\boldsymbol{\theta }}$代入式(4)，可得

(27) $ {{\boldsymbol{W}}_1} = \left[ {\begin{array}{*{20}{c}} { - 0.275\;9} \\ {0.714\;6} \\ { - 0.642\;8} \end{array}} \right] ,\;{{\boldsymbol{W}}_2} = \left[ {\begin{array}{*{20}{c}} { - 0.757\;2} \\ { - 0.116\;0} \\ { - 0.642\;8} \end{array}} \right] ,\;{{\boldsymbol{W}}_3} = \left[ {\begin{array}{*{20}{c}} {0.763\;7} \\ { - 0.059\;7} \\ { - 0.642\;8} \end{array}} \right]. $

根据CSPM动平台任意2个姿态向量之间的几何约束关系，存在

(28) $ {{\boldsymbol{V}}_i} \cdot {{\boldsymbol{V}}_j} = \cos\; {\alpha _3};\;i,j = 1,2,3,{\text{ }}i \ne j. $

式中： $ {\varphi }_{3} $为相邻姿态向量 $ {{\boldsymbol{V}}}_{i} $之间的夹角，有 $ {\alpha }_{3}= 2{{\rm{arcsin}}}\left({\rm{cos}}\;\beta {\rm{cos}}\;({{\text{π}}}/{6})\right) $ .

将式(5)、(27)代入式(7)，结合式(27)、(28)可得

(29) $ \left.\begin{array}{l} -0.275\;9{v}_{1x}+0.714\;6{v}_{1y}-0.642\;8{v}_{1z}-0.535\;8=0,\\ -0.752\;7{v}_{2x}-0.116\;0{v}_{2y}-0.642\;8{v}_{2z}-0.535\;8=0,\\ 0.763\;7{v}_{3x}-0.059\;7{v}_{3y}-0.642\;8{v}_{3z}-0.535\;8=0, \\ {v}_{1x}{v}_{2x}+{v}_{1y}{v}_{2y}+{v}_{1z}{v}_{2z}+0.398\;7=0,\\ {v}_{1x}{v}_{3x}+{v}_{1y}{v}_{3y}+{v}_{1z}{v}_{3z}+0.398\;7=0,\\ {v}_{2x}{v}_{3x}+{v}_{2y}{v}_{3y}+{v}_{2z}{v}_{3z}+0.398\;7=0, \\ {v}_{1x}^{2}+{v}_{1y}^{2}+{v}_{1z}^{2}=1,\\ {v}_{2x}^{2}+{v}_{2y}^{2}+{v}_{2z}^{2}=1,\\ {v}_{3x}^{2}+{v}_{3y}^{2}+{v}_{3z}^{2}=1. \end{array}\right\} $

式中 ： $ {v}_{1x}\mathrm{、}{v}_{1y}\mathrm{、}{v}_{1z} $为 $ {\boldsymbol{V}}_{1} $的 $ X\mathrm{、}Y\mathrm{、}Z $轴分量.

根据机构初始姿态下的装配方式，确定方程的初值为

(30) $ {{\boldsymbol{V}}_0} = \left[ {\begin{array}{*{20}{c}} {\begin{array}{*{20}{c}} 1, & { - 1}, &1, \end{array}}&{\begin{array}{*{20}{c}} { - 1}, & { - 1}, &1, \end{array}}&{\begin{array}{*{20}{c}} { - 1}, &1, &1 \end{array}} \end{array}} \right]. $

使用matlab中的fsolve函数，求解得到

(31) $ {{\boldsymbol{V}}_{\mathbf{1}}} = \left[ {\begin{array}{*{20}{c}} {0.576\;2} \\ {0.792\;5} \\ { - 0.199\;8} \end{array}} \right] , {{\boldsymbol{V}}_2} = \left[ {\begin{array}{*{20}{c}} { - 0.942\;2} \\ {0.240\;7} \\ {0.232\;9} \end{array}} \right] , {{\boldsymbol{V}}_3} = \left[ {\begin{array}{*{20}{c}} {0.068\;5} \\ { - 0.725\;5} \\ { - 0.684\;8} \end{array}} \right]. $

将式（31）的正运动学解与式（25）所示的初始值进行对比可知，两者的误差非常小，验证了所提方法的正确性.

手腕康复装置的奇异性是指CSPM处于某些特殊姿态时会失去固有刚度且机构可控性变差. 球面并联机构运动路径应该远离奇异位形，以免机构丧失自由度，降低机构的可控性.

对式(7)求导，得到

(32) $ {\boldsymbol{\dot W}}_{i}\cdot{{\boldsymbol{V}}}_{i}+{{\boldsymbol{W}}}_{i}\cdot{\boldsymbol{\dot V}}_{i}=0. $

(33) $ {{{\boldsymbol{\dot V}}_i}} = {{{\boldsymbol{\omega}} }} \times {{\boldsymbol{V}}_i} . $

(34) $ { {\boldsymbol{\dot W}} _i} = { {{\dot \theta}} _i} \cdot\left( {{\boldsymbol{U}}_i} \times {{\boldsymbol{W}}_i} \right) \cdot $

式中： $ \boldsymbol{\omega } $为终端转速， ${{{\dot \theta}} _i}$为驱动转角速度.

将式(33)、(34)代入式(32)，化简得到

(35) $ {\dot \theta}_{i}\cdot({{\boldsymbol{U}}}_{i}\times {{\boldsymbol{W}}_i})\cdot{{\boldsymbol{V}}}_{i}=({{\boldsymbol{W}}}_{i}\times {{\boldsymbol{V}}}_{i})\cdot{\boldsymbol{\omega}}. $

(36) $ {\boldsymbol{J\omega}} ={\boldsymbol{\dot \theta}}. $

式中：J为雅可比矩阵， ${\boldsymbol{J}} = {[{{\boldsymbol{J}}_1}{\mathbf{,}}\;{{\boldsymbol{J}}_2}{\mathbf{,}}\;{{\boldsymbol{J}}_3}]^{\rm{T}}}$，

(37) $ {{\boldsymbol{J}}}_{i}=\frac{{{\boldsymbol{W}}}_{i}\times {{\boldsymbol{V}}}_{i}}{({{\boldsymbol{U}}}_{i}\times {{\boldsymbol{W}}}_{i})\cdot{{\boldsymbol{V}}}_{i}},\;i=1,2,3. $

条件系数定义为

(38) $ \zeta ({\boldsymbol{J}}) = \frac{1}{{\kappa ({\boldsymbol{J}})}}. $

(39) $ \kappa ({\boldsymbol{J}}) = \left\| {\boldsymbol{J}} \right\|\left\| {{{\boldsymbol{J}}^{{{ - 1}}}}} \right\|. $

通过 $\zeta ({\boldsymbol{J}}) $可以判断机构姿态是否处于奇异位形. 当 $\zeta ({\boldsymbol{J}}) \to 0$时，CSPM姿态处于奇异位形附近；反之，CSPM姿态远离奇异位形. 定义 $\zeta_{\min } {({\boldsymbol{J}})}$来判定机构位姿点是否位于工作空间内，参考众多研究， 将 $\zeta_{\min } {({\boldsymbol{J}})}$定义为0.25.

根据文献[18]可知，CSPM连杆碰撞的原因是相邻驱动杆之间的夹角过大或者过小，本文采用欧式距离判断碰撞是否发生. 如图5所示，假设近端连杆和远端连杆的连接点为 $ {B}_{i}(i=\mathrm{1,2},3) $，动平台与上端连杆的连接点为 $ {C}_{i}(i=\mathrm{1,2},3) $. 通过几何距离测量可知，共轴上点 $ {A}_{i}、{B}_{i}、{C}_{i} $与 $ O $点的欧式距离都为100 mm，即 $ {OA}_{i}={OB}_{i}={OC}_{i}=100\;\mathrm{m}\mathrm{m} $，计算 $ {A}_{i}、{B}_{i}、{C}_{i} $在固定坐标系下的坐标为

(40) $ \left.\begin{gathered} {{\boldsymbol{A}}_i} = {{\boldsymbol{U}}_i} \cdot OA, \\ {{\boldsymbol{B}}_i} = {{\boldsymbol{W}}_i} \cdot OB, \\ {{\boldsymbol{C}}_i} = {{\boldsymbol{V}}_i} \cdot OC. \\ \end{gathered} \right\}$

CSPM在运动过程中只能发生以下2种类型的碰撞^[18]：1)在相邻2个驱动杆与远端连杆的连接点处发生碰撞；2)动平台与底座之间发生碰撞.

对于所研究的CSPM手腕康复装置而言，第1类碰撞是近端连杆末端之间发生碰撞. 设定碰撞点为 $ D $，如图6(a)所示，此时测量 ${B}_{i}{B}_{j}=22\;\mathrm{m}\mathrm{m}, i、j=\mathrm{1,2},3,\;i\ne j$，因此碰撞边界条件为

(41) $ {B_i}{B_j} \leqslant 22\;{\text{mm; }}i、j = 1,2,3, {\text{ }}i \ne j. $

(42) $ {B_i}{B_j} = \sqrt {OB_i^2 - OB_j^2} ;\;i、j = 1,2,3,\;i \ne j \cdot $

对于所研究的CSPM手腕康复装置而言，第2类碰撞是动平台与底座边缘之间发生了碰撞，碰撞点为 $ F $，如图6(b)所示. 此时测量 $ {\mathbf{\angle }A}_{i}O{C}_{i}={22}^{\circ } $， $A{C}_{i}=38\;{\rm{mm}}（i=\mathrm{1,2},3）$，碰撞的边界条件为

(43) $ A{C_i} \leqslant 38\;{\text{mm; }}i = 1,2,3. $

(44) $ A{C}_{i}=2\cdot OA\cdot \mathrm{sin}\left(\frac{\angle AO{C}_{i}}{2}\right),\;i=1,2,3. $

对于CSPM手腕康复装置而言，设定不发生奇异和碰撞的欧拉角 ${\bf\textit{φ}}$= $ \left[{\varphi }_{1},{\varphi }_{2},{\varphi }_{3}\right] $为姿态角工作空间 $\varOmega $. 采用蒙特卡洛法进行工作空间计算，采样点数为200 000，假设绕X轴转动的姿态角 $ {\varphi }_{1} $的采样范围为 $ \left[{-90}^{\circ},{90}^{\circ}\right] $；绕Y轴转动的姿态角 $ {\varphi }_{2} $的采样范围为 $\left[{-90}^{\circ},{90}^{\circ}\right]$；绕Z轴转动的姿态角 $ {\varphi }_{3} $的范围为 $\left[{-180}^{\circ},{180}^{\circ}\right]$. 以每一个采样点姿态角对应的 $\zeta ({\boldsymbol{J}})$是否大于0.25以及是否满足第1、2碰撞条件作为判断依据，若某一采样点姿态角满足上述条件，则判定该点属于欧拉角工作空间 $\varOmega $. 由于CSPM在工作空间可以绕旋转轴无限旋转，求解得到工作空间中绕 $ Z $轴转动角度 $ {\varphi }_{3} $仍为 $ \left[{-180}^{\circ},{180}^{\circ}\right] $，如图7所示.

为了确保中风患者康复训练时的手腕运动角度小于表1所示的极限运动角度，避免患者手腕受到损伤，将CSPM手腕康复装置的实际运动空间限定在以下范围： $ {\varphi }_{3} $转动角度为 $ \left[{-60}^{\circ},{60}^{\circ}\right] $， $ {\varphi }_{2} $转动角度为 $ \left[{-15}^{\circ},{30}^{\circ}\right] $, $ {\varphi }_{1} $转动角度为 $ \left[{-40}^{\circ},{52}^{\circ}\right] $. 如图8所示.

根据CSPM手腕康复装置的设计目标，使用者只能在固定轨迹下进行康复训练. 根据欧拉角和单位四元数之间存在的转换关系，可以将手腕实际工作的空间姿态点转化为单位四元数，结合所提的逆运动学新方法与单位四元数插补，实时得到平滑的曲线.

CSPM姿态 $ \left({\varphi }_{1},{\varphi }_{2},{\varphi }_{3}\right) $可以转换为单位四元数：

(45) $ \begin{split} {\boldsymbol{q}} = \left[ {\begin{array}{*{20}{c}} {\cos \dfrac{{{\varphi _3}}}{2}} \\ 0 \\ 0 \\ {\sin \dfrac{{{\varphi _3}}}{2}} \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {\cos \dfrac{{{\varphi _2}}}{2}} \\ 0 \\ {\sin \dfrac{{{\varphi _2}}}{2}} \\ 0 \end{array}} \right]\left[ {\begin{array}{*{20}{c}} {\cos \dfrac{{{\varphi _1}}}{2}} \\ {\sin \dfrac{{{\varphi _1}}}{2}} \\ 0 \\ 0 \end{array}} \right] = \\ \left[ {\begin{array}{*{20}{c}} {\cos \dfrac{{{\varphi _1}}}{2}\cos \dfrac{{{\varphi _2}}}{2}\cos \dfrac{{{\varphi _3}}}{2}+\sin \dfrac{{{\varphi _1}}}{2}\sin \dfrac{{{\varphi _2}}}{2}\sin \dfrac{{{\varphi _3}}}{2}} \\ {\sin \dfrac{{{\varphi _1}}}{2}\cos \dfrac{{{\varphi _2}}}{2}\cos \dfrac{{{\varphi _3}}}{2} - \cos \dfrac{{{\varphi _1}}}{2}\sin \dfrac{{{\varphi _2}}}{2}\sin \dfrac{{{\varphi _3}}}{2}} \\ {\cos \dfrac{{{\varphi _1}}}{2}\sin \dfrac{{{\varphi _2}}}{2}\cos \dfrac{{{\varphi _3}}}{2}+\sin \dfrac{{{\varphi _1}}}{2}\cos \dfrac{{{\varphi _2}}}{2}\sin \dfrac{{{\varphi _3}}}{2}} \\ {\cos \dfrac{{{\varphi _1}}}{2}\cos \dfrac{{{\varphi _2}}}{2}\sin \dfrac{{{\varphi _3}}}{2} - \sin \dfrac{{{\varphi _1}}}{2}\sin \dfrac{{{\varphi _2}}}{2}\cos \dfrac{{{\varphi _3}}}{2}} \end{array}} \right]. \end{split} $

$ \boldsymbol{q} $可以表示为

(46) $ {\boldsymbol{q}} = {\left[ {\begin{array}{*{20}{c}} {{q_{\rm{s}}}},&{q_{{\rm{v}}1}},&{q_{{\rm{v}}2}},&{{q_{{\rm{v}}3}}} \end{array}} \right]^{\rm{T}}} . $

单位四元数 $ \boldsymbol{q} $可以转化为欧拉角：

(47) $ \left[ {\begin{array}{*{20}{c}} {{\varphi _1}} \\ {{\varphi _2}} \\ {{\varphi _3}} \end{array}} \right] {\text{ = }} \left[ {\begin{array}{*{20}{c}} {{\rm{a}}\tan 2\;(2({q_{\rm{s}}}{q_{{\rm{v}}1}}+{q_{{\rm{v}}2}}{q_{{\rm{v}}3}}),1 - 2(q_{{\rm{v}}1}^2+q_{{\rm{v}}2}^2))} \\ {\arcsin\; (2({q_{\rm{s}}}{q_{{\rm{v}}2}} - {q_{{\rm{v}}1}}{q_{{\rm{v}}3}}))} \\ {{\rm{a}}\tan 2\;(2({q_{\rm{s}}}{q_{{\rm{v}}3}}+{q_{{\rm{v}}1}}{q_{{\rm{v}}2}}),1 - 2(q_{{\rm{v}}2}^2+q_{{\rm{v}}3}^2))} \end{array}} \right]. $

设定机构起始姿态为初始位姿，经过旋转变换 $ {\boldsymbol{q}}_{1} $，机构姿态为“姿态1”；在初始位形下经过旋转变换 $ {\boldsymbol{q}}_{2} $， $ {\boldsymbol{q}}_{1}\ne {\boldsymbol{q}}_{2} $机构姿态为 “姿态2”，如图9所示.

从图9可知，“姿态1”与“姿态2”可以通过单位四元数 $ \Delta \boldsymbol{q} $表示：

(48) $ \Delta {\boldsymbol{q}} = {{\boldsymbol{q}}_2}{\boldsymbol{q}}_1^{ - 1}. $

根据单位四元数性质 $ {\boldsymbol{q}}^{-1}={\boldsymbol{q}}^{*} $，得到“姿态1”到“姿态2”的插补表达式：

(49) $ {{\boldsymbol{q}}_t} = {({{\boldsymbol{q}}_2}{\boldsymbol{q}}_1^*)^t}{{\boldsymbol{q}}_1};\;t \in [0,1.0]. $

式中： $ {\boldsymbol{q}}_{1}^{*} $ 为 $ {\boldsymbol{q}}_{1} $的共轭四元数，t为时间.

将插值点 $ {\boldsymbol{q}}_{t} $使用式(47)转化为欧拉角，利用求解运动学唯一逆运动学解的方法对全部插值点进行求解，可以得到路径点的关节转角.

搭建如图10所示的实验装置. 选择初始位姿 $ {{\boldsymbol{\gamma}} }_{1}=\left[{0}^{\circ},{0}^{\circ},{0}^{\circ}\right] $为“姿态1”，在CSPM手腕康复机构实际工作空间中选取 $ {{\boldsymbol{\gamma}} }_{2}=\left[{5}^{\circ},{15}^{\circ},{30}^{\circ}\right] $为“姿态2”，代入式(45)转换为单位四元数：

(50) $ \left.\begin{gathered} {{\boldsymbol{q}}_1} = \left[ {\begin{array}{*{20}{c}} 1&0&0&0 \end{array}} \right] , \\ {{\boldsymbol{q}}_2} = \left[ {\begin{array}{*{20}{c}} {0.958\;2}&{ - 0.008\;0}&{0.137\;2}&{0.250\;9} \end{array}} \right] \cdot \\ \end{gathered}\right\} $

通过LabVIEW编写上位机程序，arduino编写下位机程序，下位机接收上位机发送的初始点和终点欧拉角度，使用所提的方法进行运动学逆求解及驱动步进电机. 通过查询相关资料，将实验运动时间设定为4 s，采用等值压缩方法，将运动时间点映射到 $ [0,1.0] $，并代入式(47)获得路径点姿态. 为了保证实验过程中速度更加平滑，引入T型曲线规划算法规划速度. 使用无线位姿传感器测量动平台姿态，反馈给上位机.

实验结果和理论结果如图11所示. 可知，理论计算结果与试验测试结果均为光滑曲线， $ {\varphi }_{1} $、 $ {\varphi }_{2} $的最大误差约为1°， $ {\varphi }_{3}\mathrm{的} $最大误差约为2.5°，误差均在合理范围内. 上述结果表明，所提的逆运动解方案是可行、可信的. 该实验过程中的误差产生原因如下. 1)实验为开环控制策略. 2)结构安装误差. 转动副在安装过程中未能消除间隙，导致了动平台姿态误差的产生. 3)零件弹性变形. 驱动杆和远端连杆发生了微小的弹性变形，造成了动平台的姿态误差. 图11所示的实验曲线比较光滑，表明CSPM手腕康复机构在运动过程中姿态变化平滑，为后续手腕康复机构的实时控制器设计提供了良好基础.

(1) 提出并验证了共轴3RRR并联机构逆运动学解的新方法，解决了球面手腕康复装置逆运动学解完备性和无解析解的问题.

(2)基于逆运动学解的新方法，以无连杆碰撞和无奇位形为约束条件，得到手腕康复装置欧拉角姿态空间. 通过分析手腕运动的特点，获得手腕康复装置实际工作空间： $ {\varphi }_{1}=\left[{-40}^{\circ},{52}^{\circ}\right] ，$ $ {\varphi }_{2}= \left[{-15}^{\circ},{30}^{\circ}\right] $， $ {\varphi }_{3}=\left[{-60}^{\circ},{60}^{\circ}\right] $.

(3) 提出并验证了基于四元数插补的实时轨迹规划方案， $ {\varphi }_{1} $、 $ {\varphi }_{2} $的最大误差为1°， $ {\varphi }_{3} $的最大误差为2.5°，为手腕康复机构的控制器设计提供了基础.