Modeling and PDC fuzzy control of planar parallel robot

Abstract

Many works in the literature have studied the kinematical and dynamical issues of parallel robots. But it is still difficult to extend the vast control strategies to parallel mechanisms due to the complexity of the model-based control. This complexity is mainly caused by the presence of multiple closed kinematic chains, making the system naturally described by a set of differential–algebraic equations. The aim of this work is to control a two-degree-of-freedom parallel manipulator. A mechanical model based on differential–algebraic equations is given. The goal is to use the structural characteristics of the mechanical system to reduce the complexity of the nonlinear model. Therefore, a trajectory tracking control is achieved using the Takagi-Sugeno fuzzy model derived from the differential–algebraic equation forms and its linear matrix inequality constraints formulation. Simulation results show that the proposed approach based on differential–algebraic equations and Takagi-Sugeno fuzzy modeling leads to a better robustness against the structural uncertainties.

Introduction

In modern societies, manipulators have become a key tool in industry due to their great versatility in repetitive task. The mechanical architectures traditionally used in automated production line are based on serial robots. These robots are made of a sequence of rigid links serially assembled through revolute and/or prismatic active joints. ¹ In other words, the end effector is connected to the base via a single open-loop kinematic chain which gives them a large work space and high dexterity. ² Despite these advantages, serial robot suffer from high inertia due to the position of actuators, which are located on the moving part, they also suffer from a low payload to weight ratio and a low precision due to the cumulative joint errors and link deflection. Owing to the drawback of serial robot, parallel manipulators have taken a great interest in many applications, such as high-speed machining, assembly and packaging task, flight simulators, and various medical and space applications. ^{3 –7} Parallel robots are generally characterized by their nonanthropomorphic shape, they are composed of an end effector connected to a base by at least two separate and independent kinematic chain. ⁸ This multiple kinematic chain provides a high rigidity and agility with a high payload to weight ratio due to the deportation of the actuators to the base and the distribution of the load between different chain. ⁹ However, they suffer from some drawbacks such as limited work space, abundance of singularities, ¹⁰ and complex dynamic model caused by the presence of multiple closed kinematic chains (CKCs). ¹¹

Many works discussed the kinematic or the dynamic issues of CKCs. ^9,12 But it is still relatively difficult to extend the vast control theory developed for serial manipulator to the parallel one, and this difficulty is due to the complexity of the model-based control. ¹³ This complexity comes from the fact that parallel robot are described by a set of differential–algebraic equations (DAEs) of index-3. ¹⁴ This differential index represents the number of time differentiating the constraint equations to obtain a set of ordinary differential equations (ODEs). ¹⁵ The main difficulty here is the high differential index (index > 2). ^16,17 This makes the DAEs not explicitly expressed in the state–space representation which is a suitable form for most of the control strategies.

In the case where the model of the robot is given on the ODEs form, the simplest way for the control is the famous linear single-axis proportional integral derivative (PID) controller in the joints space. ^{8,18 –22} This method is well-known for its effectiveness under the assumption of locally linear dynamics. It is verified only for low-speed control and cannot be efficient over the whole work space with a single tuning. ²³ An equivalent version of linear single-axis control in the Cartesian space is given in the works by Callegari and colleagues. ^{24 –27} In this case many simplifications are assumed that leads to a lack of accuracy and stability. ²⁷ Another famous approach based on feedback linearization is the computed torque control (CTC). ²⁸ This technique is widely spread in serial robotics ^8,29,30 and has been implemented on several parallel platforms. ^23,31 From this, the literature on the control of parallel robot based on ODE models is very large and mainly depends on the complexity of the studied manipulator and the targeted control strategy such as the sliding mode control, ^32,33 the H^∞ control, ^34,35 the robust TS control, ³ the adaptive control, ^36,37 the passivity-based approach, ³⁸ and other variances. ^{10,23,39 –43} Despite this nonexhaustive scope, models of CKCs are still complex and highly nonlinear in most cases which makes the guarantee of stability in the Lyapunov sense quite difficult.

The aim of this work is to deal with the parallel robot on its DAE form, the difficulty here is that the DAE forms are not expressed explicitly in state–space representation. ⁴⁴ From this, in the second section a DAE model is presented, then the ODE model is deduced to show how complex the ODE model could be. After that, a technique based on the input–output linearization (IOL) is used to write the DAEs model on the state–space representations with an algebraic equation helping to compute the nonlinear parameters that represent the internal forces of the robots. The originality here is that by computing the internal forces the robot can be fully described by two decoupled submodel and then, depending on the complexity of those models, the robot can be controlled either in the joints space or in the operational space. In the third section, a brief description of the Takagi-Sugeno (TS) modeling and control is given. The TS modeling is applied on the legs model of 2-degree-of-freedom (DOF) parallel manipulator by transforming the internal forces as premises for the control law. The interest of bonding the internal forces is to make the system naturally robust against the structured uncertainties just by solving some simple linear matrix inequality (LMI). In the last sections, a comparison between PDC and CTC control laws is presented and some concluding remarks are given.

Mechanical modeling of a 2-DOF parallel manipulator

DAE model of Biglide

The Biglide is a 2-DOF planar parallel manipulator ^31,45 (Figure 1). This parallel kinematic machine consists of two bars mounted through a joint on two separate sliding blocks (active prismatic joints), the extremity of each bar is connected to an end effector, thereby forming a CKC. This configuration allows the positioning of the end effector (operational coordinates: $X={\left[\begin{array}{cc}x& y\end{array}\right]}^T$ ) at a specified point in the operational space by controlling the position $q={\left[\begin{array}{cc}{q}_1& q_2\end{array}\right]}^T$ (active joint coordinates) of the prismatic joints.

Figure 1.

Kinematic diagram of Biglide robot.

The kinematic analysis of the Biglide leads to the following constraint of loop closure

$\phi \left(q,X\right)=0\text{\hspace{1em}with\hspace{1em}}\phi \left(q,X\right)=\left(\begin{array}{l}{\left(x-q_1\right)}^2+y^2-a^2\\ {\left(q_2-x\right)}^2+y^2-a^2\end{array}\right)$

1

The inverse kinematics can thereby be expressed as

$q=g\left(X\right)\text{\hspace{1em}with\hspace{1em}}g\left(X\right)=\left(\begin{array}{l}x-C\left(y\right)\\ x+C\left(y\right)\end{array}\right),C\left(y\right)=\sqrt{a^2-y^2}$

2

Thanks to its simple mechanism, the forward kinematics of the Biglide can be easily derived from equation (2)

${\left[\begin{array}{cc}x& y\end{array}\right]}^T={\left[\frac{q_1+q_2}{2}\sqrt{a^2-\frac{{\left(q_1-q_2\right)}^2}{4}}\right]}^T$

3

Under the assumption that the mechanical system satisfies the loop-closure constraint, the kinematic relationship between end effector velocities and active joints velocities is computed by differentiating the constraint (1) with respect to time, and this relationship is conveniently described by two Jacobian matrices $J_x(X,q)\text{ and }{J}_q(X,q)$ ⁹ as

$J_x\left(X,q\right)\dot{X}=J_q\left(X,q\right)\dot{q}$

4

with $J_x\left(X,q\right)=\left[\begin{array}{cc}x-q_1& y\\ x-q_2& y\end{array}\right],J_q\left(X,q\right)=\left[\begin{array}{cc}x-q_1& 0\\ 0& x-q_2\end{array}\right]$

In order to derive motion equations of the Biglide, the structural characteristics of the robot are exploited to write a model based on two independent subsystems: the end effector dynamics and the open-kinematic chains dynamics. The actuator forces (robot inputs) are denoted by $\Gamma ={\left[\begin{array}{cc}{\Gamma }_1& {\Gamma }_2\end{array}\right]}^T$ (Figure 1). The passive joint coordinates are denoted by $\theta ={\left[\begin{array}{cc}{\theta }_1& {\theta }_2\end{array}\right]}^T$ . The robot is subject to gravity $\overrightarrow{g}=-g\overrightarrow{j}\text{ with }\parallel \overrightarrow{j}\parallel =1$ . The active and passive joint coordinates of each link are denoted by ${\rho }_i={\left[\begin{array}{cc}{q}_i& {\theta }_i\end{array}\right]}^T,i=1,2$ . The dynamics of each leg are denoted by H_i , while the dynamic of the end effector is denoted byF_P . Thereby, the dynamics can be easily obtained from different formalisms such as Newton–Euler equations, Lagrange equations, or virtual work ⁸ as

$H_i=M_{{\rho }_i}\left({\rho }_i\right){\ddot{\rho }}_i+N_{{\rho }_i}\left({\rho }_i,{\dot{\rho }}_i\right){\dot{\rho }}_i+G_{{\rho }_i}\left({\rho }_i\right){\rho }_i$

5

$F_p=M_x\left(X\right)\ddot{X}+G_x\left(X\right)X$

6

where $M_{{\rho }_i}({\rho }_i)=\left[\begin{array}{cc}{m}_i& -m{s}_i\text{sin}{\theta }_i\\ -m{s}_i\text{sin}{\theta }_i& J_i\end{array}\right]\text{ and }{M}_x(X)=\left[\begin{array}{cc}m& 0\\ 0& m\end{array}\right]$ represent the inertia matrices

$N_{{\rho }_i}({\rho }_i,{\dot{\rho }}_i)=\left[\begin{array}{cc}{b}_i& m{s}_i\cos {\theta }_i{\dot{\theta }}_i\\ 0& 0\end{array}\right]$ represents the centrifugal and Coriolis matrix

$G_{{\rho }_i}({\rho }_i)=\left[\begin{array}{cc}0& 0\\ 0& m{s}_{i}g\cos {\theta }_i/{\theta }_i\end{array}\right]\text{ and }{G}_x(X)=\left[\begin{array}{cc}0& 0\\ 0& mg/y\end{array}\right]$ represent the gravity matrices.

With m, m_i , ms_i , J_i , i = 1,2 are the inertial parameters of the robot, which are the end effector mass, the mass, and the first and second moments of link i related to revolute joint axis on the sliding blocks 1 and 2, respectively.

From the kinematic constraint it is easy to write passive joints coordinate functions of operational coordinates.

$\theta =\mu \left(X\right)\text{\hspace{1em}where\hspace{1em}}\mu \left(X\right)={\text{tan}}^{-1}\left(\frac{C\left(y\right)}{y}\right){\left[\begin{array}{cc}1& -1\end{array}\right]}^T$

7

Thereby, the velocities of passive joints can be expressed as

$\dot{\theta }=K_x\left(X\right)\dot{X}\text{\hspace{1em}where\hspace{1em}}{K}_x\left(X\right)=\left[\begin{array}{cc}0& \frac{1}{C\left(y\right)}\\ 0& \frac{-1}{C\left(y\right)}\end{array}\right]$

8

Using the equations (7) and(8), the dynamic of the open-kinematic chain H_i , i = {1,2} given by equation (5) are reordered into active part H_q and passive part H_θ such that

$\left[\begin{array}{c}{H}_q\\ H_{\theta }\end{array}\right]=\left[\begin{array}{cc}{M}_{11}& M_{12}\\ M_{21}& M_{22}\end{array}\right]\left[\begin{array}{c}\ddot{q}\\ \ddot{\theta }\end{array}\right]+\left[\begin{array}{cc}{N}_{11}& N_{12}\\ 0& 0\end{array}\right]\left[\begin{array}{c}\dot{q}\\ \dot{\theta }\end{array}\right]+\left[\begin{array}{cc}0& 0\\ 0& G_{22}\end{array}\right]\left[\begin{array}{c}q\\ \theta \end{array}\right]$

9

with M ₁₁ = diag(m ₁, m ₂), $M_{12}=M_{21}=diag(-m{s}_{1}sin{\theta }_1,-m{s}_{2}sin{\theta }_2),M_{22}=diag(J_1,J_2),N_{11}=diag(b_1,b_2),N_{12}=diag(m{s}_{1}cos{\theta }_1{\dot{\theta }}_1,m{s}_{2}cos{\theta }_2{\dot{\theta }}_2),G_{22}=diag(m{s}_{1}g\cos {\theta }_1/{\theta }_1,m{s}_{2}g\cos {\theta }_2/{\theta }_2)$ .

The virtual power $P_r^{*}$ of the overall dynamic of the robot given by the equations 6 and 9 ^11,46 is first written with a dependent virtual velocity field $({\dot{q}}^{*T},{\dot{\theta }}^{*T},{\dot{X}}^{*T})$ such that

$\begin{array}{l}{P}_r^{*}={\dot{q}}^{*T}\left(H_q-\Gamma \right)+{\dot{X}}^{*T}{F}_p+{\dot{\theta }}^{*T}{H}_{\theta }=0\hfill \end{array}$

10

Using equation (8), the virtual passive velocities is eliminated with ${\dot{\theta }}^{*T}={\dot{X}}^{*T}{K}_x^T$ , leading to the following form

$\begin{array}{l}{P}_r^{*}={\dot{q}}^{*T}\left(H_q-\Gamma \right)+{\dot{X}}^{*T}\left(F_p+K_x^T{H}_{\theta }\right)=0\hfill \end{array}$

11

The new virtual velocities field $({\dot{q}}^{*T},{\dot{X}}^{*T})$ is still dependent. Therefore, we introduce the Lagrange multipliers $\lambda ={\left[\begin{array}{cc}{\lambda }_1& {\lambda }_2\end{array}\right]}^T$ such that, the virtual power of the constraint (1) is expressed with equation (4) as follows

$P_c^{*}=\left({\dot{q}}^{*T}{J}_q^T-{\dot{X}}^{*T}{J}_x^T\right)\lambda =0,\forall \lambda ={\left[\begin{array}{cc}{\lambda }_1& {\lambda }_2\end{array}\right]}^T\in {ℝ}^2$

12

Finally, the virtual power with the Lagrange multipliers is expressed as follows

$P^{*}=P_r^{*}+P_c^{*}={\dot{q}}^{*T}(H_q+J_q^T\lambda -\Gamma )+{\dot{X}}^{*T}(F_p+K_x^T{H}_{\theta }-J_x^T\lambda )=0$

13

Given that the Lagrange multipliers are unknown, the virtual velocities field $({\dot{q}}^{*T},{\dot{X}}^{*T})$ may be considered as independent and then the index-1 DAEs model of the constrained system in terms of (q, X) is derived from equation (13) as

$P^{*}=0\forall {\left[{\dot{q}}^{*T}{\dot{X}}^{*T}\right]}^T\in {ℝ}^4\Rightarrow \{\begin{array}{l}{M\prime }_{11}\ddot{q}+{M\prime }_{12}\ddot{X}+{N\prime }_{11}\dot{q}+{N\prime }_{12}\dot{X}+J_q^T\lambda =\Gamma \hfill \\ {M\prime }_{21}\ddot{q}+{M\prime }_{22}\ddot{X}+{N\prime }_{22}\dot{X}+{G\prime }_{22}X-J_x^T\lambda =0\hfill \\ J_q\dot{q}-J_x\dot{X}=0\hfill \end{array}$

14

with ${M\prime }_{11}=M_{11},{M\prime }_{12}=M_{12}{K}_x,{M\prime }_{21}=K_x^T{M}_{21},{M\prime }_{22}=M_x+K_x^T{M}_{22}{K}_x$ ,

${N\prime }_{11}=N_{11},{N\prime }_{12}=M_{12}{\dot{K}}_x+N_{12}{K}_x,{N\prime }_{22}=K_x^T{M}_{22}{\dot{K}}_x+K_x^T{N}_{22}{K}_x,{G\prime }_{22}=G_x+K_x^T{G}_{22}\mu$ .

Input–output linearization

The DAEs model given by equation (14), presents the coupled dynamics between the legs and the end effector. The algebraic variables (Lagrange multipliers λ_i) represent the internal forces between the legs and the end effector. The assembly of the two subsystems is ensured through the algebraic constraint. This model is singular, the difficulty here is that DAEs are not expressed explicitly in state–space representation. ⁴⁴ From this, two main approaches have been developed. The first one is the approximation of the DAE model to a singularly perturbed model. ⁴⁷ The idea is to substitute the constraint equation with a fast dynamics representing the violation of the constraint. ^{44,48 –50} From a practical point of view, this approach is completely justified because the connections between joints and links are elastic. ⁵¹ This method helps to relax the complexity of the nonlinear model through the additional artificial dynamics, but at the same time it extends the state vector. This leads to an increase in the complexity of the control applications based on the resolution of the LMI problem. The second technique is the IOL, ⁵² and this approach consists of solving the algebraic equation by differentiating the constraint until the appearance of the algebraic variables. The solution is used in the differential equation of the DAE system and a state–space model-based control is given. The differentiation of the constraint (14) is written as

$J_q\ddot{q}-J_x\ddot{X}+{\dot{J}}_q\dot{q}-{\dot{J}}_x\dot{X}=0$

15

By combining (15) with the DAEs model (14), a redundant ODEs model is obtained:

$\{\begin{array}{l}{M}_q\ddot{q}+N_q\dot{q}+J_q^T\lambda =\Gamma \hfill \\ {M\prime }_x\ddot{X}+{N\prime }_x\dot{X}+{G\prime }_{x}X-J_x^T\lambda =0\hfill \\ J_q\ddot{q}-J_x\ddot{X}+{\dot{J}}_q\dot{q}-{\dot{J}}_x\dot{X}=0\hfill \end{array}$

16

with $M_q={M\prime }_{11}+{M\prime }_{12}{J}_x^{-1}{J}_q,{M\prime }_x={M\prime }_{22}+{M\prime }_{21}{J}_q^{-1}{J}_x$ ,

$N_q={N\prime }_{11}+{N\prime }_{12}{J}_x^{-1}{J}_q+{M\prime }_{12}{J}_x^{-1}{\dot{J}}_q-{M\prime }_{12}{J}_x^{-1}{\dot{J}}_x{J}_x^{-1}{J}_q,{N\prime }_x={N\prime }_{22}+{M\prime }_{21}{J}_q^{-1}{\dot{J}}_x-{M\prime }_{21}{J}_q^{-1}{\dot{J}}_q{J}_q^{-1}{J}_x,{G\prime }_x={G\prime }_{22}$ .

From system of equations (16) and by placing the differential equations into the algebraic one, the Biglide parallel robot can be fully described by two different subsystems as follows

$\text{Legs\hspace{0.17em}model\hspace{0.05em}}\{\begin{array}{l}{M}_q\ddot{q}+N_q\dot{q}+J_q^T\lambda =\Gamma \hfill \\ \mathcal{P}\lambda +\mathcal{L}\dot{q}-\mathcal{Q}q-\mathcal{W}\Gamma =0\hfill \end{array}$

17

$\text{Endeffector\hspace{0.17em}model\hspace{0.05em}}\{\begin{array}{ll}{M\prime }_x\ddot{X}+{N\prime }_x\dot{X}+{G\prime }_{x}X-J_x^T\lambda =0\hfill & \hfill \\ \mathcal{P}\lambda +\mathcal{L}\prime \dot{X}-\mathcal{Q}\prime X-\mathcal{W}\Gamma =0\hfill & \hfill \end{array}$

18

with $\mathcal{P}=J_q{M}_q^{-1}{J}_q^T+J_x{M\prime }_x{}^{-1}{J}_x^T,\mathcal{L}=J_q{M}_q^{-1}{N}_q-{\dot{J}}_q+{\dot{J}}_x{J}_x^{-1}{J}_q-J_x{M}_x^{-1}{N\prime }_x{J}_x^{-1}{J}_q,\mathcal{L}\prime =\mathcal{L}{J}_q^{-1}{J}_x,\mathcal{Q}\prime =J_x{M}_x^{-1}{G\prime }_x,Q=Q\prime g^{-1},\mathcal{W}=J_q{M}_q^{-1}$ .

This is the first main results of this article, the Biglide robot can be controlled by two different approaches depending on the targeted control law. The first approach is to only use the legs model equation (17). In this case, the internal forces λ are seen as nonlinear terms to be compensated by the controller. The second approach is to only consider the end effector model (18). In this case, the internal forces of the system are seen as a new input control defined by the algebraic equation.

ODEs model of Biglide

In order to get the ODEs model of Biglide in operational space, the virtual active joint velocities are eliminated from the virtual power equation (11) using the equation (4), ${\dot{q}}^{*T}={\dot{X}}^{*T}{J}_x^T{J}_q^{-T}$ , the virtual power can be rewritten with the independent velocities ${\dot{X}}^{*}$ as

$\begin{array}{l}{P}_r^{*}={\dot{X}}^{*T}\left(J_x^T{J}_q^{-T}\left(H_q-\Gamma \right)+F_p+K_x^T{H}_{\theta }\right)=0\forall {\dot{X}}^{*}\in {\mathbb{R}}^2\hfill \end{array}$

19

Therefore, the ODEs model can be obtained from equation (19) as $\Gamma =H_q+J_q^T{J}_x^{-T}(F_p+K_x^T{H}_{\theta })$ . Then to eliminate $\ddot{\theta }\text{ and }\ddot{q}$ , the first derivative of equations (19) and (8) are used. Finally the ODE model of Biglide in operational space is given as

$\Gamma =M\left(X\right)\ddot{X}+N\left(X,\dot{X}\right)\dot{X}+G\left(X\right)X$

20

with $M(X)=\left[\begin{array}{cc}{m}_1+\frac{1}{2}(m-{\eta }_1+{\eta }_2)& f_1(X)\\ m_2+\frac{1}{2}(m-{\eta }_2+{\eta }_1)& f_2(X)\end{array}\right],N(X,\dot{X})=\left[\begin{array}{cc}{n}_{11}& n_{12}\\ n_{21}& n_{22}\end{array}\right],G(X)=\left[\begin{array}{cc}0& gC(y)(m+{\eta }_1+{\eta }_2)/2y^2\\ 0& -gC(y)(m+{\eta }_1+{\eta }_2)/2y^2\end{array}\right],f_1(X)=\left[(2m_1-3{\eta }_1-{\eta }_2)y^2+mC{(y)}^2+J_1+J_2\right]/2yC(y),f_2(X)=-\left[(2m_2-3{\eta }_2-{\eta }_1)y^2+mC{(y)}^2+J_1+J_2\right]/2yC(y),{\eta }_{1,2}=m{s}_{1,2}/a$ , $n_{11}=b_1,n_{21}=b_2,n_{12}=b_{1}y/C(y)-\left[(2m_1-3{\eta }_1-{\eta }_2)y^2+(2m_1-3{\eta }_1-{\eta }_2)C{(y)}^2+J_1+J_2\right]\dot{y}/(2C{(y)}^3,n_{22}=-b_{2}y/C(y)+\left[(2m_2-3{\eta }_2-{\eta }_1)y^2+(2m_2-3{\eta }_2-{\eta }_1)C{(y)}^2+J_1+J_2\right]\dot{y}/(2C{(y)}^3$ .

In a more general case, it is not obvious to get the ODEs model of parallel robot in operational space. To get the model in the active joint space, the same logic is followed by eliminating the virtual operational velocities. Note that the ODEs model is given here just for a comparison purpose with the DAEs one. Unlike the DAEs model, where most of the nonlinear terms are hidden into the algebraic equation, the ODEs model presents a large number of nonlinear terms. For this reason, we thought that it will be better to use the DAEs model to design nonlinear controllers.

PDC fuzzy control

TS modeling

The TS fuzzy model is a mathematical representation of systems, it belongs to the quasi LPV family. ⁵³ Inside a compact set of state variables, TS fuzzy model can represent exactly a nonlinear system by a collection of linear models weighted together by a nonlinear function issued from nonlinearities of the system. ⁵⁴ The conditions of stability and stabilization of these models are generally based on the Lyapunov theory. ⁵⁵ The advantage of this representation is that it provides a systematic framework for designing control laws through the LMI constrain formulation. ⁵⁶

Let’s consider the following nonlinear plant

$\{\begin{array}{ll}\dot{x}\left(t\right)=Ax\left(t\right)+Bu\left(t\right)\hfill & \hfill \\ y\left(t\right)=Cx\left(t\right)\hfill & \hfill \end{array}$

21

where $x(t)\in {ℝ}^n$ is the state vector, $u(t)\in {ℝ}^m$ is the control input vector A, B, and C are nonlinear matrices of proper dimension.

Let r be the number of nonlinear terms of the system (21). A TS model is viewed as a convex sum of linear models via membership functions (MFs), it is defined by «If…Then» rules representing the local linear input–output relations of the plant. ⁵³ From this, the nonlinear system (21) can be written on continuous-time TS fuzzy model as follows ⁵⁶

$\{\begin{array}{l}\dot{x}\left(t\right)={\displaystyle \sum _{i=1}^r}{h}_i(z\left(t\right))\left[A_{i}x\left(t\right)+B_{i}u\left(t\right)\right]\hfill \\ y\left(t\right)={\displaystyle \sum _{i=1}^r}{h}_i(z\left(t\right))C_{i}x\left(t\right)\hfill \end{array}$

22

where z(t) is the premise vector (vector of nonlinear terms of the model). It may depend on the state, input, exogenous parameters, or time, on measurable and/or unmeasurable variables. Matrices (A_i , B_i ), $i\in \{1,2,\mathrm{...},r\}$ represent the ith linear model of the TS form (22). The MFs $h_i(z(t))\ge 0,i\in \{1,\mathrm{...},r\}$ are nonlinear functions and satisfy the convex property as follows ${\displaystyle \sum _{i=1}^r}{h}_i(z\left(t\right))=1$

TS model of the legs

For the following, let’s consider the dynamical model of the Biglide given by equation (17). In order to decompose the internal forces λ into a convex sum, the term $J_q^T\lambda$ is written as

$J_q^T\lambda =\psi \left(\lambda \right)q\text{\hspace{1em}with\hspace{1em}}\psi \left(\lambda \right)=\left[\begin{array}{cc}\frac{-{\lambda }_1}{2}& \frac{{\lambda }_1}{2}\\ \frac{{\lambda }_2}{2}& \frac{-{\lambda }_2}{2}\end{array}\right]$

23

Then, the state–space representation of the Biglide robot in term of joint variables $(q,\dot{q})$ is given as

$\left[\begin{array}{c}{\dot{q}}_1\\ {\dot{q}}_2\\ {\ddot{q}}_1\\ {\ddot{q}}_2\end{array}\right]=\left[\begin{array}{cc}I& 0\\ 0& M_q^{-1}\end{array}\right]\left(\left[\begin{array}{cccc}0& 0& 1& 0\\ 0& 0& 0& 1\\ \frac{{\lambda }_1}{2}& \frac{-{\lambda }_1}{2}& -b_1& 0\\ \frac{-{\lambda }_2}{2}& \frac{{\lambda }_2}{2}& 0& -b_2\end{array}\right]\left[\begin{array}{c}{q}_1\\ q_2\\ {\dot{q}}_1\\ {\dot{q}}_2\end{array}\right]+\left[\begin{array}{cc}0& 0\\ 0& 0\\ 1& 0\\ 0& 1\end{array}\right]\Gamma \right)$

24

where $I\in {ℝ}^{2\times 2}$ is an identity matrix.

From the state representation (24), the TS fuzzy model of the Biglide can be obtained by considering the following two nonlinearities ${\lambda }_1(\cdot )\text{ and }{\lambda }_2(\cdot )$ . From this it yields a 4-rule (2 ^r ) TS fuzzy model. Thanks to the DAE model, the number of local models have been reduced from 16-rule TS descriptor model ³ to a set of 4-rule classical TS models. Note that, ${\lambda }_i(\cdot )$ , i = {1, 2} are computed via the algebraic equation (17). A small delay ε_t is introduced to break the algebraic loop. This delay can be justified by supposing that the dynamics of the system are slow enough to not react between t and (t + ε_t ). It is also known that the control inputs are subject to saturation $\left|{\Gamma }_{1,2}\right|⩽{\Gamma }_{\text{max}}$ , the working space is restrained with the range of $q_i\in [{\underset{\_}{q}}_i,{\overline{q}}_i]$ , i = {1, 2}, and the velocities are clamped $\left|{\dot{q}}_{1,2}\right|⩽{\dot{q}}_{\text{max}}$ . From this consideration, the premises ${\lambda }_i(\cdot )$ , i = {1, 2} can be bounded as ${\lambda }_i(\cdot )\in [\underset{\_}{{\lambda }_i}\overline{{\lambda }_i}]$ . The linear submodels are obtained using the sector nonlinearity transformation ⁵⁶ by writing λ as

$\lambda \left(\mathrm{.}\right)=w_0^i\left(\mathrm{.}\right)\cdot \underset{\_}{{\lambda }_i}+w_1^i\left(\mathrm{.}\right)\cdot \overline{{\lambda }_i}$

25

with $w_0^i(\mathrm{.})=(\overline{{\lambda }_i}-\lambda (\mathrm{.}))/(\overline{{\lambda }_i}-\underset{\_}{{\lambda }_i})\mathrm{, }{w}_1^i(\mathrm{.})=(\lambda (\mathrm{.})-\underset{\_}{{\lambda }_i})/(\overline{{\lambda }_i}-\underset{\_}{{\lambda }_i}),w_1^i(\mathrm{.})=(\lambda (\mathrm{.})-\underset{\_}{{\lambda }_i})/(\overline{{\lambda }_i}-\underset{\_}{{\lambda }_i})$ . Note that the functions $h_i(z(t))⩾0$ are expressed as products between the $w_j^i,i\in \{1,\mathrm{...},r\},j\in \{1,0\}$ .

Controller design

For the stabilization of system (22), a parallel distributed compensation (PDC) control law is proposed ⁵⁷

$u\left(t\right)={\displaystyle \sum _{i=1}^r}{h}_i\left(z\left(t\right)\right)F_{i}x\left(t\right)$

26

where the matrices $F_i\in {ℝ}^{m\times n}$ need to be fixed to ensure the stabilization of the system. ⁵⁸

Theorem 1: The fuzzy TS model (22) with the PDC control law (26) is asymptotically stable if there exist a matrix P = P^T >0 such as ^56,59

$\{\begin{array}{cc}{Y}_{ii}<0& \forall i\in \{1,\mathrm{...},r\}\\ \frac{1}{r-1}{Y}_{ii}+\frac{1}{2}\left(Y_{ij}+Y_{ji}\right)<0& \forall i,j\in \{1,\mathrm{...},r\}i<j\end{array}$

27

with $Y_{ij}=A_{i}P+P{A}_i^T-B_i{M}_j-M_j^T{B}_i^T$ . If the LMI problem is feasible, the PDC controller gain matrices are given by $F_i=M_i{P}^{-1}$ .□

An integral part is added (Figure 2) in order to compensate the stationary error. Consider the extended state vector ${\overline{x}}^T=\left[\begin{array}{cc}{x}^T& x_l^T\end{array}\right]$ , then ³

Figure 2.

Extended PDC control scheme.

${\dot{x}}_l=y_d-{\displaystyle \sum _{i=1}^r}{h}_i(z\left(t\right))C_i\text{\hspace{0.17em}x}\left(t\right)$

28

where y_d is the desired input vector. The extended system is expressed as

$\{\begin{array}{l}\overline{\dot{x}}\left(t\right)={\displaystyle \sum _{i=1}^r}{h}_i(z\left(t\right))\left[{\overline{A}}_i\overline{x}\left(t\right)+{\overline{B}}_{i}u\left(t\right)\right]+B_0{y}_d(t)\hfill \\ y\left(t\right)={\displaystyle \sum _{i=1}^r}{h}_i(z\left(t\right)){\overline{C}}_i\overline{x}\left(t\right)\hfill \end{array}$

29

with ${\overline{A}}_i=\left[\begin{array}{cc}{A}_i& 0\\ -C_i& 0\end{array}\right],{\overline{B}}_i=\left[\begin{array}{c}{B}_i\\ 0\end{array}\right],{\overline{C}}_i=\left[\begin{array}{cc}{C}_i& 0\end{array}\right],B_0=\left[\begin{array}{c}0\\ I\end{array}\right]$ ; where 0 and I are matrices of proper dimension. Therefore, the extended PDC control law can be written as

$u\left(t\right)={\displaystyle \sum _{i=1}^r}{h}_i(z\left(t\right)){\overline{F}}_i\overline{x}\left(t\right)$

30

with the extended gain vector ${\overline{F}}_i=\left[\begin{array}{cc}{F}_i& L_i\end{array}\right]$ .

Simulation results

The model of the parallel robot used for numerical simulations includes structured and unstructured uncertainties (Appendix 1). The structured uncertainties represent the variation in the end-effector mass, this variation is given by $\Delta m\in \left[00.816\right]kg$ . To simulate a realistic behavior, The numerical model includes unmodeled dynamics such as elasticity between joints ⁶⁰ and the Stribeck friction ⁶¹ (Appendix 2). Furthermore, two resonant modes are added to simulate the elastic joints, the lower value of the resonant frequency is $\omega \simeq 29rad/s$ . To compare, a CTC controller is designed. It is a well-known approach for serial robotics. ^8,29,30 It has also been implemented on several platforms of parallel kinematic machines. ^23,31,62 Usually, the tuning of the CTC controller with a PID control action uses a pole placement technique for robot manipulators. ²³ The gains are adjusted in order to get a negative real triple pole with a frequency less than the half of the lower resonant frequency. ⁸ The suitable value is used for simulation. For the TS approach, the PDC controller was tuned in the same manner as the CTC with a convex optimization algorithm via the LMI solver of Matlab LMI Toolbox. ⁶³

In order to get the TS numerical model, the joint positions were bounded in a manner to avoid singularities and to cover a large area of the working space $q_1\in \left[54.5,89.5\right]\text{\hspace{0.17em}mm}\text{ and }{q}_2\in \left[129.5,164.5\right]\text{\hspace{0.17em}mm}$ . The velocity of each joint was bounded as $\left|{\dot{q}}_{1,2}\right|\le 428\text{\hspace{0.17em}mm}/\text{s}$ . The control inputs was saturated as $\left|{\Gamma }_{1,2}\right|\le 50\text{\hspace{0.17em}N}$ . The internal forces λ were bounded using the algebraic equation (17). To guarantee the robustness against the structured uncertainty, the internal forces were bounded in the worst case (where the payload is equal to the global mass of the robot).

For the simulation, we chose a circular trajectory in the operational space. The response of both controllers, CTC and TS, are depicted in Figures 3 and 4 for the model without payload $(\Delta m=0\text{\hspace{0.17em}kg})$ , Figures 5 and 6 for a payload representing 50% of the global mass $(\Delta m=0.816\text{\hspace{0.17em}kg})$ , and Figures 7 and 8 for a payload representing 100% of the global mass $(\Delta m=1.632\text{\hspace{0.17em}kg})$ .

Figure 3.

Simulation results for $\Delta m=0kg$ : (a) velocity profile of $\dot{x}$ , (b) velocity profile of $\dot{y}$ , (c) torque profile of joint 1, (d) torque profile of joint 2, (e) internal force profile of leg 1, and (f) internal force profile of leg 2.

Figure 4.

Simulation results for $\Delta m=0kg$ : (a) trajectory profile and (b) radial error.

Figure 5.

Simulation results for $\Delta m=0.816kg$ : (a) velocity profile of x, (b) velocity profile of y, (c) torque profile of joint 1, (d) torque profile of joint 2, (e) internal force profile of leg 1, and (f) internal force profile of leg 2.

Figure 6.

Simulation results for $\Delta m=0.816kg$ : (a) trajectory profile and (b) radial error.

Figure 7.

Simulation results for $\Delta m=1.632kg$ : (a) velocity profile of x, (b) velocity profile of y, (c) torque profile of joint 1, (d) torque profile of joint 2, (e) internal force profile of leg 1, and (f) internal force profile of leg 2.

Figure 8.

Simulation results for $\Delta m=1.632kg$ : (a) trajectory profile and (b) radial error.

In the case where the robot does not carry any payload, both control laws ensure a good tracking with a slight advantage to the TS control law, so we can say that the TS and CTC are robust against the unmodeled dynamics (stick/slip friction). Moreover, from Figure 3 (c and d) we can note that the control inputs are almost similar, this indicates that the two controllers were tuned in a manner to get similar performances in the nominal case. In the second case, $(\Delta m=0.816\text{\hspace{0.17em}kg})$ , the TS control law ensures good tracking performances, while the CTC controller loses the desired path but the robot remains stable. In the latter case, $(\Delta m=1.632\text{\hspace{0.17em}kg})$ , the TS-based control law presents the best results exhibiting its robustness according to structured uncertainties. The CTC, as expected, presents an unstable trajectory.

Finally, two well-known criteria are computed over the simulation time (T = 0.5 s) in order to quantify the behavior of both controllers (Table 1): the first criterion is the integral of absolute error: $J_{IAE}={\displaystyle {\int }_0^T}\left|\epsilon \right|dt$ . where ε(t) is the radial error. The second criterion is the integral of square value of the control: $J_{ISV}={\displaystyle {\int }_0^T}{\displaystyle \sum _{i=1}^2}{\Gamma }_i{(t)}^{2}dt$ .

Table 1.

ISV and IAE criteria.

	J ISV		J IAE
	TS	CTC	TS	CTC
Δm=0 kg	50.05	55.2	0.4287	28.06
Δm=0.816 kg	58.84	58.23	0.6281	253.7
Δm=1.632 kg	277.33	446.6	1.458	943.6

ISV: integral of square value of the control; IAE: integral of absolute error; TS: Takagi-Sugeno; CTC: computed torque control.

Conclusion

This article proposes a novel approach for the control of parallel robots. Those systems are particularly difficult to control due to the presence of multiple CKCs. These mechanical loops lead to a model naturally described by a set of DAEs of index-3. Designing a controller based on the DAE model is not an easy task, it requires a particular knowledge on the control theory of singular systems. For this reason, most of the time we transform the DAE model into an ODE one by differentiating the constraint twice. This manipulation is appropriate for mechanical systems because if the system is well designed, the initial conditions should respect the algebraic constraint. Unfortunately, the obtained ODE model is most of the time complex and highly nonlinear. This makes the design of nonlinear controllers that requires the use of convex optimization techniques such as TS modeling and LMIs quite difficult. Because, in this case, the complexity of the controller is an exponential function of the number of nonlinear terms. For this reason, this article proposes to keep the DAE form of the system in order to decouple it into two independent subsystems. The idea here is to collect most of the nonlinear terms inside the Lagrange multipliers ${\lambda }_i(\cdot )$ . According to our way of decomposing, the Lagrange multipliers represent the internal forces of the system. Knowing that all the physical parameters are bounded, the internal forces were bounded according to the range of the state variables and the saturation of the control inputs. From those bounds, the internal forces were decomposed into a convex sum and used as premises to establish the TS model. A PDC control law was designed. An extended integral part was added in order to compensate the stationary errors. Finally a comparison with the well-known CTC control law has been presented. The robustness of the PDC controller has been well checked against the structural uncertainties.