2. Impedance control of SISO mechanical systems

In this section, we would like to consider the impedance control of a 1-DOF linear mechanical system for simplicity. Let m, b, and k be the mass, damping, and stiffness of the system, respectively, and the equation of motion can be represented as

(1) \begin{equation} m\ddot{x}+b\dot{x}+kx=u-f_{\mathrm{ext}} \end{equation}

where $x\in \Re$ is the position, $u\in \Re$ is the control input and $f_{\mathrm{ext}}\in \Re$ is the external force from the environment. A target impedance is specified as the desired behavior of the system during the interaction of the system and the environment [Reference Huang and Chien3–Reference Al-Shuka, Leonhardt, Zhu, Song, Ding and Li7]

(2) \begin{equation} M\!\left(\ddot{x}-\ddot{x}_{d}\right)+B\!\left(\dot{x}-\dot{x}_{d}\right)+K\!\left(x-x_{d}\right)=-f_{\mathrm{ext}} \end{equation}

where M, B, and K are positive constants to give proper dynamics, and $x_{d}\in \Re$ is the desired trajectory for the motion of the mechanical system. To design controller u in (1) so that the closed-loop system behaves like (2), we may firstly find $\ddot{x}$ from (2) as

(3) \begin{equation} \ddot{x}=\ddot{x}_{d}-M^{-1}B\!\left(\dot{x}-\dot{x}_{d}\right)-M^{-1}K\!\left(x-x_{d}\right)-M^{-1}f_{\mathrm{ext}} \end{equation}

Substitute (3) into (1), and Hogan’s impedance controller can be found to be

(4a) \begin{equation} \begin{array}{l} u=\underset{\text{known}}{\underbrace{m\ddot{x}_{d}+mM^{-1}B\dot{x}_{d}+mM^{-1}Kx_{d}}}\\\\ \qquad +\underset{\text{state feedback}}{\underbrace{\left(b-mM^{-1}B\right)\!\dot{x}+\left(k-mM^{-1}K\right)\!x}}+\underset{\text{force feedback}}{\underbrace{\left(1-mM^{-1}\right)f_{\mathrm{ext}}}} \end{array} \end{equation}

Thus, to implement this controller, we need to know all plant parameters $(m,b,k)$ , target impedance coefficients $(M,B,K)$ , system states $(x,\dot{x})$ , and the external force $f_{\mathrm{ext}}$ . Due to the difficulties in the force feedback as stated in the previous section, a force sensorless algorithm is desirable. Most of the force sensorless designs [Reference Murakami, Yu and Ohnishi8–Reference Roveda, Shahid, Iannacci and Piga25] replaced force sensors by force estimators or observers that are generally with complex structures and are involved in their mathematical justifications. It was mentioned in ref. [Reference Al-Shuka, Leonhardt, Zhu, Song, Ding and Li7] that for the simple SISO case here if it is possible to select $M=m$ , then the force feedback term $(1-mM^{-1})f_{\mathrm{ext}}$ in (4a) vanishes and a force sensorless impedance controller is obtained as

(4b) \begin{equation} \begin{array}{l} u=m\ddot{x}_{d}+mM^{-1}B\dot{x}_{d}+mM^{-1}Kx_{d}\\[5pt] \qquad +\left(b-mM^{-1}B\right)\dot{x}+\left(k-mM^{-1}K\right)x \end{array} \end{equation}

With this controller, the mechanical system (1) will behave like the target impedance (2) without using the force sensor, and it is no need for constructing any force estimator. So, the implementation of this impedance controller is much simplified. However, the condition $M=m$ is not always achievable even for the simple SISO case here simply because selections of M relate to the eigen-structure of the target impedance and the closed-loop dynamics as well. Therefore, the presence of the force feedback term in (4a) is unavoidable.

Let us consider the case when $M\neq m$ , and we would like to use some intuitive way to scale the coefficient of (2) to try to get an equivalent effect of $M=m$ . This is done by firstly rewrite (2) as a monic equation

\begin{equation*} \left(\ddot{x}-\ddot{x}_{d}\right)+M^{-1}B\!\left(\dot{x}-\dot{x}_{d}\right)+M^{-1}K\!\left(x-x_{d}\right)=-M^{-1}f_{\mathrm{ext}} \end{equation*}

Then, multiplying m to the both sides to have

(5) \begin{equation} m\!\left(\ddot{x}-\ddot{x}_{d}\right)+mM^{-1}B\!\left(\dot{x}-\dot{x}_{d}\right)+mM^{-1}K\!\left(x-x_{d}\right)=-mM^{-1}f_{\mathrm{ext}} \end{equation}

It seems that the inertia of the target impedance (5) becomes m which is exactly the inertia of the system equation (1). Can the force feedback term be canceled this way? Let us find $\ddot{x}$ from (5) first to compute the impedance controller as we had done in (3) and (4a), but the result is again (3), and the impedance controller goes back to (4a). As a result, we still have the same force feedback term, and the force feedback is unavoidable. Therefore, scaling of the coefficients in the target impedance is infeasible in giving the condition $M=m$ to work properly. This is simply because (2) and (5) are the same equation, and hence, we still obtain the same impedance controller (4a) which requires the force feedback.

It is seen that the leading coefficient of (5) is the same as the system inertia, but the impedance controller for (1) based on the target impedance (5) will still go back to (4a) which needs force feedback. Let us look at the differences in (2) and (5) to find any possible way to cancel out the force feedback terms in the controller. It is obvious that the coefficient sets $(M,B,K)$ in (2) and $(m,mM^{-1}B,mM^{-1}K)$ in (5) will give identical eigenvalues, and hence, the two target impedances share the same dynamics. The only difference now is the coefficients in the external force terms in (2) and (5): one is with $-1$ and the other with $-mM^{-1}$ . With $-1$ , it will result in cancelation of the external force term, and with $-mM^{-1}$ , it will go back to (4a) which requires force feedback. Hence, let us try to modify the coefficient of the external force term in (5) to $-1$ as

(6a) \begin{equation} m\!\left(\ddot{x}-\ddot{x}_{d}\right)+mM^{-1}B\!\left(\dot{x}-\dot{x}_{d}\right)+mM^{-1}K\!\left(x-x_{d}\right)=-f_{\mathrm{ext}} \end{equation}

Now, we would like to find the impedance controller of (1) based on the target impedance (6a). The result is exactly (4b) which is a force sensorless impedance controller even though the condition $M=m$ is not satisfied. Hence, a force sensorless impedance controller without using force estimators for the SISO case can be designed in two steps:

Step 1: Scaling the coefficients of the target impedance

Step 2: Modify the coefficient of the force term to $-1$

After this modification, the target impedance becomes (6a), and the force feedback term in the impedance controller is canceled.

Further observations can be done by writing (6a) into

(6b) \begin{equation} M\!\left(\ddot{x}-\ddot{x}_{d}\right)+B\!\left(\dot{x}-\dot{x}_{d}\right)+K\!\left(x-x_{d}\right)=-Mm^{-1}f_{\mathrm{ext}} \end{equation}

It is seen that (2) and (6b) share the coefficient set $(M,B,K)$ , and hence, the eigen-structure is identical which implies that the dynamics will be the same even though one will result in an impedance controller with force feedback and the other one without. Hence, we have proved that the cancelation of the force feedback term in the impedance controller this way will not destroy closed-loop stability. During the free space tracking phase (i.e., $f_{\mathrm{ext}}=0$ ), (4a) becomes (4b) and the tracking performance will be the same. In the compliant motion phase (i.e., $f_{\mathrm{ext}}\neq 0$ ), the coefficient $Mm^{-1}$ in (6b) will change the contact force. However, by proper selection of M, the contact force can be adjusted, and this may be useful in implementation.

3. Impedance control of robot manipulators

The new method developed in the previous section for SISO mechanical systems can be extended to the impedance control of robot manipulators to have a force sensorless design. Let us consider a robot manipulator described by the model

(7) \begin{equation} \mathbf{D}\!\left(\mathbf{q}\right)\!\ddot{\mathbf{q}}+\mathbf{C}\!\left(\mathbf{q},\dot{\mathbf{q}}\right)\!\dot{\mathbf{q}}+\mathbf{g}\!\left(\mathbf{q}\right)=\boldsymbol{\tau }-\mathbf{J}^{T}\mathbf{F}_{\mathrm{ext}} \end{equation}

where $\mathbf{q}\in \Re ^{n}$ is a vector of generalized coordinates, $\mathbf{D}(\mathbf{q})\in \Re ^{n\times n}$ is a positive definite inertia matrix, $\mathbf{C}(\mathbf{q},\dot{\mathbf{q}})\dot{\mathbf{q}}\in \Re ^{n}$ is a vector of centrifugal and Coriolis forces, $\mathbf{g}(\mathbf{q})\in \Re ^{n}$ is a vector of gravitational forces, $\boldsymbol{\tau }\in \Re ^{n}$ is the control torque vector, $\mathbf{J}(\mathbf{q})\in \Re ^{n\times n}$ is the non-singular Jacobian matrix, and $\mathbf{F}_{\mathrm{ext}}\in \Re ^{n}$ is the external force vector at the end-effector. Let $\mathbf{x}\in \Re ^{n}$ be the position vector of (7) in the Cartesian space, and we may use the Jacobian matrix to relate the velocity vector $\dot{\mathbf{x}}$ and $\dot{\mathbf{q}}$ in the joint space by using the relationship $\dot{\mathbf{x}}=\mathbf{J}\dot{\mathbf{q}}$ . Define the matrices $\mathbf{D}_{x}(\mathbf{x})=\mathbf{J}^{-T}\mathbf{D}\mathbf{J}^{-1}, \mathbf{C}_{x}(\mathbf{x},\dot{\mathbf{x}})=\mathbf{J}^{-T}(\mathbf{C}-\mathbf{D}\mathbf{J}^{-1}\dot{\mathbf{J}})\mathbf{J}^{-1}$ , and $\mathbf{g}_{x}(\mathbf{x})=\mathbf{J}^{-T}\mathbf{g}$ , and the system model in the Cartesian space becomes [Reference Craig2, Reference Huang and Chien3, Reference Spong, Hutchinson and Vidyasagar5]

(8) \begin{equation} \mathbf{D}_{x}\!\left(\mathbf{x}\right)\ddot{\mathbf{x}}+\mathbf{C}_{x}\!\left(\mathbf{x},\dot{\mathbf{x}}\right)\!\dot{\mathbf{x}}+\mathbf{g}_{x}\!\left(\mathbf{x}\right)=\mathbf{J}^{-T}\boldsymbol{\tau }-\mathbf{F}_{\mathrm{ext}} \end{equation}

We would like to design an impedance controller to (8) such that the closed-loop dynamics follows the target impedance

(9) \begin{equation} \mathbf{M}\!\left(\ddot{\mathbf{x}}-\ddot{\mathbf{x}}_{d}\right)+\mathbf{B}\!\left(\dot{\mathbf{x}}-\dot{\mathbf{x}}_{d}\right)+\mathbf{K}\!\left(\mathbf{x}-\mathbf{x}_{d}\right)=-\mathbf{F}_{\mathrm{ext}} \end{equation}

where $\mathbf{M},\mathbf{B},\mathbf{K}\in \Re ^{n\times n}$ are diagonal matrices representing the desired apparent inertia, damping, and stiffness, respectively, and $\mathbf{x}_{d}\in \Re ^{n}$ is the desired trajectory of the end-effector. If all parameters in (8) are available, then Hogan’s impedance controller can be designed to be

(10) \begin{align} \boldsymbol{\tau }&=\underset{\text{known}}{\underbrace{\mathbf{J}^{T}\mathbf{D}_{x}\ddot{\mathbf{x}}_{d}+\mathbf{J}^{T}\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{B}\dot{\mathbf{x}}_{d}+\mathbf{J}^{T}\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{K}\mathbf{x}_{d}}}\nonumber\\[5pt] &\quad +\underset{\text{state feedback}}{\underbrace{\mathbf{J}^{T}\!\left(\mathbf{C}_{x}-\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{B}\right)\dot{\mathbf{x}}+\mathbf{J}^{T}\mathbf{g}_{x}-\mathbf{J}^{T}\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{Kx}}}\\[5pt] &\quad +\underset{\text{force feedback}}{\underbrace{\mathbf{J}^{T}\!\left(\mathbf{I}-\mathbf{D}_{x}\mathbf{M}^{-1}\right)\mathbf{F}_{\mathrm{ext}}}} \nonumber \end{align}

Similar to the SISO case in the previous section, implementation of this controller requires the knowledge of all plant parameters in $(\mathbf{D}_{x},\mathbf{C}_{x},\mathbf{g}_{x})$ , target impedance coefficient matrices $(\mathbf{M},\mathbf{B},\mathbf{K})$ , system states $(\mathbf{x},\dot{\mathbf{x}})$ , and the external force vector $\mathbf{F}_{\mathrm{ext}}$ . Again, if we may select M such that $\mathbf{D}_{x}\mathbf{M}^{-1}=\mathbf{I}$ , then the force feedback term in (10) will be canceled, and it becomes a force sensorless impedance controller. Since $\mathbf{D}_{x}$ is a function of the robot states and M is a constant matrix, the condition $\mathbf{D}_{x}\mathbf{M}^{-1}=\mathbf{I}$ cannot be satisfied in general. Therefore, a force sensorless impedance control may not be achieved simply by selecting M to cancel the force feedback term in (10) directly. In the previous section, a force sensorless impedance controller was proposed for SISO mechanical systems in two steps. Let us extend the design to the control of robot manipulators as

Step 1: Scaling the coefficients of the target impedance (9) into

\begin{equation*} \mathbf{D}_{x}\!\left(\ddot{\mathbf{x}}-\ddot{\mathbf{x}}_{d}\right)+\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{B}\!\left(\dot{\mathbf{x}}-\dot{\mathbf{x}}_{d}\right)+\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{K}\!\left(\mathbf{x}-\mathbf{x}_{d}\right)=-\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{F}_{\mathrm{ext}} \end{equation*}

Step 2: Modify the coefficient of the external force term to $-\mathbf{I}$

(11) \begin{equation} \mathbf{D}_{x}\!\left(\ddot{\mathbf{x}}-\ddot{\mathbf{x}}_{d}\right)+\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{B}\!\left(\dot{\mathbf{x}}-\dot{\mathbf{x}}_{d}\right)+\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{K}\!\left(\mathbf{x}-\mathbf{x}_{d}\right)=-\mathbf{F}_{\mathrm{ext}} \end{equation}

Now, the leading coefficient matrix is the same as the inertia matrix in (8), and the right-hand sides of (11) and (9) are identical. Then, the impedance controller of (8) with the given target impedance (11) is derived to be exactly (10) but without the force feedback term as

(12) \begin{equation} \begin{array}{l} \boldsymbol{\tau }=\underset{\text{known}}{\underbrace{\mathbf{J}^{T}\mathbf{D}_{x}\ddot{\mathbf{x}}_{d}+\mathbf{J}^{T}\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{B}\dot{\mathbf{x}}_{d}+\mathbf{J}^{T}\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{K}\mathbf{x}_{d}}}\\\\ \qquad +\underset{\text{state feedback}}{\;\underbrace{\mathbf{J}^{T}\!\left(\mathbf{C}_{x}-\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{B}\right)\dot{\mathbf{x}}+\mathbf{J}^{T}\mathbf{g}_{x}-\mathbf{J}^{T}\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{Kx}}} \end{array} \end{equation}

This is a force sensorless impedance controller without using the force estimator. With this impedance controller, the system in (8) will behave like the target impedance in (11) which can also be written as

(13) \begin{equation} \mathbf{M}\!\left(\ddot{\mathbf{x}}-\ddot{\mathbf{x}}_{d}\right)+\mathbf{B}\!\left(\dot{\mathbf{x}}-\dot{\mathbf{x}}_{d}\right)+\mathbf{K}\!\left(\mathbf{x}-\mathbf{x}_{d}\right)=-\mathbf{M}\mathbf{D}_{x}^{-1}\mathbf{F}_{\mathrm{ext}} \end{equation}

It shares the same dynamics as the one in (9), since they have the same eigen-structure. However, their right-hand sides are different that will induce different contact forces. During the free space tracking phase, there is no external force induced, and hence, (10) and (12) are identical and the two impedance controllers will give the same closed-loop dynamics which implies the same free space tracking performance. On the other hand, during the compliant motion phase, the external force will be induced with different magnitudes due to the coefficient I in (9) and $\mathbf{M}\mathbf{D}_{x}^{-1}$ in (13). By proper selection of M, it is possible to make the force induced by the proposed controller be larger or smaller than the traditional one. It is seen that the closed-loop stability is not destroyed in both the free space tracking and compliant motion phases, and there is no need for any force estimators.

4. Adaptive impedance control of robot manipulators

In the previous section, we have proposed a force sensorless impedance controller without using force estimators under the assumption that all the system parameters are available. If there are uncertainties in (8), we would like to develop a force sensorless adaptive impedance controller in this section. An intuitive design of the adaptive version of the impedance controller is to replace matrices $(\mathbf{D}_{x},\mathbf{C}_{x},\mathbf{g}_{x})$ in (12) by using their estimates $(\hat{\mathbf{D}}_{x},\hat{\mathbf{C}}_{x},\hat{\mathbf{g}}_{x})$ as

\begin{equation*} \boldsymbol{\tau }=\mathbf{J}^{T}\hat{\mathbf{D}}_{x}\ddot{\mathbf{x}}_{d}+\mathbf{J}^{T}\hat{\mathbf{D}}_{x}\mathbf{M}^{-1}\mathbf{B}\dot{\mathbf{x}}_{d}+\mathbf{J}^{T}\hat{\mathbf{D}}_{x}\mathbf{M}^{-1}\mathbf{K}\mathbf{x}_{d}+\mathbf{J}^{T}\!\left(\hat{\mathbf{C}}_{x}-\hat{\mathbf{D}}_{x}\mathbf{M}^{-1}\mathbf{B}\right)\!\dot{\mathbf{x}}+\mathbf{J}^{T}\hat{\mathbf{g}}_{x}-\mathbf{J}^{T}\hat{\mathbf{D}}_{x}\mathbf{M}^{-1}\mathbf{Kx} \end{equation*}

This controller can then be substituted into (8) to get the error dynamics, and the update laws are able to be derived based on the Lyapunov stability theory. However, the adaptive impedance controller thus obtained requires the feedback of the joint accelerations, and there will be singularity problem during the estimation of $\hat{\mathbf{D}}_{x}$ [Reference Huang and Chien3]. The Slotine and Li’s modification can be applied with the support of the model reference design [Reference Huang and Chien3] to avoid the feedback of joint accelerations as well as the singularity problem, but the force sensorless structure will be destroyed in the computation of the reference model states. Hence, the traditional adaptive designs are not feasible here.

Let us rewrite (12) into

\begin{equation*} \boldsymbol{\tau }=\mathbf{J}^{T}\left\{\left(\mathbf{D}_{x}\ddot{\mathbf{x}}_{d}+\mathbf{C}_{x}\dot{\mathbf{x}}+\mathbf{g}_{x}\right)-\left[\underset{=-\mathbf{D}_{x}\!\left(\ddot{\mathbf{x}}-\ddot{\mathbf{x}}_{d}\right)-\mathbf{F}_{ext}\;\text{due to }\left(11\right)}{\underbrace{\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{B}\!\left(\dot{\mathbf{x}}-\dot{\mathbf{x}}_{d}\right)+\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{K}\!\left(\mathbf{x}-\mathbf{x}_{d}\right)}}\right]\right\} \end{equation*}

The term $\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{B}(\dot{\mathbf{x}}-\dot{\mathbf{x}}_{d})+\mathbf{D}_{x}\mathbf{M}^{-1}\mathbf{K}(\mathbf{x}-\mathbf{x}_{d})$ above can be replaced by the term $-\mathbf{D}_{x}(\ddot{\mathbf{x}}-\ddot{\mathbf{x}}_{d})-\mathbf{F}_{\mathrm{ext}}$ due to (11) to give

\begin{align*} \boldsymbol{\tau } & = \mathbf{J}^{T}\left\{\left(\mathbf{D}_{x}\ddot{\mathbf{x}}_{d}+\mathbf{C}_{x}\dot{\mathbf{x}}+\mathbf{g}_{x}\right)+\mathbf{D}_{x}\!\left(\ddot{\mathbf{x}}-\ddot{\mathbf{x}}_{d}\right)+\mathbf{F}_{\mathrm{ext}}\right\}\\[5pt] & =\mathbf{J}^{T}\left\{\left(\mathbf{D}_{x}\ddot{\mathbf{x}}+\mathbf{C}_{x}\dot{\mathbf{x}}+\mathbf{g}_{x}+\mathbf{F}_{\mathrm{ext}}\right)-\left[\underset{=0\;\text{due to }\left(13\right)}{\underbrace{\mathbf{M}\!\left(\ddot{\mathbf{x}}-\ddot{\mathbf{x}}_{d}\right)+\mathbf{B}\!\left(\dot{\mathbf{x}}-\dot{\mathbf{x}}_{d}\right)+\mathbf{K}\!\left(\mathbf{x}-\mathbf{x}_{d}\right)+\mathbf{M}\mathbf{D}_{x}^{-1}\mathbf{F}_{\mathrm{ext}}}}\right]\right\} \end{align*}

This can be rearranged into

\begin{equation*} \boldsymbol{\tau }=\mathbf{J}^{T}\left[\underset{\equiv \mathbf{h}\left(t\right)}{\underbrace{\mathbf{D}_{x}\ddot{\mathbf{x}}+\mathbf{C}_{x}\dot{\mathbf{x}}+\mathbf{g}_{x}+\mathbf{F}_{\mathrm{ext}}-\mathbf{M}\ddot{\mathbf{x}}-\mathbf{M}\mathbf{D}_{x}^{-1}\mathbf{F}_{\mathrm{ext}}}}+\mathbf{M}\ddot{\mathbf{x}}_{d}-\mathbf{B}\!\left(\dot{\mathbf{x}}-\dot{\mathbf{x}}_{d}\right)-\mathbf{K}\!\left(\mathbf{x}-\mathbf{x}_{d}\right)\right] \end{equation*}

Define $\mathbf{h}(t)=\mathbf{D}_{x}\ddot{\mathbf{x}}+\mathbf{C}_{x}\dot{\mathbf{x}}+\mathbf{g}_{x}+\mathbf{F}_{\mathrm{ext}}-\mathbf{M}\ddot{\mathbf{x}}-\mathbf{M}\mathbf{D}_{x}^{-1}\mathbf{F}_{\mathrm{ext}}$ as a vector containing uncertain plant parameters, acceleration feedback terms, and force feedback terms. This uncertain vector is time-varying without knowing its variation bounds, and hence, the traditional adaptive strategies as well as the conventional robust designs are infeasible to give proper treatment. Here, we would like to propose a function approximation [Reference Huang and Chien3]-based adaptive controller

(14) \begin{equation} \boldsymbol{\tau }=\mathbf{J}^{T}\left[\hat{\mathbf{h}}+\mathbf{M}\ddot{\mathbf{x}}_{d}-\mathbf{B}\!\left(\dot{\mathbf{x}}-\dot{\mathbf{x}}_{d}\right)-\mathbf{K}\!\left(\mathbf{x}-\mathbf{x}_{d}\right)\right] \end{equation}

where $\hat{\mathbf{h}}$ is an estimate of h. It is seen that the realization of (14) does not require the knowledge of system parameters, and it is a force sensorless algorithm. Substituting (14) into (8) and the closed-loop system becomes in the form

(15) \begin{equation} \mathbf{M}\!\left(\ddot{\mathbf{x}}-\ddot{\mathbf{x}}_{d}\right)+\mathbf{B}\!\left(\dot{\mathbf{x}}-\dot{\mathbf{x}}_{d}\right)+\mathbf{K}\!\left(\mathbf{x}-\mathbf{x}_{d}\right)=\left(\hat{\mathbf{h}}-\mathbf{h}\right)-\mathbf{M}\mathbf{D}_{x}^{-1}\mathbf{F}_{\mathrm{ext}} \end{equation}

Therefore, if we may find a proper update law for $\hat{\mathbf{h}}$ to have $\hat{\mathbf{h}}\rightarrow \mathbf{h}$ , then (15) will converge to (13) asymptotically and the force cancelation is achievable. To find the update law, let us define $\mathbf{e}=\mathbf{x}-\mathbf{x}_{d}$ as the trajectory tracking error, and (15) can be rewritten as

(16) \begin{equation} \ddot{\mathbf{e}}+\mathbf{M}^{-1}\mathbf{B}\dot{\mathbf{e}}+\mathbf{M}^{-1}\mathbf{Ke}=\mathbf{M}^{-1}\left(\hat{\mathbf{h}}-\mathbf{h}\right)-\mathbf{D}_{x}^{-1}\mathbf{F}_{\mathrm{ext}} \end{equation}

Let $\mathbf{x}_{e}=[\mathbf{e}^{T}\,\dot{\mathbf{e}}^{T}]^{T}\in \Re ^{2n}$ be the state vector, and then, (16) can be represented in the state space form

(17) \begin{equation} \dot{\mathbf{x}}_{e}=\mathbf{A}_{e}\mathbf{x}_{e}-\mathbf{B}_{e}\!\left(\hat{\mathbf{h}}-\mathbf{h}\right)-\mathbf{B}_{e}\mathbf{M}\mathbf{D}_{x}^{-1}\mathbf{F}_{\mathrm{ext}} \end{equation}

where

\begin{equation*} \mathbf{A}_{e}=\left[\begin{array}{c@{\quad}c} \mathbf{0} & \mathbf{I}_{n}\\[5pt] -\mathbf{M}^{-1}\mathbf{K} & -\mathbf{M}^{-1}\mathbf{B} \end{array}\right]\in \Re ^{2n\times 2n} \qquad \mathbf{B}_{e}=\left[\begin{array}{c} \mathbf{0}\\[5pt] \mathbf{M}^{-1} \end{array}\right]\in \Re ^{2n\times n} \end{equation*}

The vectors h and $\hat{\mathbf{h}}$ can be approximated by using the function approximation technique [Reference Huang and Chien3, Reference Chien and Huang26] as

\begin{equation*} \mathbf{h}=\mathbf{W}^{T}\mathbf{z}+\boldsymbol{\varepsilon } \end{equation*}

\begin{equation*} \hat{\mathbf{h}}=\hat{\mathbf{W}}^{T}\mathbf{z} \end{equation*}

where $\mathbf{W}\in \Re ^{n\beta \times n}$ is the weighting matrix, $\hat{\mathbf{W}}\in \Re ^{n\beta \times n}$ is its estimate, $\mathbf{z}(t)\in \Re ^{n\beta \times 1}$ is a vector of basis functions, and $\boldsymbol{\varepsilon }\in \Re ^{n}$ is the approximation error vector. Suppose a large enough number $\beta$ is used and the approximation error can be ignored, then (17) may be further written into

(18) \begin{equation} \dot{\mathbf{x}}_{e}=\mathbf{A}_{e}\mathbf{x}_{e}-\mathbf{B}_{e}{\tilde{\boldsymbol{{W}}}}{}^{T}\mathbf{z}-\mathbf{B}_{e}\mathbf{M}\mathbf{D}_{x}^{-1}\mathbf{F}_{\mathrm{ext}} \end{equation}

where ${\tilde{\boldsymbol{{W}}}}=\hat{\mathbf{W}}-\mathbf{W}$ is the approximation error of W. Consider the Lyapunov function like candidate

\begin{equation*} V=\frac{1}{2}\mathbf{x}_{e}^{T}\mathbf{P}\mathbf{x}_{e}+\frac{1}{2}Tr\!\left({\tilde{\boldsymbol{{W}}}}{}^{T}\boldsymbol{\Gamma }{\tilde{\boldsymbol{{W}}}}\right) \end{equation*}

where $Tr(.)$ is the trace operation, $\boldsymbol{\Gamma }\in \Re ^{n\beta \times n\beta }$ is a positive definite weighting matrix, and $\mathbf{P}=\mathbf{P}^{T}\in \Re ^{2n\times 2n}$ is a positive definite matrix satisfying the Lyapunov equation $\mathbf{A}_{e}^{T}\mathbf{P}+\mathbf{P}\mathbf{A}_{e}=-\mathbf{Q}$ for a given positive definite and symmetric matrix $\mathbf{Q}\in \Re ^{2n\times 2n}$ . Taking the time derivative of $V$ along the trajectory of (18), we have

\begin{align*} \dot{V} & =\frac{1}{2}\mathbf{x}_{e}^{T}\!\left(\mathbf{A}_{e}^{T}\mathbf{P}+\mathbf{P}\mathbf{A}_{e}\right)\mathbf{x}_{e}-\mathbf{x}_{e}^{T}\mathbf{P}\mathbf{B}_{e}{\tilde{\boldsymbol{{W}}}}{}^{T}\mathbf{z}-Tr\!\left({\tilde{\boldsymbol{{W}}}}{}^{T}\boldsymbol{\Gamma }\dot{\hat{\mathbf{W}}}\right)-\mathbf{x}_{e}^{T}\mathbf{P}\mathbf{B}_{e}\mathbf{M}\mathbf{D}_{x}^{-1}\mathbf{F}_{\mathrm{ext}}\\[5pt] & =-\frac{1}{2}\mathbf{x}_{e}^{T}\mathbf{Q}\mathbf{x}_{e}-Tr\left[{\tilde{\boldsymbol{{W}}}}{}^{T}\!\left(\mathbf{z}\mathbf{x}_{e}^{T}\mathbf{P}\mathbf{B}_{e}+\boldsymbol{\Gamma }\dot{\hat{\mathbf{W}}}\right)\right]-\mathbf{x}_{e}^{T}\mathbf{P}\mathbf{B}_{e}\mathbf{M}\mathbf{D}_{x}^{-1}\mathbf{F}_{\mathrm{ext}} \end{align*}

Therefore, the update law can be selected to be

(19) \begin{equation} \dot{\hat{\mathbf{W}}}=-\boldsymbol{\Gamma }^{-1}\mathbf{z}\mathbf{x}_{e}^{T}\mathbf{P}\mathbf{B}_{e} \end{equation}

then we may have

(20) \begin{equation} \dot{V}=-\frac{1}{2}\mathbf{x}_{e}^{T}\mathbf{Q}\mathbf{x}_{e}-\mathbf{x}_{e}^{T}\mathbf{P}\mathbf{B}_{e}\mathbf{M}\mathbf{D}_{x}^{-1}\mathbf{F}_{\mathrm{ext}} \end{equation}

Hence, in the free space tracking phase, there will be no contact force and (20) implies uniform boundedness of $\mathbf{x}_{e}$ and ${\tilde{\boldsymbol{{W}}}}$ . Square integrability of $\mathbf{x}_{e}$ can also be proved by the derivation

\begin{align*} \int _{0}^{\infty }\!\left(\sqrt{\mathbf{Q}}\mathbf{x}_{e}\right)^{T}\!\left(\sqrt{\mathbf{Q}}\mathbf{x}_{e}\right)dt & =\int _{0}^{\infty }\mathbf{x}_{e}^{T}\mathbf{Q}\mathbf{x}_{e}dt\\[5pt] & =-2\int _{0}^{\infty }\dot{V}dt=-2\left[V\!\left(\infty \right)-V\!\left(0\right)\right]\lt \infty \end{align*}

Together with boundedness of $\dot{\mathbf{x}}_{e}$ ensured in (18), we may have asymptotic convergence of the tracking error during the free space tracking phase by following Barbalat’s lemma [Reference Huang and Chien3]. On the other hand, in the compliant motion phase ( $\mathbf{F}_{\mathrm{ext}}\neq \mathbf{0}$ ), we may not ensure convergence of the tracking error by using (20). Let us assume that there exists some unknown $\alpha \gt 0$ such that $\| \mathbf{F}_{\mathrm{ext}}\| \lt \alpha$ , and we would like to prove boundedness of $\mathbf{e}=\mathbf{x}-\mathbf{x}_{d}$ . Rewrite (20) into

(21) \begin{equation} \dot{V}\leq \left(-\frac{1}{2}\lambda _{\min }\!\left(\mathbf{Q}\right)\left\| \mathbf{x}_{e}\right\| +\alpha \left\| \mathbf{P}\mathbf{B}_{e}\mathbf{M}\mathbf{D}_{x}^{-1}\right\| \right)\!\left\| \mathbf{x}_{e}\right\| \end{equation}

This implies $\dot{V}\leq 0$ if $\| \mathbf{x}_{e}\| \geq \beta$ where

\begin{equation*} \beta =\frac{2\alpha }{\lambda _{\min }\!\left(\mathbf{Q}\right)}\left\| \mathbf{P}\mathbf{B}_{e}\mathbf{M}\mathbf{D}_{x}^{-1}\right\| \end{equation*}

is a positive constant which is adjustable by proper selection of the pair $(\mathbf{P},\mathbf{Q})$ . Therefore, boundedness of the tracking error during complaint motion can be ensured.

In summary, for a robot manipulator containing parametric uncertainties, the adaptive controller (14) and update law (19) can give impedance control performance without using force sensors or force estimators. In the free space tracking phase, the robot motion will converge to the desired trajectory. In the compliant motion phase, the robot end-effector is going to touch the environment, and robot motion will deviate from the desired trajectory. However, the tracking error will still be bounded.

Figure 1. A two-link robot manipulator.

5. Computer simulations

Simulation cases are given in this section to verify the effectiveness of the proposed designs. Let us consider a two-link rigid robot as shown in Fig. 1 where $m_{1}=m_{2}=0.5(\mathrm{kg})$ are link masses, $l_{1}=l_{2}=0.75(\mathrm{m})$ are link lengths, $l_{c1}=l_{c2}=0.375(\mathrm{m})$ are distances of the center of gravity to the actuator, and $I_{1}=I_{2}=0.0234(\mathrm{kg}\cdot \mathrm{m}^{2})$ are moments of inertia of the links. The model of the robot can be derived to be in the form of (7) where [Reference Huang and Chien3, Reference Spong, Hutchinson and Vidyasagar5].

\begin{equation*} \mathbf{D}\!\left(\mathbf{q}\right)=\left[\begin{array}{c@{\quad}c} m_{1}l_{c1}^{2}+m_{2}d_{1}+I_{1}+I_{2} & m_{2}\!\left(l_{c2}^{2}+l_{1}l_{c2}\cos q_{2}\right)+I_{2}\\[5pt] m_{2}\!\left(l_{c2}^{2}+l_{1}l_{c2}\cos q_{2}\right)+I_{2} & m_{2}l_{c2}^{2}+I_{2} \end{array}\right] \end{equation*}

\begin{equation*} \mathbf{C}\!\left(\mathbf{q},\dot{\mathbf{q}}\right)=\left[\begin{array}{c@{\quad}c} -m_{2}l_{1}l_{c2}\sin q_{2}\dot{q}_{2} & -m_{2}l_{1}l_{c2}\sin q_{2}\!\left(\dot{q}_{1}+\dot{q}_{2}\right)\\[5pt] m_{2}l_{1}l_{c2}\sin q_{2}\dot{q}_{1} & 0 \end{array}\right] \end{equation*}

\begin{equation*} \mathbf{g}\!\left(\mathbf{q}\right)=\left[\begin{array}{c} \left(m_{1}l_{c1}+m_{2}l_{1}\right)\!g\cos q_{1}+m_{2}l_{c2}g\cos \!\left(q_{1}+q_{2}\right)\\[5pt] m_{2}l_{c2}g\cos \!\left(q_{1}+q_{2}\right) \end{array}\right] \end{equation*}

with $d_{1}=l_{1}^{2}+l_{c2}^{2}+2l_{1}l_{c2}\cos q_{2}$ . The motion of the robot manipulator is in the vertical plane, and a vertical smooth constraint surface is at the position $x=0.95\,\mathrm{m}$ whose stiffness is assumed to be $k=5\times 10^{4}\,\mathrm{N}/\mathrm{m}$ . The desired trajectory begins from $(0.8m,0.8m)$ which is a circle centering at $(0.8m,1m)$ with radius $0.2(m)$ . The end-effector starts from $(0.8m,0.75m)$ to track the desired trajectory in 10 s. The initial condition of the robot in the joint space is then $\mathbf{q}(0)=[\begin{array}{l@{\quad}l@{\quad}l@{\quad}l} 0.0022 & 1.5019 & 0 & 0 \end{array}]^{T}$ . There is a 0.05m initial error in the Y direction which can be used to observe the transient response of the free space tracking.

There are two simulation cases. The first case is for the performance comparisons of Hogan’s impedance controller and the proposed force sensorless impedance controller under the assumption that all system parameters are known. For case 2, the system parameters are not available and the proposed adaptive controller is employed in the simulation.

CASE 1: In this case, all system parameters are assumed to be given, and Hogan’s impedance controller (10) and the proposed force cancelation impedance controller (12) are used in the simulation for performance comparisons. The target impedance coefficient matrices are selected as $\mathbf{M}=0.15\mathbf{I}_{2}, \mathbf{B}=100\mathbf{I}_{2}$ , and $\mathbf{K}=3000\mathbf{I}_{2}$ for both the controllers. These values are selected to have a natural frequency $140\,\mathrm{rad}/\mathrm{s}$ and damping ratio 2.375 in the closed-loop dynamics so that the impact during contact can be alleviated. The matrices M, B, and K are selected to be diagonal so that the target impedance becomes a decoupled MIMO system. Selections of the diagonal elements are equivalent to the selections of the dynamics of the corresponding second-order linear subsystems. The behavior of these subsystems represents the performance of both the free space tracking and compliant motion phase. The simulation results are shown in Figs. 2–5. The trajectory tracking performance in the Cartesian space is presented in Fig. 2 where the dash-dot line is the desired trajectory, the solid line is the trajectory of the end-effector based on the proposed controller, and the dashed line is the trajectory of Hogan’s controller. The two control results are almost identical in the free space tracking phase. When the end-effector touches the constraint surface, there are no obvious impact activities for the two controllers. During the compliant motion phase, both trajectories slide on the surface smoothly. When leaving the surface, the two controllers drive the robot manipulator effectively to catch up with the desired trajectories nicely. Figure 3 shows the tracking performance in the joint space. It is seen that both the transients end in about 0.1 s, and the tracking errors are small afterward in the free space. The deviations of the trajectories from the desired path are in the compliant motion phases. Figure 4 gives the control efforts of the two controllers. We see that they are identical during the free space tracking phase, while the control effort of the new controller during the compliant motion is slightly smaller. Figure 5 is the time histories of the contact forces. Since the constraint surface is smooth, it can only exhibit a normal force. The simulation setup is specifically designed to make the contact force on the second joint be always zero. The peak value of the contact force for the traditional impedance control is about 140(N), and it is about 160(N) for the new controller.

Figure 2. Tracking performance in Cartesian space in Case 1. The dash-dot line is the desired trajectory, the solid line is the trajectory of the end-effector under the control of the proposed controller, and the dashed line is the trajectory of the end-effector under the control of Hogan’s impedance control.

Figure 3. Tracking performance in joint space forces in Case 1.

Figure 4. Time histories of control efforts in Case 1.

Figure 5. Time histories of contact forces in Case 1.

Figure 6. Tracking performance in Cartesian space in Case 2.

Figure 7. Tracking performance in joint space in Case 2.

Figure 8. Time history of control efforts in Case 2.

Figure 9. Time history of contact force in Case 2.

Figure 10. Function approximation results.

The simulation results show that the proposed controller gives almost identical motion trajectory as the Hogan’s design but induces slightly different contact force under the assumption that all system parameters are available. Therefore, the force cancelation technique is feasible to avoid the installation of force sensors in the impedance control paradigm.

CASE 2: The previous case verifies that the proposed force cancelation controller can give similar response as Hogan’s controller under the assumption that the system parameters are available even though the force sensor information is not used. Here, in this case, we would like to assume that the parameters in $(\mathbf{D}_{x},\mathbf{C}_{x},\mathbf{g}_{x})$ are not available and the proposed adaptive controller (14) with the update law (19) is implemented. The target impedance coefficient matrices are the same as those in the previous case, and the controller parameters are selected to be $\boldsymbol{\Gamma }=0.1^{\ast}\mathbf{I}_{21\times 21}$ and

\begin{equation*} \mathbf{Q}=\left[\begin{array}{c@{\quad}c@{\quad}c@{\quad}c} 20000 & 0 & 3000 & 0\\[5pt] 0 & 140000 & 0 & 1200\\[5pt] 3000 & 0 & 600 & 0\\[5pt] 0 & 1200 & 0 & 300 \end{array}\right]. \end{equation*}

A 21-term Fourier series is utilized as the basis functions, and the initial estimate is selected to be $\hat{\mathbf{h}}(0)=[\begin{array}{l@{\quad}l} -3 & 5 \end{array}]$ . The initial weighting vectors for function approximation are $\hat{\mathbf{w}}_{1}(0)=[\!-\!30\cdots 0]^{T}\in \Re ^{21x1}$ and $\hat{\mathbf{w}}_{2}(0)=[50\cdots 0]^{T}\in \Re ^{21x1}$ . The simulation results are shown in Figs. 6–10. The tracking performance in the Cartesian space given in Fig. 6 is like the one in case 1. The transient trajectory is smooth even though the system parameters are unavailable. Figure 7 shows that the tracking performance in the joint space is nice where the transient is within 0.1 s. Figure 8 gives the control efforts to the joints whose magnitudes are realizable. Figure 9 is the time histories of the contact forces whose peak value is the same as the one in Case 1. Figure 10 is the parameter estimation time history where all values are bounded as proved in Section 4. The discrepancies in the function approximation do not invalidate the convergence of the free space tracking as seen in Figs. 6 and 7. In summary, the adaptive control version of the proposed force sensorless impedance controller can give almost the same performance as the new controller for the known system parameters case.