2.1.1. Inverse kinematics

In this study, simulation and experimentation are conducted on a 6-DoF PM. It consists of a fixed base (lower) and a moving (upper or end-effector) platform (frames) connected by six extendable legs. The legs are connected to the fixed base via spherical or universal joints and to the moving platform via universal joints. 6-DoF robot involves linear position and orientations as represented in Figure 1. Linear motions encompass movements along the longitudinal (surge), lateral (sway), and vertical (heave) axes. Angular movements, on the other hand, are described by rotations about x- (roll), y- (pitch), and z-axes (yaw), respectively. The inverse kinematics problem of PM robot focuses on determining the required leg lengths corresponding to a desired position and orientation of the moving platform, which is essential for control of moving platform and ensuring the robot’s ability to achieve complex tasks with required accuracy.

Figure 1.

Kinematic configurations of the 6-DoF parallel manipulator.

To derive the equations of inverse kinematics, we establish two reference frames, $B$ and $P$ , fixed to the base and moving platform, respectively, with origins $O_B$ and $O_P$ . The connection points of the lower and upper plate joints are designated as $B_i$ and $P_i$ , respectively, where $i=1,\ldots, 6$ . The angles between $B_1$ to $B_2$ and $B_3$ to $B_5$ are set to $120^\circ$ to symmetrically distribute connection points on the base plate, a symmetry also maintained for the upper platform. The angles between adjacent joints on the upper platform are denoted as $\theta _P$ , while those on the lower platform are denoted as $\theta _B$ as shown in Figure 1. The position vectors of the centers of spherical joints in frame $P$ are denoted by $ \mathbf {P}_i$ and can be expressed as follows:

(1) \begin{align} \mathbf {P_i} = \begin{bmatrix} r_P \cos\!(\lambda _{P,i}) \\[3pt] r_P \sin\!(\lambda _{P,i}) \\[3pt] P_{iz} \end{bmatrix} = \begin{bmatrix} P_{ix} \\[3pt] P_{iy} \\[3pt] P_{iz} \end{bmatrix} \end{align}

where $\lambda _P$ is the upper frame joints’ position angles relative to the moving frame $O_P$ axis and $r_P$ is the radius of the moving platform passing through the connecting points.

Similarly, the position vectors of the centers of spherical joints in frame $B$ are denoted by $ \mathbf {B}_i$ and can be expressed as follows:

(2) \begin{align} \mathbf {B_i} = \begin{bmatrix} r_B \cos\!(\lambda _{B,i}) \\[3pt] r_B \sin\!(\lambda _{B,i}) \\[3pt] B_{iz} \end{bmatrix} = \begin{bmatrix} B_{ix} \\[3pt] B_{iy} \\[3pt] B_{iz} \end{bmatrix} \end{align}

where $\lambda _B$ is the lower frame joints’ position angles relative to the fixed frame $O_B$ axis and $r_B$ is the radius of the fixed platform passing through the connecting points.

The 6-DoF of the upper frame including three translations and three rotations are defined within the following vector:

(3) \begin{align} \mathbf {X} = \begin{bmatrix} x & \quad y & \quad z & \quad \alpha & \quad \beta & \gamma \end{bmatrix}^T \end{align}

In this context, $ x$ , $ y$ , and $ z$ denote surge, sway, and heave translations, respectively, and $ \alpha$ , $ \beta$ , and $ \gamma$ represent roll, pitch, and yaw rotations of the moving platform relative to the fixed frame $O_B$ . The inverse kinematics of a PM typically provide a unique solution. Specifically, when the position and orientation of the upper platform are specified, the lengths of the legs corresponding to these configurations can be precisely determined using the inverse kinematics model. This method employs closed-loop vectors based on geometric relationships and trigonometry within the system. Subsequently, the leg vector pointing from attachment point $ B_i$ to attachment point $ P_i$ can be expressed as follows:

(4) \begin{align} \mathbf {l}_{i}^B = \mathbf t_{p}^B + \mathbf R_{\alpha \beta \gamma }.\mathbf p_{i}^P - \mathbf b_{i}^B = \mathbf t_{p}^B + \underbrace {\mathbf R_{z}(\gamma ) . \mathbf R_{y}(\beta ) . \mathbf R_{x}(\alpha )}_{Rotation\,Matrices} . \mathbf p_{i}^P - \mathbf b_{i}^B \end{align}

where $\mathbf t_{p}^B = \begin{bmatrix} p_x & p_y & p_z \end{bmatrix}^T$ is the moving platform translation vector, with respect to the base origin $O_B$ , and $\mathbf p_{i}^P = \begin{bmatrix} P_{ix} & P_{iy} & P_{iz} \end{bmatrix}^T$ and $\mathbf b_{i}^B = \begin{bmatrix} B_{ix} & B_{iy} & B_{iz} \end{bmatrix}^T$ are the translational offset vectors of the $i^{th}$ leg on the moving and base platforms, respectively. The matrix $R_{\alpha \beta \gamma }$ representing the orientation of the frame $P$ with respect to frame $B$ is expressed as follows:

(5) \begin{align} \mathbf R_{\alpha \beta \gamma } = \begin{bmatrix} \mathbf R_{z}(\gamma ) . \mathbf R_{y}(\beta ) . \mathbf R_{\alpha }\end{bmatrix} = \left[\begin{array}{c@{\quad}c@{\quad}c} c\alpha c\beta & -s\alpha c\gamma + c\alpha s\beta s\gamma & s\alpha s\gamma + c\alpha s\beta c\gamma \\[3pt] s\alpha c\beta & c\alpha c\gamma + s\alpha s\beta s\gamma & -c\alpha s\gamma + s\alpha s\beta c\gamma \\[3pt] -s\beta & c\beta s\gamma & c\beta c\gamma \end{array}\right] = \left[\begin{array}{c@{\quad}c@{\quad}c} r_{11} & r_{12} & r_{13} \\[3pt] r_{21} & r_{22} & r_{23} \\[3pt] r_{31} & r_{32} & r_{33} \end{array}\right] \end{align}

where $c$ and $s$ denote, respectively $cos(.)$ and $sin(.)$ . From Eq. (4), the length of the $i^{th}$ leg can be calculated by $||l_{i}||$ , which solves the six leg lengths of the PM robot. Hence the complete solution for the inverse kinematics can be expressed as follows:

(6) \begin{align} ||l_{i}|| = \sqrt {\begin{bmatrix} (t_{p}^B + R_{pqr}.p_{i}^P - b_{i}^B)^T & (t_{p}^B + R_{\alpha \beta \gamma }.p_{i}^P - b_{i}^B) \end{bmatrix}} \end{align}

for i = 1, 2, …, 6

2.1.2. Kinematic Jacobian

For parallel robots, the direct kinematic problem (finding the end-effector position given the joint positions) is often difficult to solve due to their closed-loop structure. Therefore, the inverse kinematic problem (finding the joint positions given the end-effector position) is typically solved first. Once the inverse kinematics are known, the inverse kinematic Jacobian can be derived. It’s used to connect joint velocities and forces to the end-effector velocities, torque, and forces. Moreover, the characteristics of the Jacobian are frequently employed in examining kinematic singularities. The leg positions being functions of end-effector positions are defined as follows:

(7) \begin{align} l_i = f_i(x,y,z,\alpha, \beta, \gamma ) = f_i(\mathbf X) \quad \text {for} \quad i = 1, 2, \ldots, 6 \end{align}

The leg velocities, in terms of the upper platform velocities, are obtained by applying the chain rule to the partial derivatives of the functions $f_i$ . This results in the derivatives of $l_i$ being expressed as functions of the derivatives of $x$ , $y$ , z, $\alpha$ , $\beta$ , $\gamma$ .

(8) \begin{align} \delta l_i = \frac {\partial f_i}{\partial x} \delta x + \frac {\partial f_i}{\partial y} \delta y + \frac {\partial f_i}{\partial z} \delta z + \frac {\partial f_i}{\partial \alpha } \delta \alpha + \frac {\partial f_i}{\partial \beta } \delta \beta + \frac {\partial f_i}{\partial \gamma } \delta \gamma \end{align}

Dividing both sides by the differential time element $ \delta t$ and arranging in matrix form yields:

(9) \begin{align} \begin{bmatrix} \dot {l}_1 \\[3pt] \dot {l}_2 \\[3pt] \vdots \\[3pt] \dot {l}_6 \end{bmatrix} = \left[\begin{array}{l@{\quad}l@{\quad}l@{\quad}l@{\quad}l@{\quad}l} \dfrac {\partial f_1}{\partial x} & \dfrac {\partial f_1}{\partial y} & \dfrac {\partial f_1}{\partial z} & \dfrac {\partial f_1}{\partial \alpha } & \dfrac {\partial f_1}{\partial \beta } & \dfrac {\partial f_1}{\partial \gamma } \\[6pt] \dfrac {\partial f_2}{\partial x} & \frac {\partial f_2}{\partial y} & \dfrac {\partial f_2}{\partial z} & \dfrac {\partial f_2}{\partial \alpha } & \dfrac {\partial f_2}{\partial \beta } & \dfrac {\partial f_2}{\partial \gamma } \\[6pt] \vdots & \vdots & \vdots & \vdots & \vdots & \vdots \\[6pt] \dfrac {\partial f_6}{\partial x} & \dfrac {\partial f_6}{\partial y} & \dfrac {\partial f_6}{\partial z} & \dfrac {\partial f_6}{\partial \alpha } & \dfrac {\partial f_6}{\partial \beta } & \dfrac {\partial f_6}{\partial \gamma } \end{array}\right] \begin{bmatrix} \dot {x} \\[6pt]\dot {y} \\[6pt] \dot {z} \\[6pt]\dot {\alpha } \\[6pt] \dot {\beta } \\[6pt] \dot {\gamma } \end{bmatrix} \end{align}

Let $ J^{-1} \in \mathbb {R}^{6 \times 6}$ denote the inverse Jacobian (kinematic Jacobian) matrix composed of these partial derivatives. The forward and inverse kinematics are then expressed as follows:

(10) \begin{align} \dot {\mathbf {X}} = J \dot {\mathbf {l}} \\[-24pt] \nonumber \end{align}

(11) \begin{align} \dot {\mathbf {l}} = J^{-1} \dot {\mathbf {X}} \end{align}

Here, $ \mathbf {X} \in \mathbb {R}^6$ represents the task-space end-effector positions, while $ \dot {\mathbf {X}}$ and $ \ddot {\mathbf {X}}$ are the corresponding velocities and accelerations. $ J \in \mathbb {R}^{6 \times 6}$ is the Jacobian matrix.

2.1.3. Forward kinematics

The forward kinematics determines the position and orientation of the upper platform based on the provided leg lengths. The forward kinematics of a PM involve a set of complex nonlinear differential equations, making them challenging to solve. As a result, iterative solvers such as the Newton–Raphson method are commonly used to obtain approximate solutions. The general form of the function used for forward kinematics can be given as follows:

(12) \begin{align} \begin{split} F_i(x,y,z,\alpha, \beta, \gamma ) = & x^2 + y^2 + (z+z_0)^2 + r_p^2 + r_b^2 \\[3pt] & + 2(r_{11} P_{ix} + r_{12} P_{iy})(p_x - B_{ix}) \\[3pt] & + 2(r_{21} P_{ix} + r_{22} P_{iy})(p_y - B_{iy}) \\[3pt] & + 2(r_{31} P_{ix} + r_{32} P_{iy})(p_z) \\[3pt] & - 2(P_x B_{ix} + P_y B_{iy})-||l_{i}||^2 \quad \text {for} \quad i = 1, 2, \ldots, 6 \end{split} \end{align}

3.1.1. Description of the controller

AbC or acceleration control (ACC) is a control strategy that directly employs joint accelerations to command robotic manipulators. By utilizing feedforward acceleration reference, AbC can potentially achieve faster behavior and more accurate tracking performance compared to traditional control methods. However, disturbances such as external forces, friction, or modeling errors can affect the system’s performance. To mitigate the effects of disturbances, a DOB is integrated into the AbC architecture as described in ref. [Reference Sariyildiz, Oboe and Ohnishi40]. While the AbC with DOB offers improved tracking performance and disturbance rejection, it may still struggle with highly nonlinear systems and random disturbances. In this study, a sliding mode controller is introduced alongside the AbC with DOB, constituting a hybrid controller. This hybridization allows the system to leverage the strengths of each individual method while compensating for their respective weaknesses. Let us define the joint position errors $\mathbf {e} \in \mathbb {R}^6$ of the joints as follows:

(17) \begin{align} \mathbf e(t) = \mathbf q_r(t) - \mathbf q(t) \end{align}

where $\mathbf q_r \in \mathbb {R}^6$ is the vector of reference actuators (joints) trajectories. Let us define the sliding surface as follows:

(18) \begin{align} \mathbf {s}(t) = \dot {\mathbf {q}}_r(t) - \dot {\mathbf {q}}(t) + \lambda (\mathbf {q}_r(t) - \mathbf {q}(t)) = \dot {\mathbf {e}}(t) + \lambda \mathbf {e}(t) \end{align}

where $ \mathbf {s}(t) \in \mathbb {R}^6$ represents the sliding surface, $\dot {\mathbf {e}}(t) \in \mathbb {R}^6$ is the joint tracking rate error, and $\dot {\mathbf {q}}_r(t) \in \mathbb {R}^6$ is the vector of desired velocities. $\lambda = \text {diag}(\lambda )$ with $\lambda \gt 0$ for $(i = 1,\ldots, 6)$ is the diagonal positive definite constant matrix that determines the slope of the sliding surface.

The time derivative of the sliding surface is then given as follows:

(19) \begin{align} \dot {\mathbf {s}}(t) = \ddot {\mathbf {q}}_r(t) - \ddot {\mathbf {q}}(t) + \lambda (\dot {\mathbf {q}}_r(t) - \dot {\mathbf {q}}(t)) = \ddot {\mathbf {e}}(t) + \lambda \dot {\mathbf {e}}(t) \end{align}

Based on Eq. (15), the expression of the acceleration vector is written as follows:

(20) \begin{align} \ddot {\mathbf {q}}(t) = {M}_n^{-1}(\mathbf {q}) (\tau (t) - \boldsymbol {\tau _{\text {dis}}}(t, \mathbf {q}) - V_n(\mathbf {q}, \dot {\mathbf {q}})\dot {\mathbf {q}} - G_n(\mathbf {q})) \end{align}

Inserting Eq. (19) into Eq. (18), one gets

(21) \begin{align} \dot {\mathbf {s}}(t) = \ddot {\mathbf {q}}_r(t) - {M}_n^{-1}(\mathbf {q}) \left (\tau (t) - \boldsymbol {\tau _{\text {dis}}}(t, \mathbf q) - V_n(\mathbf {q}, \dot {\mathbf {q}})\dot {\mathbf {q}} - G_n(\mathbf {q})\right ) + \lambda \dot {\mathbf {e}}(t) \end{align}

The system dynamics are designed to reach the sliding surface in a finite time. A common reaching law could be as follows:

(22) \begin{align} \dot {\mathbf s}(t) = -k_0 \text {sgn} (\mathbf s(t)) - \eta \mathbf s(t) \end{align}

where $ k_0$ and $ \eta$ are positive constants and $ \text {sgn} (\mathbf s(t))$ denotes the signum function. Based on Eq. (21) and Eq. (22), the input force or control law can be obtained as follows:

(23) \begin{align} \boldsymbol {\tau }(t) = M_n(\mathbf {q}) (\ddot {\mathbf {q}}_r(t) + k_0 \text {sgn}(\mathbf {s}(t)) + \eta \mathbf {s}(t) + \lambda \dot {\mathbf {e}}(t)) + \boldsymbol {\tau }_{\text {dis}}(t, \mathbf {q}) + V_n(\mathbf {q}, \dot {\mathbf {q}})\dot {\mathbf {q}} + G_n(\mathbf {q})) \end{align}

Here $\boldsymbol {\tau }(t)$ is control law, but $\boldsymbol \tau _{\text {dis}}(t, q)$ contains modeling mismatch, friction, and other disturbances. In order to address these generalized disturbances, an appropriate estimation law $\boldsymbol {\tau _{\text {dob}}} (t, q)$ , which is the estimated output of $\boldsymbol {\tau _{\text {dis}}}(t, q)$ , will be incorporated into the control law. Based on [Reference Sariyildiz, Oboe and Ohnishi40–Reference Mohammadi, Tavakoli, Marquez and Hashemzadeh42], the estimation of the disturbance vector can be derived in terms of an auxiliary variable and the states of the system by utilizing the known nonlinear dynamic model of PMs as follows:

(24) \begin{align} \boldsymbol {\tau _{\text {dob}}} (t, q) = \mathbf {z} - L \dot {\mathbf {q}} \end{align}

where $\boldsymbol {\tau _{\text {dob}}} (t, q)$ represents the estimation of $\boldsymbol \tau _{\text {dis}}(t, q)$ , $L \in \mathbb {R}$ represents the observer gain of DOB to be tuned, and $\mathbf {z} \in \mathbb {R}$ represents the auxiliary variable vector, which can be obtained by integrating the following equation:

(25) \begin{align} \dot {\mathbf {z}} = L{M}_n^{-1}(\mathbf {q})(\boldsymbol {\tau }(t) - \boldsymbol {\tau }_{\text {dob}}(t, \mathbf {q})) - L{M}_n^{-1}(\mathbf {q})(V_n(\mathbf {q}, \dot {\mathbf {q}})\dot {\mathbf {q}} - G_n(\mathbf {q})) \end{align}

The dynamic equation of disturbance estimation can be obtained by taking the derivative of Eq. (24) and inserting it into Eq. (25) and Eq. (15) as follows:

(26) \begin{align} \dot {\boldsymbol {\tau }}(t, \mathbf {q})_{\text {dob}} &= L{M}_n^{-1}(\mathbf {q})(\boldsymbol {\tau }(t) - \boldsymbol {\tau }_{\text {dob}}(t, \mathbf {q}) - V_n(\mathbf {q}, \dot {\mathbf {q}})\dot {\mathbf {q}} - G_n(\mathbf {q})) \nonumber \\[3pt] &\quad - L{M}_n^{-1}(\mathbf {q})(\boldsymbol {\tau }(t) - \boldsymbol {\tau }_{\text {dis}}(t, \mathbf {q}) - V_n(\mathbf {q}, \dot {\mathbf {q}})\dot {\mathbf {q}} - G_n(\mathbf {q})) \nonumber \\[3pt] &\quad \dot {\boldsymbol {\tau }}(t, \mathbf {q})_{\text {dob}} = L{M}_n^{-1}(\mathbf {q})(\boldsymbol {\tau }_{\text {dis}}(t, \mathbf {q})-\boldsymbol {\tau }_{\text {dob}}(t, \mathbf {q})) \end{align}

If $\boldsymbol {{\dot \tau }}_{dob}$ is subtracted from both side, then Eq. (27) can be formed as follows:

(27) \begin{align} \dot {\mathbf {e}}(t)_{\text {dis}} + L{M}_n^{-1}(\mathbf {q})\mathbf {e}_{\text {dis}} = \boldsymbol {{\dot \tau }}(t, \mathbf {q})_{\text {dis}} \end{align}

where $\mathbf e_{\text {dis}} = \boldsymbol {\tau _{\text {dob}}} (t, \mathbf q) - \boldsymbol {\tau _{\text {dis}}} (t, \mathbf q) \in \mathbb {R}^6$ .

Assumption 2. The change rate of the lump disturbance $\boldsymbol {{\dot \tau }}_{dis} (t, \mathbf q)$ is bounded by $\|\boldsymbol {{\dot \tau }}_{dis} (t, \mathbf q)|_{t\gt 0} \leq \delta _2$ , where $\delta _2 \in \mathbb {R}^+$ is a constant.

With Assumptions 1–2 and $L$ strictly positive definite, then it can be shown that the error of disturbance estimation is uniformly ultimately bounded, that is,

(28) \begin{align} \| e_{\text {dis}}(t) \| \leq \exp ({-}\lambda _u (t - t_0)) \| e_{\text {dis}}(t_0) \| + \frac {\delta _2}{\lambda _u} \end{align}

where $\lambda _u = L M_n^{-1}(q) \gt 0$ .

Equation (28) shows that the convergence rate and the accuracy of disturbance estimation are directly adjusted by tuning the nominal inertia matrix and the gain of the DOB.

Replacing the designed DOB $\boldsymbol {\tau _{\text {dob}}} (t, \mathbf q)$ with $\boldsymbol {\tau _{\text {dis}}} (t, \mathbf q)$ in Eq. (23), the control rule becomes as follows:

(29) \begin{align} \boldsymbol {\tau }(t) = {M}_n(\mathbf {q}) \left (\ddot {\mathbf {q}}_r(t) + k_0 \operatorname {sgn}(\mathbf s(t)) + \eta \mathbf s(t) + \lambda \dot {\mathbf {e}}(t) \right ) +{\boldsymbol {\tau }}_{\text {dob}}(t, \mathbf {q}) + V_n(\mathbf {q}, \dot {\mathbf {q}})\dot {\mathbf {q}} + G_n(\mathbf {q})) \end{align}

where the term $\ddot {\mathbf {q}}_r(t) + \eta \mathbf s(t) + \lambda \dot {\mathbf {e}}(t)$ is called the reference feedforward acceleration, and it can be represented as follows:

(30) \begin{align} \ddot {\mathbf {q}}_{d}(t) = \ddot {\mathbf {q}}_r(t) (\eta + \lambda ) (\dot {\mathbf {e}}(t) + (\eta \lambda ) \mathbf {e}(t)) = \ddot {\mathbf {q}}_r(t) + k_v \dot {\mathbf {e}}(t) + k_p \mathbf {e}(t) \end{align}

where $k_v$ and $k_p$ $\in \mathbb {R}^{6x6}$ are the diagonal velocity and position control gain matrices. Then the controller rule can be obtained as follows:

(31) \begin{align} \mathbf {u}(t) = {M}_n(\mathbf {q})\left ( \ddot {\mathbf {q}}_d(t) + k_0 \operatorname {sgn}(\mathbf s(t)) \right )+{\boldsymbol {\tau }}_{\text {dob}}(t, \mathbf {q})+V_n(\mathbf {q}, \dot {\mathbf {q}})\dot {\mathbf {q}} + G_n(\mathbf {q})) \end{align}

As can be seen, it is actually a form of acceleration feedforward control with integration of sliding reaching term and DOB. Hence, it represents a hybrid control approach. The controller parameters $k_v$ and $k_p$ are tuned by considering the nominal system dynamics only, since the model mismatch and disturbances can be compensated by sliding face and DOB.

3.1.2. Stability analysis

In order to prove the stability of the proposed controller, let us refer to the following stability theorem:

Stability Theorem: When dealing with a nonlinear uncertain dynamical system as defined by Eq. (15), the application of the robust control law outlined by Eq. (29) ensures asymptotic robust stability of the closed-loop system in the presence of model uncertainties and disturbances.

Proof: Let us consider the following Lyapunov function candidate and its time derivative:

(32) \begin{align} V = \frac {1}{2}\mathbf {s}(t)^T \mathbf {s}(t)+\frac {1}{2}\mathbf {e}_{\text {dis}}^T \mathbf {e}_{\text {dis}} \\[-24pt] \nonumber \end{align}

(33) \begin{align} \dot {V} = \mathbf {s} \dot {\mathbf {s}} + \mathbf {e}_{\text {dis}} \dot {\mathbf {e}}_{\text {dis}} \end{align}

Inserting control law of Eq. (29) into Eq. (18), one gets sliding face as follows:

(34) \begin{align} \dot {\mathbf {s}}(t) = -k_0 \operatorname {sgn}\mathbf {s}(t) - \eta \mathbf {s}(t)-M_n^{-1}\mathbf {e}_{\text {dis}} \end{align}

Using Eq. (30) of DOB error dynamic with Assumption2 of the change rate of the lump disturbance $\boldsymbol {{\dot \tau }}_{dis} (t, \mathbf q)$ is bounded by $\|\|\boldsymbol {{\dot \tau }}_{dis} (t, \mathbf q)\| \leq \delta _2$ , and taking $\delta _2 = 0$ , and inserting into Eq. (34) into Eq. (33), Eq. (33) yields as follows:

(35) \begin{align} \dot {V} = \mathbf {s}(t) ({-}k_0 \operatorname {sgn}\mathbf {s}(t) - \eta \mathbf {s}(t)-{M}_n^{-1}\mathbf {e}_{\text {dis}})- L{M}_n^{-1} ||\mathbf {e}_{\text {dis}}||^{2} \end{align}

Expending Eq. (35)

(36) \begin{align} \dot {V} = -k_0 ||\mathbf {s}(t)|| - \eta ||\mathbf {s}(t)||^{2}-{M}_n^{-1}\mathbf {s}(t)\mathbf {e}_{\text {dis}}- L{M}_n^{-1} ||\mathbf {e}_{\text {dis}}||^{2} \end{align}

where $k_0 \geq M_n^{-1} ||\mathbf {e}_{\text {dis}}||_{max}$ .

Consider $||\mathbf {e}_{\text {dis}}(0)|| = ||\mathbf {e}_{\text {dis}}||_{max}$ , we have $k_0 \geq M_n^{-1} ||\mathbf {e}_{\text {dis}}(0)||$ , then

(37) \begin{align} \dot {V} \leq - \eta ||\mathbf {s}(t)||^{2}-L M_n^{-1} ||\mathbf {e}_{\text {dis}}||^{2} \leq -k_1 \left (\frac {1}{2}\mathbf {s}(t)^T \mathbf {s}(t)+\frac {1}{2}\mathbf {e}_{\text {dis}}^T \mathbf {e}_{\text {dis}}\right ) = -k_1 V \end{align}

where $k_1 = 2\min \{ \eta, LM_n^{-1} \}$ .

The following lemma can be used to the solution of the inequality of Eq. (37):

Lemma. [Reference Liu36]–[Reference Ioannou and Sun37 ]

Let $f, V: [0, \infty ) \in \mathbb {R}$ , then $\dot {V} \leq -\alpha V + f$ for $\forall t \geq t_0 \geq 0$ implies that

(38) \begin{align} V(t) \leq e^{2\alpha (t - t_0)} V(t_0) + \int _{t_0}^{t} e^{2\alpha (t - \tau )} f(\tau ) d\tau \end{align}

for any finite constant $ \alpha$ .

According to [Reference Atassi and Khalil39], we have the proof as follows:

Let $ \omega (t) \triangleq \dot {V} - \alpha V - f$ , we have $ \omega (t) \leq 0$ , and

(39) \begin{align} \dot {V} = -\alpha V + f + \omega \end{align}

implies that

(40) \begin{align} V(t) = e^{-\alpha (t - t_0)} V(t_0) + \int _{t_0}^{t} e^{-\alpha (t - \tau )} f(\tau ) d\tau + \int _{t_0}^{t} e^{-\alpha (t - \tau )} \omega (\tau ) d\tau \end{align}

Since $ \omega (t) \lt 0$ and $\forall t \geq t_0 \geq 0$ , one has the following inequality:

(41) \begin{align} V(t) \leq e^{-\alpha (t - t_0)} V(t_0) + \int _{t_0}^{t} e^{-\alpha (t - \tau )} f(\tau ) d\tau . \end{align}

Moreover, setting $ f = 0$ , one has $ \dot {V} \leq -\alpha V$ , implying the following inequality:

(42) \begin{align} V(t) \leq e^{-\alpha (t - t_0)} V(t_0) \end{align}

If $ \alpha$ has a positive constant value, then $ V(t)$ will tend to zero exponentially.

Using above the lemma, Eq. (37) can be obtained as follows:

(43) \begin{align} \dot {V} \leq -k_1 V \\[-24pt] \nonumber \end{align}

(44) \begin{align} V(t) \leq e^{-k_1(t-t_0)} V(t_0) \end{align}

Notice that the right side of inequality decreases linearly as time $t$ increases, and so the function will converge to zero exponentially, showing the stability of the proposed controller.

It is worth noting that the sliding mode part of the control rule leads to oscillations in the controlled system due to the discontinuous nature of the sign function used in the control law. This oscillatory behavior, known as chattering, is also present in the proposed controller described in Eq. (29). To mitigate chattering, a saturation function can be employed, leading to a modification of the control law as follows:

(45) \begin{align} \text {sat}(\dot {\mathbf {e}} + \lambda \mathbf {e}, \phi ) = \begin{cases} \text {sign}(\dot {\mathbf {e}} + \lambda \mathbf {e}) & \text {if } |\dot {\mathbf {e}} + \lambda \mathbf {e}| \leq \phi \cr \dfrac {\dot {\mathbf {e}} + \lambda \mathbf {e}}{\phi } & \text {if } |\dot {\mathbf {e}} + \lambda \mathbf {e}| \gt \phi \end{cases} \end{align}

The final tracking error $e$ can be constrained within a specified precision $\varepsilon$ , known as the boundary layer width as specified in ref. [Reference Slotine and Li20]. The tracking error $e$ can be calculated as follows:

(46) \begin{align} \varepsilon = \frac {\phi }{\lambda } \end{align}

Equation (46) indicates that the final maximum tracking error can be regulated by appropriately selecting the boundary layer $\phi$ and the slope constant $\lambda$ . Then, the final form of the hybrid joint-space control rule to be tuned is as follows:

(47) \begin{align} \mathbf {u}(t) = \mathbf {M}_n(\mathbf {q})\left ( \ddot {\mathbf {q}}_d(t) + k_0 \text {sat}(\dot {\mathbf {e}} + \lambda \mathbf {e}, \phi ) \right )+{\boldsymbol {\tau }}_{\text {dob}}(t, \mathbf {q}) \end{align}

3.2.1. Stability analysis of combined control method

Let us consider the following Lyapunov function candidate $V_2$ and its time derivative:

(55) \begin{align} V_2 = \frac {1}{2}\mathbf {s}(t)^T \mathbf {s}(t)+\frac {1}{2}\mathbf {e}_{\text {dis}}^T \mathbf {e}_{\text {dis}} \end{align}

(56) \begin{align} \dot V_2 = \mathbf {s} \dot {\mathbf {s}} + \mathbf {e}_{\text {dis}} \dot {\mathbf {e}}_{\text {dis}} \end{align}

Substituting the combined control law inside the dynamic model in joint space, one gets the following closed-loop acceleration dynamic model:

(57) \begin{align} \begin{aligned} \ddot {\mathbf {q}}(t) = &\ \ddot {\mathbf {q}}_r(t) + k_0 \operatorname {sgn}(\mathbf {s}(t)) + \eta \mathbf {s}(t) + \lambda \dot {\mathbf {e}}(t) \\[3pt] & + M_n(\mathbf {q})^{-1} \left [\boldsymbol {\tau }_{\text {dob}}(t, \mathbf {q}) - \boldsymbol {\tau }_{\text {dis}}(t, \mathbf {q}) \right ] \\[3pt] & + M_n(\mathbf {q})^{-1} J^T \mathbf {r}(t)) \end{aligned} \end{align}

The time derivative of the sliding surface is given as follows:

(58) \begin{align} \dot {\mathbf {s}}(t) = \ddot {\mathbf {q}}_r(t) - \ddot {\mathbf {q}}(t) + \lambda (\dot {\mathbf {q}}_r(t) - \dot {\mathbf {q}}(t)) \end{align}

Substituting Eq. (57) into (58) yields as follows:

(59) \begin{align} \begin{aligned} \dot {\mathbf {s}}(t) = &\ - k_0 \operatorname {sgn}(\mathbf {s}(t)) - \eta \mathbf {s}(t) - M_n(\mathbf {q})^{-1} \mathbf {e}_{\text {dis}}(t) - M_n(\mathbf {q})^{-1} J^T \mathbf {r}(t)) \end{aligned} \end{align}

Now inserting Eq. (59) into Eq. (56) results in the following expression of $\dot V_2$ :

(60) \begin{align} \dot V_2 = \mathbf {s}(t) ({-}k_0 \operatorname {sgn}\mathbf {s}(t) - \eta \mathbf {s}(t)-{M}_n^{-1}\mathbf {e}_{\text {dis}}- M_n(\mathbf {q})^{-1} J^T \mathbf {r}(t)))- L{M}_n^{-1} ||\mathbf {e}_{\text {dis}}||^{2} \end{align}

Considering $||\mathbf {e}_{\text {dis}}(0)|| = ||\mathbf {e}_{\text {dis}}||_{max}$ , and $k_0 \geq M_n^{-1} ||\mathbf {e}_{\text {dis}}(0)||$ , then

(61) \begin{align} \dot V_7 \leq - \eta ||\mathbf {s}(t)||^{2}-L M_n^{-1} ||\mathbf {e}_{\text {dis}}||^{2}-\mathbf {s}(t) M_n(\mathbf {q})^{-1} J^T \mathbf {r}(t) \end{align}

With $k =1$ and $\lambda = c$ , the sliding mode in joint-space $\mathbf {s}(t)$ and sliding mode in task-space $\mathbf {r}(t)$ have the following relationship:

(62) \begin{align} \mathbf {s}(t) = J^T \mathbf {r}(t) \end{align}

Hence, the last terms in $\dot V_2$ equation can be written as follows:

(63) \begin{align} -\mathbf {s}(t) M_n(\mathbf {q})^{-1} J^T \mathbf {r}(t) = -\mathbf {s}(t) M_n(\mathbf {q})^{-1} \mathbf {s}(t) \end{align}

This expression now represents a quadratic form involving the vector $\mathbf {s}(t)$ and the inverse of the inertia matrix $M_n(\mathbf {q})^{-1}$ . When the intention is to interpret the result as a weighted norm, we typically write

(64) \begin{align} -\mathbf {s}(t)^T M_n(\mathbf {q})^{-1} \mathbf {s}(t) = -M_n(\mathbf {q})^{-1} ||\mathbf {s}(t)||^2 \end{align}

Hence, inserting Eq. (64) into Eq. (60) yields the following:

(65) \begin{align} \dot V_2 \leq - \eta ||\mathbf {s}(t)||^{2}-L M_n^{-1} ||\mathbf {e}_{\text {dis}}||^{2}-M_n(\mathbf {q})^{-1} ||\mathbf {s}(t)||^2 \end{align}

Factoring out for the same terms:

(66) \begin{align} \dot V_2 \leq - \left (\eta + M_n^{-1}\right ) ||\mathbf {s}(t)||^2 - L M_n^{-1} ||\mathbf {e}_{\text {dis}}||^{2} \end{align}

(67) \begin{align} \dot V_2 \leq -k_3 \left (\frac {1}{2}\mathbf {s}(t)^T \mathbf {s}(t)+\frac {1}{2}\mathbf {e}_{\text {dis}}^T \mathbf {e}_{\text {dis}} \right ) \end{align}

(68) \begin{align} \dot V_2 \leq -k_3 V_2 \end{align}

where $k_3 = 2\min \{ (\eta + M_n^{-1}), LM_n^{-1}\}$ .

Based on the same Lemma [Reference Ioannou and Sun37] again, one gets

(69) \begin{align} V_2(t) \leq e^{-k_3 (t - t_0)} V_2(t_0) \end{align}

If $k_3$ has a positive constant value, then $ V_2(t)$ will tend to zero exponentially. Specifically, the right side of inequality decreases linearly as time $t$ increases, and so the function will converge to zero exponentially, showing the stability of the proposed controller.