Design and single-parameter adaptive fuzzy control of pneumatic lower limb exoskeleton with full state constraints

Abstract

With the excellent characteristic of intrinsic compliance, pneumatic artificial muscle can improve the interaction comfort of wearable robotic devices. This paper resolves the safety tracking control problem of a pneumatically actuated lower limb exoskeleton system. A single-parameter adaptive fuzzy control strategy is proposed with high control precision and full state constraints for the safe gait training tasks. Based on the barrier Lyapunov function, all signals in the closed-loop system can be bounded in finite time, which guarantees the deviation of the exoskeleton’s moving trajectory within a bounded range. Furthermore, with the proposed single-parameter adaptive law, the computational burden and the complexity of the control are reduced significantly. Finally, numerical simulations, no-load tracking experiments, and passive and active gait training experiments with healthy subjects validate the effectiveness of the proposed method.

1. Introduction

The pneumatic artificial muscle (PAM) has been extensively investigated in the application of rehabilitation robotics for its unique characteristics of high power-to-weight ratio and intrinsically compliant movement in human-robot interaction tasks [1]. Numerical research on PAM-driven exoskeletons has emerged in the past several years. Ferris et al. designed a pneumatically actuated ankle-foot orthosis for the assistance of ankle flexion [2]. He et al. proposed a pneumatic upper limb rehabilitation robot with four degrees of freedom (DOF), and each DOF is powered by the antagonistic PAMs [3]. Recently, Cao et al. designed a lower limb exoskeleton with the antagonistic PAMs actuating both the hip and knee joints to assist passive gait training tasks and proposed a single-layer learning control strategy to achieve high precision tracking control [4, 5]. Since the PAM only provides unidirectional actuation force, antagonistic configurations are commonly utilized in the existing PAM-driven exoskeleton design. Furthermore, in most previous research, the PAM actuators are positioned parallel to the human limb to drive the joint rotation, which increases the weight of these pneumatic exoskeletons and also makes the high precision control a challenging task. In this paper, a pneumatically actuated lower limb exoskeleton is newly designed. The actuation force of the PAM is transmitted through the Bowden cable to separate the PAM actuator module from the wearable exoskeleton module. With an elastic spring installed at the exoskeleton generating antagonistic directional tension, the exoskeleton can drive the movement of a hip joint with only one PAM. The proposed exoskeleton design significantly reduces the weight of the wearable parts of the exoskeleton system while maintaining the advantage of the PAM’s intrinsic compliance.

Despite the advantages mentioned above, the PAM actuator possesses highly nonlinear dynamics with hysteresis and time-varying model parameters, which makes the control of PAM-driven systems a challenging problem. Various control strategies have been investigated to achieve satisfactory control performances in the past two decades, such as PID-based control [6, 7], proxy-based sliding mode control [8], fuzzy control [9], and neural network-based control [10]. In ref. [11], Cao et al. proposed a predictive control strategy based on the echo state network and successfully reduced the relative error to $3.0\%$ . In ref. [12], Zhu et al. proposed an adaptive control strategy to tackle the high-frequency motion control problem of the PAM and achieved a $10\%$ relative error at $2\;Hz$ motion frequency. With the techniques of the neural network, robust control, adaptive control, etc., the high tracking control precision of PAM-driven systems has been achieved.

However, it is worth noting that few studies have paid attention to the safety of such PAM-driven systems theoretically. During the last decade, researchers have considered introducing the modeling and control of the PAM’s compliance into the trajectory tracking control to enhance the mechanism safety of PAM-driven systems. In ref. [13], a sliding mode control was proposed to track the desired trajectory and simultaneously vary the average pressure of the antagonistic PAMs. The compliance of the mechanism can be adjusted by varying the average pressure inside the PAMs. In ref. [14], Cao et al. proposed a multi-input-multi-output sliding mode controller to simultaneously control angular trajectory and compliance of a pneumatically actuated knee joint rehabilitation exoskeleton. In ref. [15], Choi et al. proposed a control strategy for a pneumatic manipulator to control the joint compliance and associated positions independently, which enhanced human safety during collisions. These studies could achieve accurate control of the position and compliance of the PAM-driven exoskeleton, laying the foundation for further assist-as-needed training. However, as pointed out in ref. [14], what kind of compliance or average pressure can ensure the exoskeleton is safe is unknown, and no reference compliance or average pressure could be provided.

Actually, only the desired trajectory is known as the ideal signal in the exoskeleton control. Therefore, constraining the states of the exoskeleton system to keep the tracking trajectories within a predefined range is a possible way to theoretically achieve the safety control of the exoskeleton system. If the tracking error $e$ between the actual movement trajectory of the exoskeleton and the desired trajectory can be ensured within a prescribed range, the overshoot of the actual trajectory can be limited, and thus, the safety of the exoskeleton system can be guaranteed. If the derivative of the tracking error $\dot{e}$ can be restricted to a specified range, then the oscillations of the exoskeleton can be limited to a certain extent, thus improving the safety of the exoskeleton system. In ref. [16], Li et al. proposed an adaptive finite time control strategy for the nonlinear strict-feedback system with full state constrained and dead-zone. In ref. [17], He et al. proposed an adaptive neural network control strategy based on the barrier Lyapunov function and successfully prevented the violation of the full state constraints in the control of the uncertain n-link robot. In the past several years, the barrier Lyapunov function has emerged as a powerful method to constrain the system states theoretically. Nevertheless, adaptive control based on fuzzy systems or neural networks, which are essential parts of these methods, usually requires updating the full weights of the neural network or the full consequent parameters of the fuzzy system. When the neurons or fuzzy rules are increased, the adaptive parameters are also increased, which inevitably increases the complexity of the algorithm. Therefore, this paper proposes a novel single-parameter adaptive fuzzy control (SAFC) strategy based on the barrier Lyapunov function and investigates the safe control of the PAM-driven lower limb exoskeleton. With the proposed SAFC method, the exoskeleton system states could be constrained in the predefined ranges so that the exoskeleton’s movement would not largely deviate from the desired trajectory, and thus, the safe control of the exoskeleton can be achieved. Furthermore, in the proposed SAFC, a novel single-parameter adaptive law is proposed in the fuzzy system design, which can effectively reduce the complexity of the control design. Numerical simulations and practical experiments with no-load validate that the proposed SAFC method can effectively achieve the tracking control of the designed pneumatic lower limb exoskeleton. In the control process, all the signals of the exoskeleton system are guaranteed to be bounded in a predefined range, which ensures the safety of the exoskeleton’s movement. In addition, both passive and active gait training experiments are conducted to verify the feasibility of the proposed exoskeleton system in rehabilitation.

In summary, the contributions of this paper can be concluded as follows:

1. (1) A lightweight pneumatically actuated lower limb exoskeleton with a Bowden cable transmitting the actuation force is newly designed.

2. (2) For the safety tracking control of the exoskeleton system, a novel SAFC method considering full state constraints is proposed. With only one single adaptation parameter, the algorithm complexity of the controller design is significantly reduced.

3. (3) Based on the barrier Lyapunov function, the stability of the closed-loop system is theoretically guaranteed. Numerical simulations and practical gait training experiments validate the effectiveness of the proposed SAFC method and the feasibility of the exoskeleton system.

The rest of this paper is organized as follows. The development of the lower limb exoskeleton and the dynamic modeling are presented in Section 2. The controller design is given in Section 3. In Sections 4 and 5, numerical simulations and experimental results are provided. Section 6 draws the conclusions.

2. Development of the lower limb exoskeleton system

A pneumatically actuated lower limb exoskeleton assisting the movement of the hip joint is proposed in this paper. By separating the actuating part from the exoskeleton, the weight of the entire wearable exoskeleton is reduced to about 5 kg. Using the PAMs as actuators and Bowden cables to transmit the driving force, the exoskeleton weight is reduced while the compliance of the movement is ensured. In this section, the mechanical design of the exoskeleton is described in detail. Then, the dynamic model of the exoskeleton is extracted.

2.1. Mechanical design

The mechanical design of the proposed 1-DOF hip exoskeleton is represented in Fig. 1. The main components of the exoskeleton include the Bowden cable, the elastic spring, the thigh linkage, the thigh support part, the torso support part, the sliding block, the angle sensor, and the force sensor. The pulling force $F_P$ generated by the PAMs is transmitted through the Bowden cable to drive the thigh linkage to rotate counterclockwise, and the force $F_T$ generated by the spring actuates the linkage to rotate clockwise. The Bowden cable transmission allows remote placement of the PAM actuators so as to minimize gravitational and inertial loading on the user. The sliding block connected with the torso support component and the thigh support components can be adjusted to match the height of the subject before the exoskeleton is put on, which enables the rotation axis of the robot to be aligned with the axis of the human hip joint. Furthermore, there is a sliding block at the back of the torso support component so that the width of the exoskeleton can be adjusted to match the width of the subject’s trunk.

Figure 1. The mechanical design of the lower limb exoskeleton.

The side view of the exoskeleton is represented in Fig. 2. When the exoskeleton is at the flexion state shown in Fig. 2(a), the elastic spring is at the original length, and the thigh linkage is at the end of the swing phase. When the exoskeleton is at the extension state shown in Fig. 2(b), the thigh linkage is at the beginning of the swing phase and also the ending of the standing phase.

Figure 2. The side view of the exoskeleton.

When the pneumatic actuating part acts on the thigh linkage through the Bowden cable, the thigh linkage of the exoskeleton rotates counterclockwise, and the spring is stretched to store the elastic potential energy. When the pneumatic muscle is deflated, the force in the Bowden cable is released. Thus, the spring releases its elastic potential energy to actuate the thigh linkage rotates clockwise and returns to the initial position.

2.2. Control architecture

As represented in Fig. 3, an xPC control system is constructed in this paper to realize the real-time control of the proposed lower limb exoskeleton. The control architecture consists of the following three parts: (1) the compressed air part (yellow); (2) the driving force part (green); (3) the circuit part (blue). In the compressed air part, the air compressor (Fusheng, SA08+) provides the compressed air to the proportional valve (Festo, VPPM-6L-L-1-G18-0L10H-V1N), and the proportional valve gives the compressed air with the desired air pressure to the PAMs (Festo, DMSP-40-600N-RM-RM). The air pressure of the valve is proportional to the control signal. In the driving force part, the PAMs actuated by the compressed air generate a pulling force transmitted by the Bowden cable to actuate the exoskeleton hip joints’ movement. In the control circuit part, the control strategy is programmed in the host PC with the Simulink module of Matalb2013b, and the xPC target (Advantech, IPC-610L) downloads the control strategy in the form of machine codes from the host PC. Then, the xPC target equipped with the data acquisition card (NI, PCI-6229) receives the real-time data of the angle sensors (Rep-avago,18S-1024-2MD-2-15-00E) and calculates the control signal. The control signal is then sent to the proportional valve, and the inside air pressure of the PAMs can be controlled. Thus, the real-time control of the exoskeleton’s movement is realized.

Figure 3. The real-time control architecture of the lower limb exoskeleton system.

2.3. Dynamic model

To facilitate the dynamic analysis, we construct the schematic diagram of the exoskeleton represented in Fig. 4. The motion of the thigh linkage is equated as a rigid, uniformly distributed mass connecting linkage BC rotating around the fixed axis O.

Figure 4. The schematic diagram of the exoskeleton dynamics.

The black dotted line segments AB, BC, and DE represent the positions of the elastic spring, the thigh rod, and the Bowden cable at the initial flexion state, respectively. The point O is the axis of rotation of the thigh linkage. $l_T$ denotes the length of OB, and $m$ represents the mass of the whole thigh linkage. The dotted line OZ denotes the reference line perpendicular to the ground. $\theta _0$ denotes the angle between the thigh linkage and OZ at the initial state, and $\theta$ indicates the rotation angle of the thigh linkage relative to the initial state. According to Lagrange’s theorem, the dynamics of this single-link model can be expressed as:

(1) \begin{equation} \tau _{hip}=m\ddot{\theta }(l_P+l_T)^2+I\ddot{\theta }+\Delta G, \end{equation}

where the first term is related to the linear velocity, the second term is associated with the angular velocity, and the third term is related to the gravitational potential energy. The effect of friction on the joint motion would be regarded as modeling uncertainties later.

The mass of a unit length of the linkage can be expressed as:

(2) \begin{equation} \lambda =\frac{m}{l_P+l_T}. \end{equation}

The moment of inertia $I$ can be expressed as:

(3) \begin{equation} \begin{aligned} I&=\int _{-l_P}^{l_T}\lambda x^2dx\\[5pt] &=\frac{m}{3(l_P+l_T)}(l_P^3+l_T^3) \end{aligned} \end{equation}

The gravitational potential energy term $\Delta G$ can be expressed as:

(4) \begin{equation} \Delta G=[\cos{\theta _0}-\cos{(\theta -\theta _0)}]\frac{mgl_P(l_P-l_T)}{2(l_P+l_T)}, \end{equation}

where $m_P$ and $m_T$ denote the mass of the OB part and the mass of the OC part with the axis O as the dividing line. Thus, $m_P$ and $m_T$ can be expressed as:

(5) \begin{equation} m_P=\frac{ml_P}{l_P+l_T}, \end{equation}

(6) \begin{equation} m_T=\frac{ml_T}{l_P+l_T}. \end{equation}

The simplified single-link dynamic model of the exoskeleton can be expressed as:

(7) \begin{equation} \tau _{hip}=m\ddot{\theta }(l_P+l_T)^2+\frac{m\ddot{\theta }}{3(l_P+l_T)}(l_P^3+l_T^3)+[\cos{\theta _0}-\cos{(\theta -\theta _0)}]\frac{mgl_P(l_P-l_T)}{2(l_P+l_T)}. \end{equation}

The moments produced by the Bowden cable pulling force and the elastic spring force on the joint are analyzed below. AB’, B’C’, and DE’ represent the current positions of the elastic spring, the exoskeleton thigh rod, and the Bowden cable when the PAMs are inflated and the Bowden cable actuates the exoskeleton to rotate at an angle $\theta$ . $F_P$ and $F_T$ in the red part of the figure represent the force transmitted by the Bowden cable and the force generated by the spring stretching, respectively. The red dashed lines represent the force arm $d_T$ of the force $F_P$ and the force arm $d_P$ of the spring force $F_T$ , respectively. It can be found that the force arm $d_T$ of the spring hardly changes during the rotation, so it is assumed that the force arm $d_T$ of the spring is constant and $d_T=l_T$ . The approximate length of the spring stretch is expressed as:

(8) \begin{equation} \Delta x=\theta l_T. \end{equation}

Then, the moment generated by the force of the spring can be expressed as:

(9) \begin{equation} \begin{aligned} \tau _T&= \Delta x k_T d_T\\[5pt] &=\theta k_T l_T^2. \end{aligned} \end{equation}

For the force arm $d_P$ of the Bowden cable, approximated here as $\textsf{OD}=\textsf{OE}$ , the force arm of the Bowden cable can be calculated as:

(10) \begin{equation} d_P=\frac{l_P}{2}\cos{\left(\frac{\theta _1-\theta }{2}\right)}, \end{equation}

where $\theta _1$ is the angle between OC and OD at the initial state.

The tension in the Bowden cable $F_P$ can be obtained through the three-element PAM model which has been extensively used in the previous research [18, 19]. Therefore, $F_P$ can be expressed as:

(11) \begin{equation} F_P=F(P)-B(P)\Delta \dot{x}-K(P)\Delta x, \end{equation}

where $P$ is the internal air pressure of the PAM, $B(P)$ denotes the damping coefficient, $K(P)$ denotes the elasticity coefficient, and the values of $B(P)$ and $K(P)$ are related to the states of the PAM. $\Delta x$ denotes the displacement of the PAM, which is also the contraction length of the Bowden cable. Therefore, the relationship between $\Delta x$ and the rotation angle $\theta$ can be expressed as:

(12) \begin{equation} \Delta x=l_P\left(\sin{\frac{\theta _1}{2}} -\sin{\frac{\theta _1-\theta }{2}}\right), \end{equation}

and:

(13) \begin{equation} \Delta \dot{x}= \frac{1}{2}\dot{\theta }\cos{\frac{(\theta _1-\theta )}{2}}. \end{equation}

Substituting (12) and (13) into (11), the pulling force in the Bowden cable can be expressed as:

(14) \begin{equation} F_P=F(P)-B(P)\frac{1}{2}\dot{\theta }\cos{\frac{(\theta _1-\theta )}{2}}-K(P)l_P\left(\sin{\frac{\theta _1}{2}} -\sin{\frac{\theta _1-\theta }{2}}\right). \end{equation}

Therefore, the torque generated by the Bowden cable can be expressed as:

(15) \begin{equation} \begin{aligned} \tau _P&=F_Pd_P\\[5pt] &=\frac{l_P}{2}\cos\left(\frac{\theta _1-\theta }{2}\right)\left[F(P)-B(P)\frac{1}{2}\dot{\theta }\cos{\frac{(\theta _1-\theta )}{2}}-K(P)l_P\left(\sin{\frac{\theta _1}{2}} -\sin{\frac{\theta _1-\theta }{2}}\right)\right]. \end{aligned} \end{equation}

The total torque generated by the Bowden cable and the antagonistic elastic spring can be expressed as:

(16) \begin{equation} \tau _{hip}=\tau _P-\tau _T. \end{equation}

Substituting (7), (9) and (15) into (16) yields:

(17) \begin{equation} M\dot{x_2}=H_1(x_1,x_2)+H_2(x_1,x_2)u, \end{equation}

where

(18) \begin{equation} M=\left(\frac{4}{3}l_P^2+\frac{4}{3}l_T^2+\frac{5}{3}l_Pl_T\right)m, \end{equation}

(19) \begin{equation} \begin{aligned} H_1(x_1,x_2)&=\frac{l_P}{2}\cos{\frac{\theta _1-x_1}{2}}\left[F_0-\frac{B_0x_2}{2}\cos{\frac{\theta _1-x_1}{2}}-K_0l_P\left(\sin{\frac{\theta _1}{2}}-\sin{\frac{\theta _1-x_1}{2}}\right)\right]\\[5pt] &-\frac{mgl_P(l_P-l_T)}{2(l_P+l_T)}\left[\cos{\theta _0}-\cos{(x_1-\theta _0)}\right]-x_1k_Tl_T^2, \end{aligned} \end{equation}

(20) \begin{equation} H_2(x_1,x_2)=\frac{l_P}{2}\cos{\frac{\theta _1-x_1}{2}}\left[F_1-B_1\frac{1}{2}x_2\cos{\frac{\theta _1-x_1}{2}}-K_1l_P\left(\sin{\frac{\theta _1}{2}}-\sin{\frac{\theta _1-x_1}{2}}\right)\right]. \end{equation}

In summary, the dynamic model of the exoskeleton can be expressed as follows:

(21) \begin{equation} \left \{ \begin{aligned} \dot{x}_1 &=x_2\\[5pt] \dot{x}_2 &=f(\bar{x})+g(\bar{x})u(t)+d(t)\\[5pt] y &=x_1, \end{aligned} \right. \end{equation}

(22) \begin{equation} g(\bar{x})=H_1(x_1,x_2)/M, \end{equation}

(23) \begin{equation} f(\bar{x})=H_2(x_1,x_2)/M, \end{equation}

where $y$ denotes the output of the system, $x_1$ and $x_2$ are system states and $\bar{x}=[x_1,x_2]^T$ . Due to the different sizes of the subjects, the human-robot interaction, and the errors of the system parameter identification, there exist errors in the above constructed dynamic model, and $d(t)$ denotes a lumped term including the dynamic modeling uncertainties of the exoskeleton system and the external disturbances.

3. Controller design

The target of the trajectory tracking control is to make the joints of the exoskeleton move with a desired trajectory. The desired trajectory of the joint movement is obtained from the angle data of the hip joint during normal human walking. If the tracking error is too large, that is, the exoskeleton moves at an angle that exceeds the proper angle of movement of the hip joint for normal human walking, then there is a risk of injury to the human joint. If the tracking error between the actual movement trajectory of the exoskeleton and the desired trajectory can be ensured within a prescribed range, the safety of the exoskeleton system can be guaranteed. Therefore, the target of the controller design is to make the output of the system track the desired trajectory $y_r$ guarantee that all the states of exoskeleton system (21) are constrained in predefined compact sets, that is, $|x_i|\lt k_{ci}$ , where $k_{ci}\;(i=1,2)$ are predefined positive constants in finite time. Furthermore, the lumped term $d(t)$ is approximated and compensated by a designed adaptive fuzzy system with only one adaptation parameter.

3.1. Preliminaries

Lemma 3.1. [20] For a continuous unknown function $f(x)$ defined on a compact set $\Omega$ , there exists a fuzzy system $\Theta ^T\Phi (x)$ such that for any given constant $\varepsilon \gt 0$ , the following inequality holds:

(24) \begin{equation} \mathop{sup}\limits _{x \in \Omega }|f(x)-\Theta ^T\Phi (x)|\leq \varepsilon. \end{equation}

Lemma 3.2. [21] For $z_k\in \mathbb{R}$ , $k=1,2,\ldots,n$ , $0\lt l\lt 1$ , the following inequality holds:

(25) \begin{equation} \left(\sum _{k=1}^n|z_k|\right)^l\leq \sum _{k=1}^n|z_k|^l\leq n^{1-l}\left(\sum _{k=1}^n|z_k|\right)^l. \end{equation}

Lemma 3.3. [22, 23] Considering the system $\dot{\varsigma }=f(\varsigma,u)$ , for smooth positive definite function $V(\varsigma )$ , if there exist scalar $c\gt 0$ , $0\lt l\lt 1$ and $\varrho \gt 0$ , such that:

(26) \begin{equation} \dot{V}(\varsigma )\leq -cV^l(\varsigma )+\varrho,\;\;\;t\geq 0, \end{equation}

the nonlinear system $\dot{\varsigma }=f(\varsigma,u)$ is semi-global practical finite time stable (SGPFS).

3.2. The single-parameter adaptive fuzzy controller (SAFC) design

According to Lemma 3.1, the lumped term $d(t)$ can be expressed in the form of a fuzzy system as follows:

(27) \begin{equation} d(t)=\Theta ^{*T}\Phi (\bar{x})+\epsilon, \end{equation}

where $*$ means the ideal fuzzy consequence parameter value, $\bar{x}$ denotes the vector containing the system states, and $\epsilon$ is the approximation error.

The fuzzy system approximating the actual $d(t)$ is expressed as:

(28) \begin{equation} \hat{d}=\hat{\Theta }^T\Phi (\bar{x}), \end{equation}

where $\hat{}$ denotes the approximate fuzzy consequence parameter value.

We define the single adaptation parameter as:

(29) \begin{equation} \hat{\omega }=\hat{\Theta }^T\hat{\Theta }. \end{equation}

The approximate error can be expressed as:

(30) \begin{equation} \begin{aligned} \tilde{d}(t)&=d(t)-\hat{d}\\[5pt] &=\Theta ^{*T}\Phi (\bar{x})-\hat{\Theta }^T\Phi (\bar{x})+\epsilon \\[5pt] &=\tilde{\Theta }^T\Phi (\bar{x})+\epsilon. \end{aligned} \end{equation}

The main results of the controller are summarized in the following theorem.

Theorem 3.4. For the nonlinear system (21) with full state constraints, the lumped term $d(t)$ can be approximated and compensated by a fuzzy system with one adaptation parameter $\hat{\omega }$ , all signals of the close-loop system are SGPFS in a finite time with the following virtual control signal $\alpha _1$ , controller $u(t)$ and adaptive law $\dot{\hat{\omega }}$ :

(31) \begin{equation} \alpha _1=\dot{y}_r-\frac{e_1^{2l-1}}{(k_{b1}^2-e_1^2)^{l-1}}-\frac{e_1}{2(k_{b1}^2-e_1^2)}, \end{equation}

(32) \begin{equation} u(t)=\frac{1}{g(\bar{x})}[-f(\bar{x})+\dot{\alpha }_1-\frac{e_2^{2l-1}}{(k_{b2}^2-e_2^2)^{l-1}}-\frac{e_2}{2(k_{b2}^2-e_2^2)}-\frac{\hat{\omega }\Phi (\bar{x})^T\Phi (\bar{x})e_2}{2l_1^2(k_{b2}^2-e_2^2)}-\frac{e_2}{2}(k_{b2}^2-e_2^2)], \end{equation}

(33) \begin{equation} \dot{\hat{\omega }}=\gamma [\frac{\Phi (\bar{x})^T\Phi (\bar{x})e_2^2}{2l_1^2(k_{b2}^2-e_2^2)^2}-l_2\hat{\omega }], \end{equation}

where $0\lt l\lt 1$ , $l_1\gt 0$ , $l_2\gt 0$ are parameters to be designed, and $k_{b1}\gt 0$ and $k_{b2}\gt 0$ are constraints of the errors $e_1$ and $e_2$ .

Proof. We define the error terms as follows:

(34) \begin{equation} e_1=x_1-y_r, \end{equation}

(35) \begin{equation} e_2=x_2-\alpha _1, \end{equation}

where $\alpha _1$ denotes the virtual control signal of $x_2$ .

Taking the derivative of $e_1$ yields:

(36) \begin{equation} \begin{aligned} \dot{e}_1&=x_2-\dot{y}_r. \end{aligned} \end{equation}

The first Lyapunov function with a barrier function disposing of the constraint of $e_1$ is designed as:

(37) \begin{equation} V_1=\frac{1}{2}\log{\frac{k^2_{b1}}{k^2_{b1}-e_1^2}}, \end{equation}

where $\log$ denotes the nature logarithm function, $k_{b1}\gt 0$ denotes the constraint of $e_1$ .

Taking the derivative of $V_1$ yields:

(38) \begin{equation} \begin{aligned} \dot{V}_1&=\frac{e_1(\dot{y}_r-\dot{x}_1)}{k_{b1}^2-e_1^2}\\[5pt] &=\frac{e_1(\dot{y}_r-x_2)}{k_{b1}^2-e_1^2}. \end{aligned} \end{equation}

Substituting (35) into (38) yields:

(39) \begin{equation} \dot{V}_1=\frac{e_1(\dot{y}_r-e_2-\alpha _1)}{k_{b1}^2-e_1^2}. \end{equation}

Substituting the virtual signal $\alpha _1$ into (39) yields:

(40) \begin{equation} \begin{aligned} \dot{V}_1=\frac{e_1e_2}{k_{b1}^2-e_1^2}-\left[\frac{e_1^2}{k_{b1}^2-e_1^2}\right]^l-\frac{e_1^2}{2(k_{b1}^2-e_1^2)^2}. \end{aligned} \end{equation}

According to Young’s inequality, we have:

(41) \begin{equation} \frac{e_1e_2}{k_{b1}^2-e_1^2}\leq \frac{e_1^2}{2(k_{b1}^2-e_1^2)^2}+\frac{1}{2}e_2^2. \end{equation}

Substituting (41) into (40), we get:

(42) \begin{equation} \dot{V}_1\leq -\left[\frac{e_1^2}{k_{b1}^2-e_1^2}\right]^l+ \frac{1}{2}e_2^2. \end{equation}

It is obvious that $\log \frac{k_{b1}^2}{k_{b1}^2-e_1^2}\leq \frac{e_1^2}{k_{b1}^2-e_1^2}$ when $|e_1|\lt k_{b1}$ , which proves that $[\log \frac{k_{b1}^2}{k_{b1}^2-e_1^2}]^l\leq [\frac{e_1^2}{k_{b1}^2-e_1^2}]^l$ holds for $0\lt l\lt 1$ . Then, we can obtain:

(43) \begin{equation} \dot{V}_1\leq -V_1^l+\frac{1}{2}e_2^2. \end{equation}

We design another Lyapunov function as:

(44) \begin{equation} V_2=\log{\frac{k_{b2}^2}{k_{b2}^2-e_2^2}}+\frac{1}{2\gamma }\tilde{\omega }^T\tilde{\omega }, \end{equation}

where $\omega ^*=\Theta ^{*T}\Theta ^*$ , $\tilde{\omega }=\omega ^*-\hat{\omega }$ , $\gamma \gt 0$ , and $k_{b2}\gt 0$ denotes the constraint of $e_2$ .

Taking the derivative of $V_2$ and substituting $\dot{x}_2$ into the result yield:

(45) \begin{equation} \begin{aligned} \dot{V}_2&=\frac{e_2\dot{e}_2}{k_{b2}^2-e_2^2}-\frac{1}{\gamma }\tilde{\omega }\dot{\hat{\omega }}\\[5pt] &=\frac{e_2(\dot{x}_2-\dot{\alpha }_1)}{k_{b2}^2-e_2^2}-\frac{1}{\gamma }\tilde{\omega }\dot{\hat{\omega }}\\[5pt] &=\frac{e_2(f(\bar{x})+g(\bar{x})u(t)+d(t)-\dot{\alpha }_1)}{k_{b2}^2-e_2^2}-\frac{1}{\gamma }\tilde{\omega }\dot{\hat{\omega }}. \end{aligned} \end{equation}

Substituting the control signal (32) and (27) into (45) yields:

(46) \begin{equation} \begin{aligned} \dot{V}_2&=\frac{e_2}{k_{b2}^2-e_2^2}\Bigg\{f(\bar{x})+\Bigg[-f(\bar{x})+\dot{\alpha }_1+\frac{-e_2^{2l-1}}{(k_{b2}^2-e_2^2)^{l-1}}-\frac{e_2}{2(k_{b2}^2-e_2^2)}-\frac{\hat{\omega }\Phi (\bar{x})^T\Phi (\bar{x})e_2}{2l_1^2(k_{b2}^2-e_2^2)}\\[5pt] &-\frac{e_2}{2}(k_{b2}^2-e_2^2)\Bigg]+\Theta ^{*T}\Phi (\bar{x})+\epsilon -\dot{\alpha }_1)\Bigg\}-\frac{1}{\gamma }\tilde{\omega }\dot{\hat{\omega }}\\[5pt] &=\frac{e_2}{k_{b2}^2-e_2^2}\Bigg\{\frac{-e_2^{2l-1}}{(k_{b2}^2-e_2^2)^{l-1}}-\frac{e_2}{2(k_{b2}^2-e_2^2)}-\frac{\hat{\omega }\Phi (\bar{x})^T\Phi (\bar{x})e_2}{2l_1^2(k_{b2}^2-e_2^2)}-\frac{e_2}{2}(k_{b2}^2-e_2^2)+\Theta ^{*T}\Phi (\bar{x})\\[5pt] &+\epsilon \Bigg\}-\frac{1}{\gamma }\tilde{\omega }\dot{\hat{\omega }}\\[5pt] &=-\left[\frac{e_2^2}{k_{b2}^2-e_2^2}\right]^l-\frac{e_2^2}{2(k_{b2}^2-e_2^2)^2}-\frac{\hat{\omega }\Phi (\bar{x})^T\Phi (\bar{x})e_2^2}{2l_1^2(k_{b2}^2-e_2^2)^2}-\frac{e_2^2}{2}+\frac{e_2\Theta ^{*T}\Phi (\bar{x})+e_2\epsilon }{k_{b2}^2-e_2^2}-\frac{1}{\gamma }\tilde{\omega }\dot{\hat{\omega }}. \end{aligned} \end{equation}

With Young’s inequality, we have:

(47) \begin{equation} \frac{e_2\Theta ^{*T}\Phi (\bar{x})}{k_{b2}^2-e_2^2}\leq \frac{\omega ^*\Phi (\bar{x})^T\Phi (\bar{x})e_2^2}{2l_1^2(k_{b2}^2-e_2^2)^2}+\frac{l_1^2}{2}, \end{equation}

(48) \begin{equation} \frac{e_2\epsilon }{k_{b2}^2-e_2^2}\leq \frac{e_2^2}{2(k_{b2}^2-e_2^2)^2}+\frac{\epsilon ^2}{2}. \end{equation}

Then, we can obtain:

(49) \begin{equation} \begin{aligned} \dot{V}_2\leq &-\left[\frac{e_2^2}{k_{b2}^2-e_2^2}\right]^l+\frac{\tilde{\omega }\Phi (\bar{x})^T\Phi (\bar{x})e_2^2}{2l_1^2(k_{b2}^2-e_2^2)^2}+\frac{l_1^2}{2}+\frac{\epsilon ^2}{2}-\frac{e_2^2}{2}-\frac{1}{\gamma }\tilde{\omega }\dot{\hat{\omega }}. \end{aligned} \end{equation}

Substituting the adaptive law (33) into (49) yields:

(50) \begin{equation} \dot{V}_2\leq -\left[\frac{e_2^2}{k_{b2}^2-e_2^2}\right]^l+\frac{l_1^2}{2}+\frac{\epsilon ^2}{2}-\frac{e_2^2}{2}-l_2\hat{\omega }\tilde{\omega }. \end{equation}

Using the relationship of $\hat{\omega }=\omega ^*-\tilde{\omega }$ and Young’s inequality, we have:

(51) \begin{equation} \tilde{\omega }\hat{\omega }=\tilde{\omega }(\omega ^*-\tilde{\omega })=-\tilde{\omega }^2+\tilde{\omega }\omega ^*\leq -\frac{1}{2}\tilde{\omega }^2+\frac{1}{2}\omega ^{*2}. \end{equation}

Then, we can obtain:

(52) \begin{equation} \begin{aligned} \dot{V}_2&\leq -\left[\frac{e_2^2}{k_{b2}^2-e_2^2}\right]^l+\frac{1}{2}\omega ^{*2}-\frac{1}{2}\tilde{\omega }^2+\frac{1}{2}l_1^2+\frac{1}{2}\epsilon ^2-\frac{1}{2}e_2^2\\[5pt] &\leq -V_2^l+\frac{1}{2}\omega ^{*2}+\frac{1}{2}l_1^2+\frac{1}{2}\epsilon ^2-\frac{1}{2}e_2^2. \end{aligned} \end{equation}

Define the final Lyapunov function as:

(53) \begin{equation} V=V_1+V_2. \end{equation}

Taking the derivative of $V$ and substituting (43) and (52) into it, we get:

(54) \begin{equation} \begin{aligned} \dot{V}&=\dot{V}_1+\dot{V}_2\\[5pt] &\leq -V_1^l-V_2^l+\frac{1}{2}\omega ^{*2}+\frac{1}{2}l_1^2+\frac{1}{2}\epsilon ^2. \end{aligned} \end{equation}

Applying Lemma 3.2, we obtain:

(55) \begin{equation} \begin{aligned} \dot{V}&\leq -(V_1+V_2)^l+\frac{1}{2}\omega ^{*2}+\frac{1}{2}l_1^2+\frac{1}{2}\epsilon ^2\\[5pt] &\leq -V^l+\rho, \end{aligned} \end{equation}

where $\rho =\frac{1}{2}\omega ^{*2}+\frac{1}{2}l_1^2+\frac{1}{2}\epsilon ^2$ is a positive constant.

Define $T^*$ as:

(56) \begin{equation} T^*=\frac{1}{\eta (1-l)}\left[V^{1-l}(0)-\left(\frac{\rho }{(1-\eta )^{\frac{1-l}{l}}}\right)\right], \end{equation}

where $0\lt \eta \lt 1$ , $V(0)$ is the initial value of $V(t)$ . According to Lemma 3.3, for $\forall t\gt T^*$ , $V^l\leq \frac{\rho }{(1-\eta )*1}$ , that is, all the signals of the closed-loop system are SGPFS. Furthermore, based on the definition of $V(t)$ in ref. [16], for $\forall t\gt T^*$ , we have:

(57) \begin{equation} |y-y_r|\leq k_{b1}\left[1-e^{-2\left(\frac{\rho }{(1-\eta )}\right)^{\frac{1}{l}}}\right]^{\frac{1}{2}}, \end{equation}

which implies that after the finite time $T^*$ , the tracking error remain in a small neighborhood of the origin. We assume that $|y_r|\lt B_0$ ,and thus, it is obvious that $|x_1|\leq |e_1|+|y_r|\lt k_{b1}+B_0$ . We define that $k_{b1}=k_{c1}-B_0$ , and then, it has $|x_1|\lt k_{c1}$ . Similarly, with the assumption that $\dot{y}_r\lt B_1$ , we can finally get $|x_2|\lt k_{c2}$ and $k_{c2}=k_{b2}+B_1$ . Consequently, the system states are confined to the predefined regions effectively. This completes the proof.

4. Numerical simulation results

To validate the proposed SAFC method, numerical simulations are conducted based on the developed dynamic model of the exoskeleton (21). The desired trajectory $y_r$ is generated by the recorded human hip motion data which has been successfully applied to the gait training of ref. [4].

The initial states of the exoskeleton model are set as: $x(1) = 1.0\;{\rm deg}$ , $x(2) = 1.0\;{\rm deg/s}$ . The designed parameters of the controller are set as: $k_{b1}=0.5$ , $b_{b2}=50.0$ , $l=0.5$ , $l_1=1.0$ , $l_2=1.0$ , $\omega =10.0$ , $\gamma =10.0$ . We apply a standard adaptive fuzzy controller (AFC) based on the backstepping method [24] to the same exoskeleton dynamic model to further demonstrate the effectiveness of the proposed SAFC. The standard AFC controller design with the backstepping technique is as follows:

(58) \begin{equation} \left \{ \begin{aligned} &z_1=x_1-y_r,\\[5pt] &x_{2d}=\dot{y}_r-k_1e_1,\\[5pt] &z_2=x_2-x_{2d},\\[5pt] &\hat{d}(t)=\hat{\theta }^T\Phi (\bar{x}),\\[5pt] &u(t)=\frac{1}{g}(f-z_1-k_1\dot{z}_1-k_2z_2-\hat{d}(t)+\ddot{y}_r),\\[5pt] &\dot{\hat{\theta }}=\mu (z_2\Phi (\bar{x})-\zeta \hat{\theta }), \end{aligned} \right. \end{equation}

where $x_1$ and $x_2$ denote the system states, $y_r$ denotes the desired trajectory, $x_{2d}$ denotes the virtual control signal of $x_2$ , $\hat{d}(t)$ denotes the output of the fuzzy system, $\dot{\hat{\theta }}$ denotes the adaptive law updating all the fuzzy consequence parameters, and $k_1\gt 0$ , $k_2\gt 0$ , $\mu \gt 0$ , and $\zeta \gt 0$ are designed control parameters. In the simulations, the control parameters of the standard AFC are set as $k_1=17$ , $k_2=19$ , $\mu =1$ , and $\zeta =1$ .

The tracking trajectories of these two controllers are represented in Fig. 5. It can be seen that the actual trajectories of the exoskeleton joint with the control of AFC and SAFC are nearly overlapped with the desired trajectory, which means that the proposed SAFC and the standard AFC can both effectively realize the trajectory tracking task of the lower limb exoskeleton. The tracking errors are represented in Fig. 6. The tracking errors of the proposed SAFC and the AFC are all small enough and are constrained in the predefined bound $k_{b1}$ . Therefore, the exoskeleton system state $x(1)$ would not violate the predefined constraint in movement, which validates the effectiveness of the proposed SAFC method in ensuring safety. In addition, the value of the adaptation parameter $\hat{\omega }$ and the control signal $u(t)$ is respectively represented in Figs. 7 and 8, from which we can see that these signals are ensured bounded.

Figure 5. Tracking trajectory.

Figure 6. Tracking error.

Figure 7. The adaptation parameter.

Figure 8. The control signal.

To verify the effectiveness of only one single adaptation parameter in reducing computational cost, we conduct four sets of numerical simulations to compare the running time of the proposed SAFC and the standard AFC. The number of fuzzy rules in the proposed SAFC and the standard AFC is set as $49$ . In addition, to further show the effectiveness of the proposed SAFC with a single adaptation parameter in reducing computational burden, the standard AFC with $121$ rules is added to the simulations. The simulation times are set as $50$ , $100$ , $200$ , and $300\;{\rm s}$ . The program runs on Matlab2013b, installed on a Win 10 OS equipped with 16 GB RAM and an Intel Core i7-10875H processor. Table I represents the running time results. It can be observed that the proposed SAFC with a single adaptation parameter can effectively reduce the simulation running time. As the number of fuzzy rules increases to $121$ , the adaptive parameters of the standard AFC increase, and thus, the running time of the AFC increases a lot, which validates the advantage of the proposed SAFC method in computational speed.

Table I. Time consuming comparison between the proposed SAFC and the standard AFC (UNIT: S).

5. Experimental results

5.1. Experiment platform setup

The practical exoskeleton experiment platform is represented in Fig. 9. The main components of the system include an air compressor, an xPC target, a host PC, PAM actuators, Bowden cables, a treadmill, a healthy subject wearing the exoskeleton, and a body weight support system compensating for the subject’s weight. In addition, two pairs of airbags are placed on the front and back sides of the exoskeleton leg straps part to detect the subject’s intention to move forward or backward. The mechanism and control system details have been explained in Section 2.

Figure 9. The practical exoskeleton system.

5.2. Trajectory tracking experiments with no-load

First of all, we conduct a set of experiments with no-load to validate the feasibility of the proposed SAFC method in practical implementation. The no-load experiments can protect the subjects from terrible control results. Then, the proposed SAFC method is applied to the exoskeleton control in healthy subjects’ passive and active gait training tasks, and the supplementary video file shows the experimental results. The sampling time of the real-time system is set as $0.0005\;{\rm s}$ . The desired trajectory $y_r$ of the exoskeleton is obtained by fitting the human lower hip motion data with $1.0m/s$ walking speed on the treadmill.

It is worth noting that since the proposed SAFC in this paper and the standard AFC in the above numerical simulation show similar control performances, and the proposed SAFC achieves better computational efficiency, the standard AFC is not considered in the experimental comparison. A conventional PID controller commonly utilized in practical engineering applications is adopted to compare with the SAFC to further investigate its control performance.

The designed parameters of the proposed SAFC are set as: $k_{b1}=5$ , $b_{b2}=300.0$ , $l=0.5$ , $l_1=10.0$ , $l_2=1.0$ , $\gamma =10.0$ . The initial adaptation parameter $\omega _0=10.0$ . The parameters of the PID controller are set as: $P=1050.0$ , $I=305.0$ , $D=10.0$ . For a fair comparison, we have tried to adjust the control parameters so that the controllers can perform their best in experiments.

The tracking trajectories are represented in Fig. 10, and the corresponding tracking errors are shown in Fig. 11. Only the tracking results of the right hip joint are presented here, since the left and right limbs of the exoskeleton conduct symmetrical movements. We can find that the tracking trajectory of the proposed SAFC nearly overlapped with the desired trajectory in Fig. 10, while that of PID deviates a lot from the desired trajectory, especially in the areas of peaks and bottoms of waves. From Fig. 11, we can find that the maximum tracking error of the proposed SAFC is less than $5\;{\rm deg}$ , which is bounded by the predefined boundary $k_{b1}$ , while that of PID is more than $5\;{\rm deg}$ and the maximum error is about $8\;{\rm deg}$ . Obviously, the proposed SAFC performs better than the conventional PID controller in the practical experiments with no-load. Therefore, we can conclude that the proposed SAFC method can precisely track the desired gait trajectory of the pneumatic lower limb exoskeleton to conduct the following gait training tasks, which indicates that this theoretically safe method guaranteeing the full state constraints can ensure training security of the exoskeleton in practical implementation.

Figure 10. Trajectories of experiments with no-load.

Figure 11. Tracking errors of experiments with no-load.

5.3. Passive gait training

From Figs. 10 and 11, we can find that the tracking trajectory of PID variated a lot from the desired trajectory and the maximum tracking error of PID was larger than 5 deg, which indicated that the PID was unable to constraint the tracking error to a predefined range and achieve safety control compared to the proposed SAFC. Therefore, the PID was not adopted for the poor tracking performance in the following gait tracking experiments.

To further validate the effectiveness of the proposed SAFC method in the task of passive gait training, four healthy volunteer subjects are invited to participate in our experiments. The subjects are asked to remain relaxed and walk on the treadmill with as much assistance as possible from the exoskeleton, thus achieving the effect of simulating hemiplegic patients. It is worth noting that here we consider the influence of the subject’s participation on the exoskeleton mechanism and modify the constraint $k_{b1}$ to $8$ , thus allowing a slightly larger tracking error than the no-load experiments for the successful implementation of the gait training experiments. In addition, we use the maximum absolute error(MAE) and the root mean square error (RMSE) in $10\sim 30\;s$ to evaluate the control performance of the proposed method. The RMSE is defined as:

(59) \begin{equation} RMSE=\sqrt{\frac{1}{n}\sum _{t=1}^{n}e_t^2}, \end{equation}

where $e_t$ denotes the tracking error at the $t$ th sampling time point.

Table II. The information of the healthy subjects and the gait training performances.

Figure 12. Trajectories of passive gait training.

The information of these four subjects and the details of the tracking errors are represented in Table II. The tracking results of the exoskeleton in the passive gait training are represented in Figs. 12 and 13. It can be clearly found that the tracking error is large in the first $10\;{\rm s}$ , while the tracking error can be maintained within $5\;{\rm deg}$ in the following $20\;{\rm s}$ , which indicates that after two gait periods of adjustment, the subjects can adapt well to the exoskeleton assistance to complete the passive gait training. Therefore, with the control results represented above, we can also conclude that although the four subjects vary in their heights and weights, the proposed SAFC method can still achieve effective tracking control to conduct the passive gait training task.

Figure 13. Tracking errors of passive gait training.

5.4. Active gait training

For the rehabilitation of partial patients, increasing the patient’s participation in the rehabilitation will be more beneficial than purely passive gait training [25–27]. Therefore, in this section, we conduct an active gait training experiment on one healthy subject to verify the feasibility of the proposed pneumatic lower limb exoskeleton with the SAFC method. The active control diagram is based on a zero-force control strategy, and the corresponding control diagram is represented in Fig. 14. The air bags detect the subject’s movement intention, the admittance model generates the desired trajectory $\Delta x_{ref}$ based on the interaction force $F_{int}$ and the desired force $F_d$ , and the proposed SAFC method is used to control the joint angle of the exoskeleton to achieve active control. Here we assume that the interaction force $F_{int}$ between the air bags and the subject’s leg and the air pressure difference between the air bag pairs are proportional, and the proportionality gain is $\kappa$ . We set $F_d=0$ and $x_{ref}=0$ so that the exoskeleton is completely passive and the subject is completely active. The parameters of the admittance model are set as: $m=1/500$ , $b=0.05$ , and $k=0.01$ . The subject walks on the treadmill at a speed of $1.0\,{\rm m/s}$ .

Figure 14. The block diagram of the active gait training control strategy.

The pressure difference of the air bags is represented in Fig. 15, and the trajectories are represented in Fig. 16. From Fig. 15, we can find that the air pressure difference of the air bags varied periodically, indicating that the subject’s intention of movement can be successfully detected by the air bags positioned in both the front and back of the subject’s legs. From Fig. 16, we can find that the trajectory generated by the active control strategy shown in the blue curve is close in magnitude and frequency to the reference trajectory in the passive gait training mentioned previously, which indicates that the admittance model can successfully generate the desired reference trajectory used in active training. The red curve shows the movement trajectory of the subject wearing the exoskeleton under the active control strategy. It can be found that the exoskeleton with the SAFC method can effectively track the trajectory generated by the admittance model to achieve the active gait training task.

Figure 15. Pressure difference of the air bags.

Figure 16. Trajectories of active gait training.

6. Conclusions

Incorporating the PAM actuator and the Bowden cable, this paper proposes a lightweight pneumatically actuated lower limb exoskeleton. First of all, we construct the dynamic model of the exoskeleton based on the three-element model of the PAM. Then, utilizing the barrier Lyapunov function, we propose a SAFC for the safe control of the pneumatic system considering full state constraints. Furthermore, with the proposed single-parameter adaptive law, we reduce the adaptive parameters of the fuzzy system to one, which significantly reduces the complexity of the control design. In the numerical simulations and no-load experimental results, the tracking errors of the exoskeleton with the proposed SAFC are successfully bounded in the predefined ranges. Compared with the standard AFC, the simulation time of the proposed SAFC is significantly reduced, which validates the effectiveness of the single-parameter adaptive law. Furthermore, both passive and active gait training tasks with healthy volunteer subjects validate the effectiveness of the designed exoskeleton and the proposed SAFC method in practical implementations. In the future, the function named “assist as needed” between fully active and passive gait training is the ideal function we should achieve in the exoskeleton, which will facilitate increased patient participation in rehabilitation training.