2.1. Analyze of PMA’s structure

The PMA used in the experiment is the pneumatic tendon developed by FESTO company, as shown in Fig. 1. The PMA is mainly composed of an internal rubber tube, a peripheral braided sleeve, and a fastening head at both ends. The rubber hose has good elasticity while the braid sleeve has great rigidity. When the rubber hose is pressurized with gas, its radial expansion and axial shortening will generate external tension, so as to convert the gas pressure energy into mechanical energy and drive the load.

Figure 1. FESTO pneumatic tendon. (a) Physical product. (b) Structure diagram.

Figure 2. Static test schematic diagram of single PMA.

During the process of inflatable contraction, the pneumatic muscle expands radially and contracts axially at the same time. It can be concluded that the length of a single fiber l stays constant [Reference Gu, Liu and Ma13], as shown in Fig. 1(b), because of the very large stiffness of the fiber sheath. Thus, the relationship among the geometrical parameters of fiber length, muscle length, and radial circumference can be obtained as follows:

(1) \begin{align}\left\{ {\begin{array}{*{20}{c}}{{L_0} = lcos{\theta _0}}\\[5pt] {n\pi {D_0} = lsin{\theta _0}}\end{array}} \right. \Rightarrow \left\{ {\begin{array}{*{20}{c}}{L = lcos\theta }\\[5pt] {n\pi D = lsin\theta }\end{array}} \right.\end{align}

where let ${L_0}$ and $L$ be the initial length and real-time length of the pneumatic muscle respectively, ${D_0}$ and $D$ be the initial diameter and real-time diameter of the aerodynamic muscle section, $l$ be the length of a single fiber of pneumatic muscle, $n$ be the number of winding cycles of a single fiber, $\;{\rm{and\;}}{\theta _0}$ and $\theta $ be the angle between the fiber and the axis at the initial moment and in real time, respectively.

2.2. The mathematical model of PMA’s tensile

In order to achieve accurate control of PMA, the first step is to establish the mathematical model, namely, the mathematical model between the load tension ${\rm{\;}}F$ , the shrinking percentage $\varepsilon $ , and the charging pressure $P$ . In this paper, a static examine platform (as shown in Fig. 2) is used for static test of a single pneumatic muscle, which can be divided into contraction test and tensile test. During the experiment, contraction test refers to continuously inflating air to make the pneumatic muscle contract from the original length to the limit position, while tensile test refers to continuously deflating air to restore the original length from the limit position of contraction.

MATLAB curve fitting toolbox was used to fit the test data of 110 mm pneumatic muscle in the process of contraction and stretching, the fitting degree was around 0.99, and the mathematical models obtained are as follows:

(2) \begin{align}\left\{ \begin{matrix} {\Delta P < 0:\;\;\; F = P\left( {0.0305{\varepsilon ^3} - 1.134{\varepsilon ^2} - 14.72\varepsilon + 819.9} \right) + } \\[5pt] {\left( {0.00257{\varepsilon ^4} - 0.1758{\varepsilon ^3} + 4.263{\varepsilon ^2} - 42.54\varepsilon - 0.28} \right)} \\[5pt] {\Delta P \ge 0:\;\;\; F = P\left( {0.0305{\varepsilon ^3} - 1.134{\varepsilon ^2} - 14.72\varepsilon + 819.9} \right) + } \\[5pt] {\left( {0.00257{\varepsilon ^4} - 0.1758{\varepsilon ^3} + 4.263{\varepsilon ^2} - 42.54\varepsilon + 0.21} \right)} \end{matrix} \right.\end{align}

3.1. PI robust sliding mode control

Sliding mode control is a kind of variable structure control strategy. The control system structure is not fixed but constantly changes with time, which is essentially different from the conventional control strategy. During actual motion process of flexible manipulator, many uncertainties exist in the dynamic model, where $M(q),C({q,\dot q})$ ${\rm{and\;}}G(q){\rm{\;}}$ are usually inaccurate, which can be expressed as

(3) \begin{align}\left\{ {\begin{array}{*{20}{c}}{M(q) = {M_0}(q) + {E_M}}\\[5pt] {C({q,\dot q}) = {C_0}({q,\dot q}) + {E_c}}\\[5pt] {G(q) = {G_0}(q) + {E_G}}\end{array}} \right.\end{align}

where let ${E_M},{E_c},{\rm{\;}}\textrm{and}{\rm{\;}}{E_G}$ be the modeling errors of $M(q),C({q,\dot q}),{\rm{\;and\;}}G(q),$ respectively, ${\rm{\;}}{M_0}(q),{C_0}({q,\dot q}),{\rm{\;and\;}}{G_0}(q)$ be the estimated values of the inertia coefficient, the Coriolis force, the centrifugal force term, and the gravity term, respectively, which is so-called nominal model.

Noted ${q_d}(t)$ as the ideal angle signal of the joint and $q(t)$ as the actual angle signal. So, the tracking error of the two joint angles is defined as

(4) \begin{align}e(t) = {q_d}(t) - q(t)\end{align}

The sliding mode function is defined as follows (where $\Lambda \gt 0$ ):

(5) \begin{align}r = \dot e + \Lambda e\end{align}

Noted ${\dot q_r} = r(t) + \dot q(t)$ , then ${\ddot q_r} = \dot r(t) + \ddot q(t)$ , ${\dot q_r} = {\dot q_d} + \Lambda e$ , ${\ddot q_r} = {\ddot q_d} + \Lambda \dot e$ . So, the following can be obtained from the dynamic model of flexible manipulator:

(6) \begin{align}\tau & = M(q)\ddot q + C({q,\dot q})\dot q + G(q) + {\tau _d}\nonumber \\ & = {M_0}(q){\ddot q_r} + {C_0}({q,\dot q}){\dot q_r} + {G_0}(q) - M(q)\dot r - C({q,\dot q})r + E^{\,\prime} + {\tau _d}\end{align}

where $E^{\,\prime} = {E_M}{\ddot q_r} + {E_c}{\dot q_r} + {E_G}$ .

The control rate of PI robust sliding mode control is designed as follows:

(7) \begin{align}\tau = {\tau _m} + {K_p}r + {K_i}\int rdt + {\tau _r}\end{align}

where ${K_p} \gt 0$ , ${K_i} \gt 0$ , and ${\tau _m}$ are torque vector control based on the name of the flexible manipulator model, ${K_p}r + {K_i}\smallint rdt$ is PI sliding mode surface, and ${\tau _r}$ is robust term.

According to $\left( 7 \right)$ , the following can be deduced:

(8) \begin{align}\left\{ {\begin{array}{*{20}{c}}{{\tau _m} = {M_0}(q){{\ddot q}_r} + {C_0}({q,\dot q}){{\dot q}_r} + {G_0}(q)}\\[5pt] {{\tau _r} = {K_r}sgn(r)}\end{array}} \right.;\end{align}

where $K_{r} = diag[k_{r\;ii}],\;\; k_{r\;ii}\geq |E_{i}|,\; i = 1, 2, E = E^{\prime} + \tau_{d}$ .

Further, according to (6) and (8), Eq. (9) can be deduced:

(9) \begin{align}&{M_0}(q){\ddot q_r} + {C_0}({q,\dot q}){\dot q_r} + {G_0}(q) - M(q)\dot r - C({q,\dot q})r + E^{\,\prime} + {\tau _d}\, \nonumber\\ & \quad = {M_0}(q){\ddot q_r} + {C_0}({q,\dot q}){\dot q_r} + {G_0}(q) + {K_p}r + {K_i}\int rdt + {K_r}sgn(r)\;\end{align}

Simplified both two sides of the above formula, Eq. (10) can be deduced:

(10) \begin{align}M(q)\dot r + C({q,\dot q}){r_r} + {K_i}\int rdt = - {K_p}r - {K_r}sgn(r) + E\end{align}

3.2. Stability of the PI robust sliding mode control

Lyapunov function based on integral type is designed as follows:

(11) \begin{align}V = \frac{1}{2}{r^{\rm{T}}}Mr + \frac{1}{2}{\left( {\int rd\tau } \right)^{\rm{T}}}{K_i}\left( {\int rd\tau } \right)\end{align}

Obviously, $V$ is positive definite, then:

(12) \begin{align}\dot V = {r^{\rm{T}}}\!\left({Mr + \frac{1}{2}Mr + {K_i}\int rd\tau}\right)\end{align}

Considering Skew symmetry of kinetic equation, ${r^{\rm{T}}}({M - 2C})r = 0$ , then:

(13) \begin{align}\dot V = {r^{\rm{T}}}\left( {Mr + Cr + {K_i}\int rd\tau } \right)\end{align}

Substituting Eq. (10) into Eq. (13):

(14) \begin{align}\dot V = - {r^{\rm{T}}}{K_p}r - {r^{\rm{T}}}{K_r}sgn(r) + {r^{\rm{T}}}E = - {r^{\rm{T}}}{K_p}r - \mathop \sum \limits_{i = 1}^2 {K_{rii}}{\left| r \right|_i} + {r^{\rm{T}}}E\end{align}

where ${K_{rii}} \ge \left| {{E_i}} \right|$ , so $\dot V \le - {r^{\rm{T}}}{K_p}r \le 0$ .

Figure 3. Flow chart of the fuzzy classifier.

Figure 4. Fuzzy sliding mode control.

Figure 5. Simulation model of flexible manipulator.

Figure 6. Tracking degree error of α $_{1}$ , β $_{1}$ , and γ $_{1}$ at scapula joint.

Figure 7. Tracking degree error of α $_{2}$ , β $_{2}$ , and γ $_{2}$ at humeral joint.

Because $\dot V$ is a negative definite matrix, the control system is stable.

3.3. Improvement of the sliding mode control

The PI robust sliding mode control has its inherent shortcoming that is the phenomenon high-frequency chattering, which can only be reduced, not eliminated. Therefore, a fuzzy sliding mode control strategy based on the saturation function is proposed to improve the robustness of PI robust sliding mode control. Compared with the method of eliminating high-frequency chattering directly by using saturated function [Reference Anh, Kien, Son and Nam14], we use fuzzy control strategy to modify the robust term coefficient, so as to ensure the rapidity and robustness of the sliding mode control.

From (8), the robust term of PI robust sliding mode control is as follows:

(15) \begin{align}{\tau _r} = {K_r}sgn(r)\end{align}

When ${{r}} \to 0$ , ${K_r}{\textrm{sgn}}({r})$ does not approach zero, and the system state will move in and out frequently on both sides of the sliding switching surface, which caused chattering of the control system, and the strength of chattering depends on the size of the robust term coefficient ${K_r}$ . In order to effectively reduce chattering, a fuzzy controller is adopted to adjust the size of the robust term coefficient ${K_r}$ . The flow chart of the fuzzy classifier is shown in Fig. 3.

Table I. Simulation error range of scapula joint position tracking.

Table II. Simulation error range of humeral joint location tracking.

Figure 8. Actual flexible manipulator.

In the actual control process, because the system error and external disturbance cannot be predicted accurately, angle tracking error e is bound to be generated in the two joint platforms, and the sliding mode surface r is a function of e. Therefore, according to the absolute value of r, the system error and external disturbance can be estimated to adjust the coefficient ${K_r}$ of robust term adaptively. The r as absolute value of the sliding mode surface is taken as the input variable of the fuzzy controller, and the output variable of the fuzzy controller is the variation of the robust term coefficient, and its adaptive adjustment formula is as follows:

(16) \begin{align}{K_{mr}} = K_r^0 + \Delta {K_r}\end{align}

where ${K_{mr}}$ is the robust term coefficient solved by the fuzzy controller, and $K_r^0$ is the initial value of the robust term coefficient.

Figure 9. Hardware block diagram of the flexible manipulator test platform.

Figure 10. Schematic diagram of the test platform.

Figure 11. Experimental curves of joints’ angle.

The switching process of sign function $sgn(r)$ is not continuous and smooth, which is another main reason for the high-frequency chattering phenomenon of the system. The saturation function method is used to redesign the function of the robust term, which takes the hyperbolic tangent function instead of the sign function to make the output smooth and bounded. The concrete expression of hyperbolic tangent function is as follows:

(17) \begin{align}\tanh (r) = \frac{{exp(r) - exp({- r})}}{{exp(r) + exp({- r})}}\end{align}

The concrete expression of the modified fuzzy robust term is

(18) \begin{align}{\tau _r} = {K_{mr}}{\rm{tanh}}(r)\end{align}

In summary, the control block diagram of the fuzzy sliding mode control strategy based on the saturation function is shown in Fig. 4, and its control law is

(19) \begin{align}\tau = {\tau _m} + {K_p}r + {K_i}\int rdt + {K_{mr}}{\rm{tanh}}(r)\end{align}

3.4. Simulation and verification

The manipulator arm is driven by two parallel joints, the scapula joint and the humerus joint. During the process of co-simulation, according to the desired trajectory, two joints are simultaneously controlled to simulate the motion, which can greatly reflect the interaction of the whole flexible manipulator under the actual motion. The structure of the flexible manipulator is shown in Fig. 5 as following:

Setting the desired trajectory of scapula and humerus joints as follows:

(20) \begin{align}\alpha = 10 \cdot \sin \left( {{{2\pi } \over 5}t} \right),\;\beta = 10 \cdot \sin \left( {{{2\pi } \over 5}t} \right),\;\gamma = 10 \cdot \sin \left( {{{2\pi } \over 5}t} \right)\;\;\;\;{\rm{\;\;\;}}[{{\rm{deg}}}]\end{align}

In the simulation of flexible manipulator, setting motion as 10 $\;\cdot \sin \left( {2\pi t/5} \right)$ for each direction of shoulder blade joint and humerus joints, then we analyze the simulation results, respectively, under the nominal model PI robust sliding mode control and fuzzy sliding mode control based on saturated function.

Table III. Experimental error range of scapula joint position tracking.

Table IV. Experimental error range of humerus joint position tracking.

Figure 12. Control force of each PMA of scapula joint.

Figure 13. Control force of each PMA of humerus joint.

Figure 14. Experimental images of circle trajectory tracking.

For scapula joint, the initial parameters of the two sliding mode controls are as follows:

\begin{align*}{K_{p1}} = diag[{3750,3750,3750}],\;{K_{i1}} = diag[{35,35,35}],\;{K_{r1}} = diag[{30,30,30}].\end{align*}

For humerus joint, the initial parameters of the two sliding mode controls are

\begin{align*}{K_{p2}} = diag[{3200,3200,3200}],\;{K_{i2}} = diag[{20,20,20}],\;{K_{r2}} = diag[{15,15,15}].\end{align*}

Through the cosimulation by MATLAB and ADAMS, position tracking simulation graphics of the scapula and humerus joints, respectively, under two sliding mode control strategies are shown as Figs. 6 and 7.

During the overall process of simulation, the robust parameters of PI robust sliding mode control are fixed. However, the interference of the two joints’ movement can make the system parameters change with the status of movement, and the changes are unknown that have huge influence on the control process probably. The specific control effect is shown in Table I and Table II. The tracking error range of fuzzy sliding mode was reduced. It has better dynamic quality such as high control accuracy and well adaptability.

4.1. Design of flexible manipulator test platform

In this paper, the flexible manipulator is a series-parallel mechanical structure driven by multiple pneumatic muscles, which is formed by a parallel platform of scapula joint and a parallel platform of humerus joint in series, and is jointly driven by 10 pneumatic muscles in different positions. The actual flexible manipulator is shown in Fig. 8.

The scapula joint is composed of two pneumatic muscles titled DMSP-10N-110-RM-CM and two pneumatic muscles titled DMSP-10N- 55-RM-CM. The humerus joint is driven by six pneumatic muscles titled DMSP-10N-110-RM-CM. The hardware design of the flexible manipulator test platform is mainly composed of industrial PC, filter, relief valve, proportional valve, air-pressure pump, PCI D/A controller, gyroscope, and so on. The hardware block diagram of the flexible manipulator test platform and the physical map of the test platform are shown in Figs. 9 and 10, respectively:

According to the test requirements of the flexible manipulator control system, the LABVIEW library functions provided by PCI8201 and gyroscope were called in the NI LabVIEW2013 environment, and the MATLAB script was used for complex numerical calculation. According to the above two control strategies, we controlled the movement of the executive end of flexible manipulator in order to carry out arc trajectory tracking control experiment. First, an appropriate projection surface is selected, and then the executive end of the flexible manipulator moves in an arc on the projection surface. Finally, the expected arc trajectory of the executive end is set as following:

(17) \begin{align}\left\{ {\begin{array}{*{20}{c}}{px = p{x_0}}\\[4pt]{py = p{y_0} + R\textrm{cos}({wt})}\\[4pt]{pz = p{z_0} + R\textrm{sin}({wt})}\end{array}} \right.\end{align}

In this formula, $p{x_0},p{y_0},\;$ and $p{z_0}$ are the initial positions of executive end of the flexible manipulator at the terminate frame. The radius of the arc is 40 mm, and the motion period is 21 s. The control parameters ( ${K_p}$ , ${K_i}$ , ${K_r}$ ) are consistent with the simulation.

Note that there are some errors between the actual flexible manipulator and simulation model such as assembly errors, joint frictions, and system modeling errors that cannot be estimated precisely, meaning that the actual control of the executive end of the flexible manipulator has a little deviation from the ideal condition. The observed curves of joints’ angles had slightly higher errors than simulated curves (Fig. 11).

With respect to scapula joint, tracking error range of PI robust sliding mode control for ${\alpha _1}$ , ${\beta _1}$ , and ${\gamma _1}$ were (−0.927, 0.289), (−0.444, 0.310) and (−0.210, 0.308) respectively. By contrast, tracking error range of fuzzy sliding mode control are $\left( { - 0.364,0.333} \right)$ , ${\rm{\;}}\left( { - 0.273,0.227} \right)\!,$ and $\left( { - 0.065,0.176} \right)\!,$ respectively. The tracking error range was reduced by 42.9%, 33.7%, and 53.5%, respectively, as shown in Table III. With respect to humerus joint, tracking error range of PI robust sliding mode control for ${\alpha _2}$ , ${\beta _2},$ and ${\rm{\;}}{\gamma _2}$ were $\left( { - 0.182,0.280} \right)$ , $\left( { - 0.477,0.564} \right)\!,$ ${\rm{and}}$ ${\rm{\;}}\left( { - 0.745,0.667} \right){\rm{\;}}$ respectively. By contrast, tracking error range of fuzzy sliding mode control were $(\!-\!0.106,0.124),\;(\!-\!0.378,0.377),$ and $\left( { - 0.431,0.263} \right)$ respectively. The tracking error range was reduced by 50.2%, 27.5%, and 50.8%, respectively, as shown in Table IV. These data indicate that the trajectory of joint angle under fuzzy sliding mode control is closer to the expected curve, and the tracking error is smaller than that under PI robust sliding mode control. Fuzzy sliding mode control can better realize the joint trajectory tracking control of the flexible manipulator and has higher accuracy.

During the experiment, the control force under two sliding model control strategies was recorded. Then we got the experimental curves of the control force of the pneumatic muscles of each root of the scapula and humerus joints. As can be seen from Figs. 12 and 13, the control force of each root pneumatic artificial muscle under fuzzy sliding mode control was relatively smooth, and the high-frequency shaking was effectively controlled.

Table V. Errors of circle trajectory tracking.

4.2. Deviation analysis of the trajectory curve

The spatial trajectory of the executive end of the manipulator is processed by the angle data collected by the sensor, and the tracking curve of the terminal trajectory is obtained. Fig. 14(a) shows the motion curve of the executive end of the flexible manipulator, Fig. 14(b) shows the curve of the deviation distance from the projection plane against time, Fig. 14(c) shows the motion curve of the spatial motion curve projected into the projection plane, and Fig. 14(d) shows the deviation error curve of the motion curve on the projection plane. The error calculation formula is shown as following:

(18) \begin{align}{\varepsilon _e} = \left| {\sqrt {p{y^2} + p{z^2}} - R} \right|/R\end{align}

According to Fig. 14(a) and (c), the flexible manipulator can make a planned arc motion with a certain accuracy. Compared with PI robust sliding mode control, the tracking error of fuzzy sliding mode control is further reduced, and its trajectory is closer to the expected arc, which satisfied the requirements of control better. According to Fig. 14(b), the executive end cannot move gently and stably on the projection plane, and there is a certain amount of deviation because of the strong coupling of the flexible manipulator. In the process of motion, the deviation range of PI robust sliding mode control from the projection plane was (−1.690 mm, 1.217 mm), with an average deviation of 0.276 mm; while the deviation range of fuzzy sliding mode control from the projection plane was (−1.023 mm, 0.804 mm), with an average deviation of 0.208 mm, with a reduction of 37% and 25%, respectively. According to Fig. 14(d), the maximum error of fuzzy sliding mode control was 5.266%, and the average error was 2.713%. Compared with the maximum error of 9.748% and the average error of 5.266% in PI robust sliding mode control, the errors of fuzzy sliding mode control were reduced by 35% and 48%, respectively. The specific tracking error is shown in Table V.