基于无传感器技术的人机交互变阻抗控制方法

A Variable Impedance Control Method for Human-Robot Interaction Based on Sensorless Technology

摘要: 为了提高人机交互的精确性和柔顺性，提出了一种基于无传感器技术的人机交互变阻抗控制方法。变阻抗控制的核心在于借助机械臂末端的速度和人机交互的交互力。首先，设计了一个速度控制器提升机械臂末端轨迹跟踪效果；同时，设计了一个补偿力机制来避免传统阻抗控制中交互力的测量。在此基础上，基于速度控制器和补偿力机制，构建了集成控制器，以实现精确的末端执行器轨迹跟踪，同时在线补偿相互作用力。运用李亚普诺夫理论证明，无传感器变阻抗控制下系统稳定。最后，仿真结果验证了基于无传感器技术的人机交互变阻抗控制方法的有效性。

Abstract: To improve the precision and compliance of human-robot interaction, a human-robot interaction variable impedance control method based on sensorless technology is proposed. The core of the variable impedance control is to utilize the velocity of the robot arm’s end-effector and the interaction force in human-robot interaction. First, a velocity controller is designed to enhance the trajectory tracking performance of the robot arm’s end-effector, and a compensation force mechanism is designed to avoid the measurement of the interaction force in traditional impedance control. Based on the velocity controller and the compensation force mechanism, an integrated controller is constructed to achieve accurate trajectory tracking of the end-effector and online compensation of the interaction force. The stability of the system under sensorless variable impedance control is proven using Lyapunov theory. Finally, simulation results verify the effectiveness of the human-robot interaction variable impedance control method based on sensorless technology.

1. 引言

人机交互的应用越来越广泛，例如在工业、演示和医疗康复领域[1] [2]。人机交互要求机器人能够对人的力做出响应[3]。因此，机器人需要更加人性化。为了实现人性化的交互，已经使用了各种可变阻抗控制方法来提高顺应性。在人机交互中，为了使机器人具有柔顺性，机器人的机械臂在交互过程中需要屈服于外部环境。在此基础上，有学者提出了基于速度的变阻抗自适应交互控制方法，根据力的误差来调整阻尼参数[4] [5]。Li团队[6] [7]根据机器人与人、环境之间的相互作用力和速度在线调整阻尼参数，保证系统的稳定性，使机器人与人、环境之间的交互更加柔顺。然而，上述变阻抗控制方法仅改变阻尼系数而保持惯性参数恒定，系统极点会产生右向偏移现象，导致瞬态响应过程的振荡衰减速率显著降低。

在人机交互中，保证机械臂末端位置跟踪的准确性至关重要。Gong等人[8]提出了一种结合基于能动性控制和阻抗控制的鲁棒控制方案，将挖掘机的耦合非线性动力学转化为开环港哈密顿模型，设计了改进的鲁棒扰动观测器，增强了液压挖掘机轨迹跟踪控制的鲁棒性。Foroutannia等人[9]提出了一种自适应模糊阻抗控制器和卷积神经网络结合的方法，实现机器人的期望轨迹跟踪。Kong的团队[10]、Sun的团队[11]与Jiang的团队[12]采用自适应阻抗控制来应对不确定性和外部干扰，但未考虑人机交互中的柔顺性问题。高鼎峰等人[13]设计了一种固定时间稳定的阻抗控制器，以实现障碍环境下空间机器人的柔顺控制。所设计的固定时间阻抗控制器可以实现机械臂良好的位置跟踪。然而，文献[8]-[13]不仅无法确保速度在可接受的范围内，同时不能确保人机交互环境中的顺应性[14]。

为了实现人机交互的柔顺性，力信号的采集至关重要。Huang等人[15]提出了一种刚度系数随机器人末端执行器位置变化的自适应阻抗控制方法，该方法将机器人的操作空间限制在一定的范围内，对操作者的顺应性较好。Huang等人[16]基于Q学习算法的自适应阻抗控制器，实现了阻抗自适应，保证了在人体肢体动力学完全未知的情况下，阻抗参数收敛于最优值。Xu等人[17]设计了一种利用参考轨迹学习的自适应阻抗控制，以实现位置误差和相互作用力的收敛。然而，参考文献[15]-[17]的所有系统需要力传感器，并假设力传感器是没有噪声的理想组件，这增加了系统的复杂性[18]。朱敬花等人[19]基于牛顿–欧拉法构建了机械臂由外向内的递归参数辨识方法，并运用位置阻抗控制原理对期望轨迹予以修正，进而设计出适用于参数不确定条件的非奇异终端滑模控制策略。虽然采用了无传感器技术，但是当遭遇更为复杂的动态环境，或是面临更高的控制精度需求时，需对其控制策略进行进一步的优化与完善，以提升其适应性。

为解决人机交互中轨迹跟踪精度与操作柔顺性的矛盾关系，本研究建立了免传感器的阻抗参数动态调节机制。本文的主要贡献如下：1) 结合交互力与速度，在线调整阻抗，以确保机器人能有效顺应人类的操作行为。2) 设计了一种补偿力替代变阻抗控制中的力传感器。

2. 无传感器技术的变阻抗控制的设计

2.1. 变阻抗控制的设计

在本节中，为了控制同一框架内的力和位置，设计了可变阻抗控制。然后根据力的变化调整阻抗参数，以保证机器人的顺应性和安全性。

机器人末端执行器的阻抗模型[20]为：

$M\left({\ddot{x}}_d-{\ddot{x}}_0\right)+B\left({\dot{x}}_d-{\dot{x}}_0\right)+K\left(x_d-x_0\right)=F_e$ (1)

其中，$M、B$ 和$K$ 分别为质量、阻尼和刚度。$x_0、{\dot{x}}_0$ 和${\ddot{x}}_0$ 分别是位置、速度和加速度。$x_d、{\dot{x}}_d$ 和${\ddot{x}}_d$ 分别是期望的位置、速度和加速度。$F_e$ 是人机交互力，通过第2.1节中设计的补偿模型进行重构替代，突破传统力传感器测量局限。

由于机器人在所有方向上解耦，等式(1)可以用小写符号重写为：

$m{\ddot{x}}_d+b{\dot{x}}_d=f_e$ (2)

其中，$m$ 和$b$ 分别为单方向质量和阻尼，$f_e$ 为单方向上人机交互力。则(2)的拉普拉斯变化为：

${\dot{X}}_d\left(s\right)=\frac{\frac{1}{b}}{\frac{m}{b}s+1}{F}_e\left(s\right)$ (3)

其中，${\dot{X}}_d\left(s\right)$ 和$F_e\left(s\right)$ 分别为${\dot{x}}_d$ 和$f_e$ 的拉普拉斯变换。$s$ 是拉普拉斯变量。式(3)可视为一阶系统，以驱动力为输入，速度为输出，其传递函数可表示为：

$H\left(s\right)=\frac{\frac{1}{b}}{\frac{m}{b}s+1}$ (4)

其中，$H\left(s\right)$ 是传递函数。由(4)可知，质量$m$ 影响状态的稳定性，质量和阻尼之比$m/b$ 影响系统的极点。当机器人与人交互时，阻抗参数(质量和阻尼)较大，机器人将更适合精确运动，因此需要更大的力来控制机器人。反之，阻抗参数越小，机器人越适合快速运动，因此，移动机器人所需的力越小。

根据以上分析，需要调整阻抗参数来控制机器人在不同的运动条件下，使机器人能够顺应人类的操作。为此，变阻抗是一种有效的方法。当人想要进行精确的运动时，阻抗参数需要调整得更高。相反，当人体想要进行加速运动时，则需要将阻抗参数调低。因此，阻抗参数是由人的意图来调整的。一般来说，人的意图可以分为三种模式，即加速、停止或逆转方向。由于反方向可以看作是停止后反方向的加速度，所以我们只需要关注加速度和停止。如果要进行加速运动，力的方向和速度应该是相同的；同时，阻抗参数要低。相反，如果要停止，力和速度的方向应该相反；同时，阻抗参数要高。

阻抗参数的变化取决于人的意图，即人需要加速时，阻尼减小。反之，人需要减速时，阻尼增加。因此，阻尼参数的调整策略如下：

$b=\{\begin{array}{l}{b}_0+b_1\left|f_{kl}\right|\\ b_0-b_2\left|f_{kl}\right|\end{array}$ (5)

其中，$b_0$ 为默认阻尼，$b_1$ 和$b_2$ 为可调正常数。$f_{kl}$ 为人机交互力的补偿力，将在第二节中详细设计。

由传递函数(4)可知，单一调节阻尼系数将引起极点右移，导致瞬态过程延长。为此建立质量–阻尼协同调节机制，保持$m/b$ 比值恒定，以稳定极点位置。具体实施策略为：加速阶段同步降低质量与阻尼参数，减速阶段则按比例提升二者数值。因此，质量调整如下：

$m=\frac{m_{0}b}{b_0}$ (6)

$m_0$ 为质量的默认值。由于在式(5)和式(6)中得到了质量$m$ 和阻尼$b$ 的调整方法。因此，将$m$ 和$b$ 代入式(2)中，得到期望位置和期望速度，这将在第二节中被使用。

2.2. 集成控制器的设计

接下来，将设计一个速度控制器来保证速度的稳定性，并设计一个补偿力来代替第一节(5)中的相互作用力来调节阻抗并补偿相互作用力。然后，设计了基于速度控制器和补偿力的集成控制器。

机器人的动力学为[21]：

$M_q\left(q\right)\ddot{q}+Cq\left(q,\dot{q}\right)\dot{q}+G_q\left(q\right)={\tau }_k+{\tau }_e$ (7)

其中，$M_q\left(q\right)$ 、$Cq\left(q,\dot{q}\right)$ 和$G_q\left(q\right)$ 分别为惯性矩阵、科里奥利项和离心项以及重力力矩。${\tau }_k$ 为控制力矩，${\tau }_e$ 为人机交互力矩。$q$ 、$\dot{q}$ 和$\ddot{q}$ 分别对应关节空间的位置、速度和加速度参数。

式(7)在笛卡尔空间可改写为：

$M_x\ddot{x}+C_x\dot{x}+G_x=f_k+f_e$ (8)

其中，$x$ 、$\dot{x}$ 和$\ddot{x}$ 分别是末端执行器的空间位置、速度和加速度。$M_x={\left(J{M}_q{J}^{{\rm T}}\right)}^{-1}$ 为笛卡尔空间中的矩阵，$C_x=J^{+{\rm T}}\left(C_q-M_q{J}^{+}\dot{J}\right)J^{+}$ 是笛卡尔空间中的科里奥利项，$G_x=J^{+{\rm T}}{G}_q$ 是重力。$J$ 为雅克比矩阵，$J^{+}$ 和$J^{+{\rm T}}$ 分别是雅可比矩阵的伪逆和雅可比矩阵专职的伪逆。$f_k=J^{+{\rm T}}{\tau }_k$ 为控制力，$f_e=J^{+{\rm T}}{\tau }_e$ 为人机接触力。

人–机器人相互作用力动力学模型[22]为：

$f_e=b{f}_h+K_h\left(t\right)e+B_h\left(t\right)\dot{e}$ (9)

其中，$b$ 为正常数，$f_h$ 为人的前馈力，$B_h$ 是人手的阻尼，$K_h$ 为人手的刚度。$e=x-x_d$ 和$\dot{e}=\dot{x}-{\dot{x}}_d$ 分别为位置误差和速度误差。参考速度${\dot{x}}_r$ 定义为：

${\dot{x}}_r={\dot{x}}_d-\mu e$ (10)

其中，$\mu$ 为一个常数因子，参考速度${\dot{x}}_r$ 由第一节期望速度${\dot{x}}_d$ 和位置误差$e$ 来设计的。则得到末端执行器速度$\dot{x}$ 与参考速度${\dot{x}}_r$ 之间的误差为：

$\theta =\dot{x}-{\dot{x}}_r=\dot{e}+\mu e$ (11)

为保证速度的稳定性，速度控制器可以设计为：

$f_{kv}=\underset{机器人动力学}{\underbrace{M_x{\ddot{x}}_r+C_x{\dot{x}}_r+G_x}}-\underset{参考速度误差补偿}{\underbrace{\zeta \left(1+\rho {\text{e}}^{-{\Vert {\theta }_r\Vert }^2}\right)\theta }}$ (12)

其中，$\zeta$ 和$\rho$ 均是正常数，$\text{e}$ 为指数收敛因子，$x_r$ 、${\dot{x}}_r$ 和${\ddot{x}}_r$ 分别为机器人的参考位置、速度和加速度。速度控制器(12)由两部分组成。第一部分$M_x{\ddot{x}}_r+C_x{\dot{x}}_r+G_x$ 表示机器人的动力学特性，第二部分$\zeta \left(1+\rho {\text{e}}^{-{\Vert {\theta }_r\Vert }^2}\right)\theta$ 用于补偿参考速度误差。

注意到人与机器人之间存在一种相互作用力。因此，如果仅考虑速度控制而不补偿相互作用力，系统会是不稳定的或无法达到理想控制效果。补偿力设计为：

$f_{kl}=-b{\widehat{f}}_h-{\widehat{K}}_h\left(t\right)e-{\widehat{B}}_h\dot{e}$ (13)

式中，$f_{kl}$ 为人机交互力的补偿力，${\widehat{f}}_h$ 、${\widehat{K}}_h$ 和${\widehat{B}}_h$ 分别为$f_h$ 、$K_h$ 和$B_h$ 的估计值。

在人机交互中，力信号的获取是必不可少的。然而，力传感器通常会产生测量噪声，引起控制效果。考虑倒力传感器的不完全性，利用轨迹的力来补偿相互作用力，使机器人末端执行器的速度更加稳定。

结合式(8)、(12)和(13)，可以得到：

$M_x\left(q\right)\dot{\theta }+C_x\left(q\right)\theta =-b{\tilde{f}}_h-{\tilde{K}}_h\left(t\right)e-{\tilde{B}}_h\dot{e}-\zeta \left(1+\rho {\text{e}}^{-{\Vert \theta \Vert }^2}\right)\theta$ (14)

其中，${\tilde{f}}_h={\widehat{f}}_h-f_h$ ，${\tilde{K}}_h={\widehat{K}}_h-K_h$ 和${\tilde{B}}_h={\widehat{B}}_h-B_h$ 分别表示操作者前馈力、刚度及阻尼的估计误差。

基于式(13)的补偿模型，构建参数自适应更新律为：

$\begin{array}{l}{\dot{\widehat{K}}}_h\left(t\right)={\Gamma }_K\left(\theta e^{{\rm T}}-\eta {\widehat{K}}_h\left(t\right)\right)\\ {\dot{\widehat{B}}}_h\left(t\right)={\Gamma }_B\left(\theta {\dot{e}}^{{\rm T}}-\eta {\widehat{B}}_h\left(t\right)\right)\\ {\dot{\widehat{f}}}_h\left(t\right)={\Gamma }_f\left(b\theta -\eta {\widehat{f}}_h\left(t\right)\right)\end{array}$ 15)

其中，${\Gamma }_K$ 、${\Gamma }_B$ 和${\Gamma }_f$ 为正定矩阵，$\eta =\frac{w}{\sqrt{\Vert \theta \Vert }}{\text{e}}^{-k\Vert \theta \Vert }$ ，$w$ 为正常数，$k$ 为正常数。

在$K_h$ 、$B_h$ 和$f_h$ 的自适应律(15)中$\frac{w}{\sqrt{\Vert \theta \Vert }}{\text{e}}^{-k\Vert \theta \Vert }$ ，当速度误差$\theta$ 变大时，$\eta$ 变小，反之亦然。${\text{e}}^{-k\Vert \theta \Vert }$ 使$\eta$ 对速度误差变化敏感，加快响应速度误差改变，$\sqrt{\Vert \theta \Vert }$ 避免$\eta$ 变化过于剧烈，保障系统稳定。因此，本文中$\eta$ 能更好地稳定补偿力。

结合速度控制器(12)和补偿力控制器(13)，设计了基于速度控制器和补偿力的集成控制器：

$f_k=f_{kv}+f_{kl}$ (16)

与文献[23]相比本文主要关注于通过变阻抗控制来实现机器人的柔顺性。此外，在速度控制器(12)中

加入速度补偿项$\zeta \left(1+\rho {\text{e}}^{-{\Vert {\theta }_{\text{r}}\Vert }^2}\right)\theta$ ，基于速度误差动态变化的非线性调节因子，使控制器对较大速度误差具

有更强的响应，对较小的速度误差更平滑。

3. 稳定性分析

如果机器人(8)由速度控制器(12)和补偿力控制器(13)组成的具有自适应律(15)的集成控制器(16)控制，则可以保证系统的稳定性和速度误差的有界性。

证明：为确保闭环系统的全局稳定性，需构造一个Lyapunov函数，其应同时包含速度跟踪误差和补偿力估计误差。为此，定义Lyapunov函数为：

$V\left(t\right)=V_v\left(t\right)+V_f\left(t\right)$ (17)

其中有：

$V_v\left(t\right)=\frac{1}{2}{\theta }^{{\rm T}}{M}_x\theta$ (18)

$V_v\left(t\right)=\frac{1}{2}{\theta }^{{\rm T}}{M}_x\theta$ 表征速度误差的能量项，反映机械臂末端速度跟踪的动态特性。

$V_f\left(t\right)=\frac{1}{2}tr\left({\tilde{f}}_h{}^{{\rm T}}\left(t\right){\Gamma }_f{}^{-{\rm T}}{\tilde{f}}_h\left(t\right)+{\tilde{B}}_h{}^{{\rm T}}\left(t\right){\Gamma }_B{}^{-{\rm T}}{\tilde{B}}_h\left(t\right)+{\tilde{K}}_h{}^{{\rm T}}\left(t\right){\Gamma }_K{}^{-{\rm T}}{\tilde{K}}_h\left(t\right)\right)$ (19)

其中$tr(\cdot )$ 代表矩阵的迹，$V_f\left(t\right)$ 与补偿力误差的能量项有关，目的是分析补偿力误差是否能够收敛，保证补偿力的稳定性。$V\left(t\right)=V_v\left(t\right)+V_f\left(t\right)$ 的组合形式能够同时约束动力学误差与参数估计误差，符合能量耗散原理，是典型的一类Lyapunov候选函数。式(18)的微分为：

${\dot{V}}_v\left(t\right)={\theta }^{{\rm T}}{M}_x\dot{\theta }+\frac{1}{2}{\theta }^{{\rm T}}{\dot{M}}_x\theta$ (20)

将(15)代入(20)，根据机器人动力学特性${\dot{M}}_x-2C_x$ 的斜对称性，可得${\theta }^{{\rm T}}{\dot{M}}_x\theta -2{\theta }^{{\rm T}}{C}_x\theta =0$ [24]。再结合式(14)的闭环动力学方程，进一步展开可以得到：

$\begin{array}{l}{V}_v\left(t\right)={\theta }^{{\rm T}}{M}_x\dot{\theta }+{\theta }^{{\rm T}}{C}_x\theta \\ =-{\theta }^{{\rm T}}{C}_x\theta -{\theta }^{{\rm T}}b{\tilde{f}}_h\left(t\right)-{\theta }^{{\rm T}}{\tilde{K}}_h\left(t\right)e-{\theta }^{{\rm T}}{\tilde{B}}_h\left(t\right)\dot{e}-{\theta }^{{\rm T}}\zeta \left(1+\rho {\text{e}}^{-{\Vert \theta \Vert }^2}\right)\theta +{\theta }^{{\rm T}}{C}_x\theta \\ =-b{\theta }^{{\rm T}}{\tilde{f}}_h\left(t\right)-{\theta }^{{\rm T}}{\tilde{K}}_h\left(t\right)e-{\theta }^{{\rm T}}{\tilde{B}}_h\left(t\right)\dot{e}-{\theta }^{{\rm T}}\zeta \left(1+\rho {\text{e}}^{-{\Vert \theta \Vert }^2}\right)\theta \end{array}$ (21)

式(19)的微分为：

$\begin{array}{l}{\dot{V}}_f\left(t\right)=\frac{V_f\left(t\right)-V_f\left(t-\delta \right)}{\delta }\\ =\frac{1}{2\delta }tr({\tilde{f}}_h{}^{{\rm T}}\left(t\right){\Gamma }_f^{-{\rm T}}{\tilde{f}}_h\left(t\right)+{\tilde{B}}_h{}^{{\rm T}}\left(t\right){\Gamma }_B^{-{\rm T}}{\tilde{B}}_h\left(t\right)+{\tilde{K}}_h{}^{{\rm T}}\left(t\right){\Gamma }_K^{-{\rm T}}{\tilde{K}}_h\left(t\right)-{\tilde{f}}_h{}^{{\rm T}}\left(t-\delta \right){\Gamma }_f^{-{\rm T}}{\tilde{f}}_h\left(t-\delta \right)\\ -{\tilde{B}}_h{}^{{\rm T}}\left(t-\delta \right){\Gamma }_B^{-{\rm T}}{\tilde{B}}_h\left(t-\delta \right)-{\tilde{K}}_h{}^{{\rm T}}\left(t-\delta \right){\Gamma }_K^{-{\rm T}}{\tilde{K}}_h\left(t-\delta \right))\end{array}$ (22)

其中，$\delta$ 为时间常数，${\tilde{f}}_h\left(t-\delta \right)$ 、${\tilde{B}}_h\left(t-\delta \right)$ 和${\tilde{K}}_h\left(t-\delta \right)$ 分别为延迟$\delta$ 秒的${\tilde{f}}_h\left(t\right)$ 、${\tilde{B}}_h\left(t\right)$ 和${\tilde{K}}_h\left(t\right)$ 的值。

现定义$\Delta {\tilde{f}}_h\left(t\right)={\tilde{f}}_h\left(t\right)-{\tilde{f}}_h\left(t-\delta \right)=\delta {\dot{\tilde{f}}}_h\left(t\right)$ ，则可以得到：

$\begin{array}{l}\frac{1}{2\delta }tr\left({\tilde{f}}_h{}^{{\rm T}}\left(t\right){\Gamma }_f^{-{\rm T}}{\tilde{f}}_h\left(t\right)-{\tilde{f}}_h{}^{{\rm T}}\left(t-\delta \right){\Gamma }_f^{-{\rm T}}{\tilde{f}}_h\left(t-\delta \right)\right)\\ =\frac{1}{2\delta }tr\left({\left({\tilde{f}}_h\left(t\right)-{\tilde{f}}_h\left(t-\delta \right)\right)}^{{\rm T}}{\Gamma }_f{}^{-{\rm T}}\left(2{\tilde{f}}_h\left(t\right)-{\tilde{f}}_h\left(t\right)+{\tilde{f}}_h\left(t-\delta \right)\right)\right)\\ =\frac{1}{2\delta }tr\left(2\Delta {\tilde{f}}_h{}^{{\rm T}}\left(t\right){\Gamma }_f^{-{\rm T}}{\tilde{f}}_h\left(t\right)-\Delta {\tilde{f}}_h{}^{{\rm T}}\left(t\right){\Gamma }_f^{-{\rm T}}\Delta {\tilde{f}}_h\left(t\right)\right)\end{array}$ (23)

将自适应律(15)代入式(23)：

$\begin{array}{l}\frac{1}{2\delta }tr\left(2\Delta {\tilde{f}}_h{}^{{\rm T}}\left(t\right){\Gamma }_f^{-{\rm T}}{\tilde{f}}_h\left(t\right)-\Delta {\tilde{f}}_h{}^{{\rm T}}\left(t\right){\Gamma }_f^{-{\rm T}}\Delta {\tilde{f}}_h\left(t\right)\right)\\ =b{\theta }^{{\rm T}}{\tilde{f}}_h\left(t\right)-tr\left(\eta {\widehat{f}}_h{}^{{\rm T}}\left(t\right){\tilde{f}}_h\left(t\right)+\frac{1}{2\delta }\Delta {\tilde{f}}_h{}^{{\rm T}}\left(t\right){\Gamma }_f^{-{\rm T}}\Delta {\tilde{f}}_h\left(t\right)\right)\end{array}$ (24)

同样，对于式(22)中第二项，即：$\frac{1}{2\delta }tr\left({\tilde{K}}_h{}^{{\rm T}}\left(t\right){\Gamma }_{\text{K}}^{-{\rm T}}{\tilde{K}}_h\left(t\right)-{\tilde{K}}_h{}^{{\rm T}}\left(t-\delta \right){\Gamma }_{\text{K}}^{-{\rm T}}{\tilde{K}}_h\left(t-\delta \right)\right)$，可以得到：

$\frac{1}{2\delta }tr\left({\tilde{K}}_h{}^{{\rm T}}\left(t\right){\Gamma }_{\text{K}}^{-{\rm T}}{\tilde{K}}_h\left(t\right)-{\tilde{K}}_h{}^{{\rm T}}\left(t-\delta \right){\Gamma }_{\text{K}}^{-{\rm T}}{\tilde{K}}_h\left(t-\delta \right)\right)={\theta }^{{\rm T}}{\tilde{K}}_h\left(t\right)e-tr\left(\eta {\widehat{K}}_h{}^{{\rm T}}\left(t\right){\tilde{K}}_h\left(t\right)+\frac{1}{2\delta }\Delta {\tilde{K}}_h^{{\rm T}}\left(t\right){\Gamma }_{\text{K}}^{-{\rm T}}\Delta {\tilde{K}}_h\left(t\right)\right)$(25)

对于式(22)中第三项，可得到：

$\begin{array}{l}\frac{1}{2\delta }tr\left({\tilde{B}}_h^{\text{Τ}}\left(t\right){\Gamma }_{\text{B}}^{-{\rm T}}{\tilde{B}}_h\left(t\right)-{\tilde{B}}_h^{\text{Τ}}\left(t-\delta \right){\Gamma }_{\text{B}}^{-{\rm T}}{\tilde{B}}_h\left(t-\delta \right)\right)\\ ={\theta }^{{\rm T}}{\tilde{B}}_h\left(t\right)\dot{e}-tr\left(\eta {\widehat{B}}_h^{\text{Τ}}\left(t\right){\tilde{B}}_h\left(t\right)+\frac{1}{2\delta }\Delta {\tilde{B}}_h^{{\rm T}}\left(t\right){\Gamma }_{\text{B}}^{-{\rm T}}\Delta B_h\left(t\right)\right)\end{array}$ (26)

接下来，定义$\phi \left(t\right)={\left[f_h^{{\rm T}}\left(t\right),vec{\left(B_h\left(t\right)\right)}^{{\rm T}},vec{\left(K_h\left(t\right)\right)}^{{\rm T}}\right]}^{{\rm T}}$ ，$\tilde{\phi }\left(t\right)=\widehat{\phi }\left(t\right)-\phi \left(t\right)$ ，$\Delta \phi \left(t\right)=\phi \left(t\right)-\phi \left(t-\delta \right)$ 。其中，$\widehat{\phi }\left(t\right)$ 为$\phi \left(t\right)$ 的估计，$\phi \left(t-\delta \right)$ 为延迟$\delta$ 秒的$\phi \left(t\right)$ 的值。其中，$vec(\cdot )$ 代表列向量运算符，假设$A=\left[\begin{array}{cc}{a}_{11}& a_{12}\\ a_{21}& a_{22}\end{array}\right]$ ，则$vec\left(A\right)={\left[\begin{array}{cccc}{a}_{11}& a_{21}& a_{12}& a_{22}\end{array}\right]}^{{\rm T}}$ 。

由式(22)~(26)可得到：

$\begin{array}{l}{\dot{V}}_f\left(t\right)=b{\theta }^{{\rm T}}{\tilde{f}}_h\left(t\right)+{\theta }^{{\rm T}}{\tilde{B}}_h\left(t\right)\dot{e}+{\theta }^{{\rm T}}{\tilde{K}}_h\left(t\right)e-tr\left(\eta {\widehat{\phi }}^{{\rm T}}\left(t\right)\tilde{\phi }\left(t\right)+\frac{1}{2\delta }\Delta {\tilde{\phi }}^{{\rm T}}\left(t\right){\Gamma }^{-{\rm T}}\Delta \tilde{\phi }\left(t\right)\right)\\ =b{\theta }^{{\rm T}}{\tilde{f}}_h\left(t\right)+{\theta }^{{\rm T}}{\tilde{B}}_h\left(t\right)\dot{e}+{\theta }^{{\rm T}}{\tilde{K}}_h\left(t\right)e-tr\left(\eta {\tilde{\phi }}^{{\rm T}}\left(t\right)\tilde{\phi }\left(t\right)+\eta {\phi }^{{\rm T}}\left(t\right)\tilde{\phi }\left(t\right)+\frac{1}{2\delta }\Delta {\tilde{\phi }}^{{\rm T}}\left(t\right){\Gamma }^{-{\rm T}}\Delta \tilde{\phi }\left(t\right)\right)\end{array}$ (27)

其中，${\Gamma }^{-{\rm T}}=diag\left({\Gamma }_f^{-{\rm T}},{\rm I}\otimes {\Gamma }_B^{-{\rm T}},{\rm I}\otimes {\Gamma }_K^{-{\rm T}}\right)$ ，${\rm I}$ 代表单位矩阵，$diag$ 代表对角矩阵。假设$A=\left[\begin{array}{cc}{a}_{11}& a_{12}\\ a_{21}& a_{22}\end{array}\right]$ ，$B=\left[\begin{array}{cc}{b}_{11}& b_{12}\\ b_{21}& b_{22}\end{array}\right]$ ，则$A\otimes B=\left[\begin{array}{cccc}{a}_{11}{b}_{11}& a_{11}{b}_{12}& a_{12}{b}_{11}& a_{12}{b}_{12}\\ a_{11}{b}_{21}& a_{11}{b}_{22}& a_{12}{b}_{21}& a_{12}{b}_{22}\\ a_{21}{b}_{11}& a_{21}{b}_{12}& a_{22}{b}_{11}& a_{22}{b}_{12}\\ a_{21}{b}_{21}& a_{21}{b}_{22}& a_{22}{b}_{21}& a_{22}{b}_{22}\end{array}\right]$ 。

将式(21)代入式(27)：

$\begin{array}{l}\dot{V}\left(t\right)={\dot{V}}_v\left(t\right)+{\dot{V}}_f\left(t\right)\\ =-b{\theta }^{{\rm T}}{\tilde{f}}_h\left(t\right)-{\theta }^{{\rm T}}{\tilde{K}}_h\left(t\right)e-{\theta }^{{\rm T}}{\tilde{B}}_h\left(t\right)\dot{e}-{\theta }^{{\rm T}}\zeta \left(1+\rho {\text{e}}^{-{\Vert \theta \Vert }^2}\right)\theta +b{\theta }^{{\rm T}}{\tilde{f}}_h\left(t\right)+{\theta }^{{\rm T}}{\tilde{B}}_h\left(t\right)\dot{e}\\ +{\theta }^{{\rm T}}{\tilde{K}}_h\left(t\right)e-tr\left(\eta {\tilde{\phi }}^{{\rm T}}\left(t\right)\tilde{\phi }\left(t\right)+\eta {\phi }^{{\rm T}}\left(t\right)\tilde{\phi }\left(t\right)+\frac{1}{2\delta }\Delta {\tilde{\phi }}^{{\rm T}}\left(t\right){\Gamma }^{-{\rm T}}\Delta \tilde{\phi }\left(t\right)\right)\\ =-\frac{1}{2\delta }tr\left(\Delta {\tilde{\phi }}^{{\rm T}}\left(t\right){\Gamma }^{-{\rm T}}\Delta \tilde{\phi }\left(t\right)\right)-tr\left({\theta }^{{\rm T}}\zeta \left(1+\rho {\text{e}}^{-{\Vert \theta \Vert }^2}\right)\theta +\eta {\tilde{\phi }}^{{\rm T}}\left(t\right)\tilde{\phi }\left(t\right)+\eta {\phi }^{{\rm T}}\left(t\right)\tilde{\phi }\left(t\right)\right)\end{array}$ (28)

$\dot{V}(t)<0$ 时，充分条件如下：

${\theta }^{{\rm T}}\zeta \left(1+\rho {\text{e}}^{-{\Vert \theta \Vert }^2}\right)\theta ++\eta {\tilde{\phi }}^{{\rm T}}\left(t\right)\tilde{\phi }\left(t\right)+\eta {\phi }^{{\rm T}}\left(t\right)\tilde{\phi }\left(t\right)\ge -\zeta \left(1+\rho {\text{e}}^{-{\Vert \theta \Vert }^2}\right){\Vert \theta \Vert }^2+\eta {\Vert \tilde{\phi }\left(t\right)\Vert }^2-\eta \Vert \phi \left(t\right)\Vert \Vert \tilde{\phi }\left(t\right)\Vert \ge 0$ (29)

由于$\zeta$ 、$\rho$ 和$\eta$ 均为正常数，根据式(29)，$\Vert \theta \Vert$ 和$\Vert \tilde{\phi }\left(t\right)\Vert$ 的界为：

$\begin{array}{l}\Vert \theta \Vert \le \sqrt{\frac{M_1}{\zeta \left(1+\rho \right)}}={\kappa }_1\\ \Vert \tilde{\phi }\left(t\right)\Vert \ge \frac{-M_2+\sqrt{M_2^2+4\eta M_3}}{2\eta }={\kappa }_2\end{array}$ (30)

其中，$M_1=\eta {\Vert \tilde{\phi }\left(t\right)\Vert }^2-\eta \Vert \phi \left(t\right)\Vert \Vert \tilde{\phi }\left(t\right)\Vert$ ，$M_2=-\eta \Vert \phi \left(t\right)\Vert$ 和$M_3=-\zeta \left(1+\rho {\text{e}}^{-{\Vert \theta \Vert }^2}\right){\Vert \theta \Vert }^2$ 。

从上述分析中可以看出，则$\dot{V}\left(t\right)<0$ 在以下趋于：

$\kappa =\left\{\left(\Vert \theta \Vert ,\Vert \tilde{\phi }\left(t\right)\Vert \right),{\kappa }_1\wedge {\kappa }_2\right\}$ (31)

由式(30)可知，速度误差$\theta$ 收敛，有界为${\kappa }_1$ 。补偿力误差$\tilde{\phi }\left(t\right)$ 有界为${\kappa }_2$ 。由于$\dot{V}\left(t\right)<0$ ，保证了系统的总能量不断减少，从而确保误差不会增大。同时，式(15)自适应律中的$\eta =\frac{w}{\sqrt{\Vert \theta \Vert }}{\text{e}}^{-k\Vert \theta \Vert }$ 随着速度误差变化进行自适应调节，这保证了补偿力误差不会持续增大，而是被调节到一个稳定的范围。

4. 仿真实验分析

本节中给出仿真实验，以验证所提出方法的有效性。考虑如图1所示的两连杆刚性机器。

Figure 1. Two-link robot manipulator

图1. 双连杆机器人机械臂

其机器人的动力学如下：

$M_q\left(q\right)=\left[\begin{array}{cc}{m}_1+m_2+2m_3\cos q_2& m_2+m_3\cos q_2\\ m_2+m_3\cos q_2& m_2\end{array}\right]$ ，$C_q\left(q,\dot{q}\right)=\left[\begin{array}{cc}-m_3{q}_2\sin q_2& -m_3\left(q_1+q_2\right)\sin q_2\\ m_3{q}_1\sin q_2& 0\end{array}\right]$ ，

$G_q\left(q\right)=\left[\begin{array}{c}{m}_{4}g\cos q_1+m_{5}g\cos \left(q_1+q_2\right)\\ m_{5}g\cos \left(q_1+q_2\right)\end{array}\right]$ ，$J=\left[\begin{array}{cc}-L_1\sin q_1-L_2\sin \left(q_1+q_2\right)& -L_2\sin \left(q_1+q_2\right)\\ L_1\cos q_1+L_2\cos \left(q_1+q_2\right)& L_2\cos \left(q_1+q_2\right)\end{array}\right]$ ，

其中，$q_1$ 和$q_2$ 分别为关节1与关节2的关节角度，$L_1$ 和$L_2$ 分别为关节1和关节2的长度。 $M={\left[\begin{array}{ccccc}{m}_1& m_2& m_3& m_4& m_5\end{array}\right]}^{{\rm T}}={\left[\begin{array}{ccccc}1.88& 0.8& 1.12& 3.85& 1.56\end{array}\right]}^{{\rm T}}$ ，$L_1=L_2=1$ 和$g=9.8{\text{ m/s}}^2$ 。期望轨迹$x_d$ (单位为$\text{m}$ )设置为：

$x_d={\left[\begin{array}{cc}1+0.9\sin \left(0.1t\right)& 0\end{array}\right]}^{{\rm T}}$

$f_{e}x$ 和$f_{e}y$ 分别是$x$ 方向和$y$ 方向上的交互力，$f_e={\left(\begin{array}{cc}{f}_{e}x& f_{e}y\end{array}\right)}^{{\rm T}}$ ，默认阻尼为$b_0=115\text{ Ns}\cdot {\text{m}}^{-1}$ ，默认质量为$m_0=56\text{ kg}$ 。其他相关参数设置见表1所示。

Table 1. Relative parameters

表1. 相关参数

$b_1$	$b_2$	$b$	${\Gamma }_K$	${\Gamma }_B$	${\Gamma }_f$	$\zeta$	$\rho$	$k$
$3$	$6.5$	$4$	$425I$	$38I$	$4I$	$9$	$0.6$	$0.3$

图2为阻抗控制下的交互力与轨迹跟踪效果，图3为本文方法下的补偿力与位置跟踪。图2中的交互力为传统阻抗控制下机械臂末端的交互力，图3中的补偿力为本文方法下的机械臂末端的补偿力。从中可知，在X方向上和Y方向上，本文方法的补偿力与交互力的变化趋势保持一致性。验证了无传感器方案的可行性。因此，实施时需确保补偿力矢量方向与实际交互力吻合，幅值允许存在适度偏差。

图2 (传统阻抗控制)和图3 (本文方法)分别示出了x方向和y方向上的期望位置和实际位置。通过人的意图获得期望的位置。由图2和图3可知，机器人可以在x方向上跟踪期望的位置。只有在交互力发生突变的瞬间，才会出现位置误差。由图3可知，在交互力突变时刻($2.5\text{ s}<t<8\text{ s}$ )，机器人末端实际轨迹可以在突然变化后快速跟踪所需位置。在x方向和y方向上的轨迹跟踪误差在图4和图5中示出。

Figure 2. Interaction force and position tracking under impedance control

图2. 阻抗控制下交互力与位置跟踪

Figure 3. Compensation force and position tracking under the proposed method

图3. 本文方法下的补偿力与位置跟踪

如图4与图5分别展示了X与Y方向上的轨迹跟踪误差动态响应。通过对比可得到以下结论：

在X方向上，传统阻抗控制中，在突加交互力阶段($2.5\text{ s}<t<8\text{ s}$ )，期望轨迹与实际轨迹误差峰值达到0.3359 m，且稳态阶段($t>12\text{ s}$ )误差持续波动。这表明单一阻尼调节导致系统极点右移，引发振荡衰减速率降低。本文方法轨迹误差峰值显著降低至0.08 m。误差曲线的平滑性提升(见图5所示)验证了式(12)的速度补偿机制，通过动态调整参考速度误差，有效抑制了交互力引起的超调。

在Y方向上，传统阻抗控制方法的动态跟踪过程中，误差峰值为0.3 m ($t=4\text{ s}$ )，误差波动显著。本文所提方法的稳态误差稳定在0.02 m。误差的快速收敛(见图5所示)表明力补偿机制(式(13))通过模型估计替代直接测量，有效抵消了低频环境扰动。

Figure 4. Trajectory tracking error in the X and Y directions under impedance control

图4. 阻抗控制下X与Y方向上的轨迹跟踪误差

Figure 5. Trajectory tracking error in the X and Y directions under the proposed method

图5. 本文方法下X与Y方向上的轨迹跟踪误差

图6 (传统阻抗控制)与图7 (本文方法)展示了X与Y方向上的力矩跟踪性能，进一步验证了无传感器方案的有效性：

X方向上力矩(图6传统阻抗控制方法的${\tau }_x$ 与图7本文方法的${\tau }_x$ )：传统阻抗控制方法的力矩响应呈现显著波动(峰值22.5 N∙m)。在突加交互力阶段($2.5\text{ s}<t<3.5\text{ s}$ )，力矩维持震荡，表明单一阻尼调节导致能量过冲与极点不稳定。本文方法中的力矩动态响应更加平滑(见图7)。这一优化归因于质量–阻尼保持${m_0/b}_0$ 比值恒定，从而稳定系统极点位置，避免振荡发散。

Y方向上力矩：传统阻抗控制方法的力矩稳态误差为2.1 N∙m，本文方法中力矩稳态误差为1.1 N∙m。补偿力机制通过自适应更新律(见式15所示)在线估计交互力，避免了传感器噪声的直接引入。

Figure 6. Torque in the X and Y directions under impedance control

图6. 阻抗控制下X与Y方向上的力矩

Figure 7. Torque in the X and Y directions under the proposed method

图7. 本文方法下X与Y方向上的力矩

实验数据揭示了误差与力矩的动态耦合关系：

瞬态阶段的协同抑制：在突加交互力阶段($2.5\text{ s}<t<3.5\text{ s}$ )，传统方法中力矩的剧烈波动(见图6)与误差峰值(见图4)直接相关，表明能量过冲导致跟踪失稳。而本文方法通过速度补偿的动态调整，在误差增大的瞬间同步降低力矩幅值(见图7)，实现了误差与力矩的协同抑制。例如，X方向误差峰值0.08 m对应力矩为20 N∙m，较传统方法的误差0.335 m对应的力矩22.5 N∙m有明显抑制。

稳态阶段的优化：本文方法下力矩信号的幅值与波动均显著降低(见图7)，与误差的稳态收敛(见图5)形成闭环验证。

根据实验数据可知，增加速度补偿与力补偿之后，X方向上与Y方向上的轨迹跟踪精度均有明显的提升。因此本文所提出的基于无传感器技术的人机交互变阻抗控制方法能够有效提升机械臂末端的轨迹跟踪精度，采用补偿力对力传感器进行替代是可行的。

5. 结语

在人机交互场景下，本研究致力于提升轨迹跟踪性能与拓展无传感器技术在阻抗控制中的应用。研究结果显示，速度补偿机制的引入显著提升了机械臂末端期望轨迹的跟踪效果，使机械臂在期间能更精确地契合期望轨迹，有效抑制了轨迹跟踪中的最大位置误差。与此同时，本研究还探讨了力补偿在无传感器阻抗控制的可行性。传统阻抗控制方法通常依赖力传感器获取外部力信息，但力传感器的使用存在成本高、安装复杂等问题。通过理论分析与实验验证，证实了力补偿可有效模拟力传感器的功能，实现对外部力的估计与补偿。

基金项目

国家重点研发计划资助项目(2020YFC2007502)。