An Adaptive Controller Based on IDA-PBC for Underactuated Mechanical Systems: Application to the Ball and Beam System

1. Introduction

The port-controlled Hamiltonian (PCH) model, which is regarded as another alternative model for the Euler-Lagrange model, is widely used to describe dynamic equations for nonlinear systems. The system described by PCH structure has many advantages: a number of natural physical systems are covered and significant structural properties are preserved. One of effective technologies used to control physical systems is interconnection and damping assignment passivity-based control (IDA-PBC) [1], which has been resoundingly used to solve the stabilization problems of various underactuated systems described by PCH framework. And this technology has been extensively used in induction machine [2], power converters [3], flexible spacecrafts [4] and aircrafts [5] and so on.

However, one of the main shortcomings of IDA-PBC method is that a set of partial differential equations (PDEs) need to be solved. In order to simplify this problem, outstanding contributions have been made by a large number of researchers. For instance, in [6], by parameterizing the expected inertia matrix, the potential PDE was enormously simplified, and this approach was extended to separable and nonseparable PCH systems. In order to ensure the solvability of PDEs, some conditions were added to the expected structure matrices $J_d$ and $R_d$, which were allowed to depend on the control input [7]. The good performance of this technique was demonstrated by the well-known boost power converter. In addition, some constructive solutions had also been proposed to simplify PDEs of underactuated mechanical systems (UMS) in [5,8,9,10] .

Many theoretical extensions and practical researches of the IDA-PBC approach had been reported in literature. In [11], two design methods of IDA-PBC were proposed in view of the existence of physical damping in hamiltonian frame. By combining the data sampling method with IDA-PBC, a sampling data controller [12] was designed, and the target dynamics was stabilized to the equilibrium point. In order to tolerate the limitation of actuator faults, IDA-PBC method with fault tolerance was improved in [13], and a high-gain adaptive IDA-PBC scheme was proposed. The effectiveness of the improved control law was verified by the experiment of a hexarotor UAV. .

Furthermore, the robustness of IDA-PBC strategies to disturbances has also been a hot topic in recent years. As reported in [14], an outer loop controller was designed to solve the matched disturbance suppression problem of UMS. In [15,16], a new IDA-PBC control law was constructed by combining model reference adaptive control method with IDA-PBC, which could more effectively compensate for disturbance compared to the standard IDA-PBC in [1]. In [17,18], a method of adding integral effects to IDA-PBC was presented for a kind of UMS with constant disturbances. In order to solve the problem of matched and unmatched disturbance suppression, specific coordinate changes were added to the damping term in [19]. As far as UMS are concerned, external interference is also abundant, which can not be ignored during system modeling. In [20], IDA-PBC approach was applied to the inertial wheel inverted pendulum, and the results showed that it had good robustness to external interference. Considering that the inertia matrix depends on non-actuated coordinates for underactuated systems, an integral effect with specific coordinate transformations was added to the outer-loop of the IDA-PBC scheme in [21]. The designed control scheme was applied to a UAV, which proved its effectiveness. Besides, the influence of viscous friction was studied by using the controlled Lagrangians method [22], and the closed-loop system was more stable. In [23], the IDA-PBC strategy was used to analyze the continuous friction.

Considering the above situation, an IDA-PBC control scheme based on adaptive method is proposed in this paper in view of the unknown frictions in UMS and uncertainties in the modeling process, which are better compensated. Only the matched input disturbances were considered in [14,24], only external frictions of the system were compensated in [20,25]. Finally, the uncertainties in friction and potential energy were handled respectively in [26]. Compared with the above, the uncertainties in external frictions, kinetic energy M, and potential energy V are estimated adaptively in this paper, which expand the research scope.

The main contributions can be summarized as

(1) An adaptive controller is designed for UMS.

(2) The estimate values of the unknown terms are placed in the damping injection controller $u_{di}$ instead of the energy shaping controller $u_{es}$, which simplifies the solution of partial differential equations.

(3) By using the LaSalle’s invariance principle and approximate linearization, the locally asymptotic stability of the state of the ball and beam system is achieved.

The rest of the paper is organized as follows: In Section 2, the design steps of IDA-PBC are briefly reviewed, and the problems to be solved are formulated. A new adaptive controller is proposed and the stability analysis is given in Section 3. In Section 4, the new control scheme is applied to the ball and beam system, and numerical simulation results are provided. Finally, the summary is presented in Section 5.

2. Problem Statement

In this section, the standard IDA-PBC method [1] for UMS [27] is briefly looked back. The various possible uncertainties are discussed and the PCH system with uncertainties is presented.

2.1. Review on IDA-PBC design

Consider a mechanical system defined by

$$M\left(q\right)\ddot{q}+\left(M\left(q\right)-\frac{1}{2}\frac{\partial {\dot{q}}^{T}M\left(q\right)}{\partial q}\right)\dot{q}+{\nabla }_{q}V\left(q\right)=G\left(q\right)u,$$

(1)

where $M\left(q\right)=M^T\left(q\right)>0$ is the inertia matrix and $V\left(q\right)$ is the potential energy function, $q\in {\mathbb{R}}^n$ is the generalized position, $u\in {\mathbb{R}}^m$, $m\le n$ is the control input. The matrix $G\left(q\right)\in {\mathbb{R}}^{n\times m}$ is an input matrix. The system is called fully actuated when $m=n$ and $rank\left(G\right)=m=n$, whereas it is called underactuated when $m<n$ and $rank\left(G\right)=m<n$. ${\nabla }_{q}V$ is the gradient of $V\left(q\right)$, i.e. ${\nabla }_{q}V\left(q\right)=\frac{\partial V\left(q\right)}{\partial q}$.

The Hamiltonian function H, which is defined as the sum of the kinetic energy and the potential energy, is the total energy of the system. It can be written as

$$H\left(q,p\right)=\frac{1}{2}{p}^T{M}^{-1}\left(q\right)p+V\left(q\right),$$

(2)

where $p=M\dot{q}$ is momenta. Then, the dynamic equation (1) can be represented as the following PCH form

$$\begin{array}{ccc}\hfill \left[\begin{array}{c}\dot{q}\\ \dot{p}\end{array}\right]& =& \left[\begin{array}{cc}0& I_n\\ -I_n& 0\end{array}\right]\left[\begin{array}{c}{\nabla }_{q}H\\ {\nabla }_{p}H\end{array}\right]+\left[\begin{array}{c}0\\ G\left(q\right)\end{array}\right]u\hfill \\ \hfill y& =& G^T\left(q\right){\nabla }_{p}H,\hfill \end{array}$$

(3)

where $y\in {\mathbb{R}}^m$ is the output. $I_n$ represents the $n\times n$ identity matrix.

The IDA-PBC method is composed of two parts, namely energy shaping and damping injection, i.e.

1) Energy Shaping:The state feedback controller $u_{es}$ should be designed so that the closed-loop system takes the following form

$$\left[\begin{array}{cc}0& I_n\\ -I_n& 0\end{array}\right]\left[\begin{array}{c}{\nabla }_{q}H\\ {\nabla }_{p}H\end{array}\right]+\left[\begin{array}{c}0\\ G\left(q\right)\end{array}\right]u_{es}=\left[\begin{array}{cc}0& M^{-1}{M}_d\\ -M_d{M}^{-1}& J_2\left(q,p\right)\end{array}\right]\left[\begin{array}{c}{\nabla }_q{H}_d\\ {\nabla }_p{H}_d\end{array}\right].$$

(4)

where $J_2=-J_2^T=\left[\begin{array}{cc}0& j\\ -j& 0\end{array}\right]$ is a free parameter and $H_d$ is the desired Hamiltonian function, which is defined by

$$H_d\left(q,p\right)=\frac{1}{2}{p}^T{M}_d^{-1}\left(q\right)p+V_d\left(q\right),$$

(5)

with $M_d\left(q\right)=M_d^T\left(q\right)>0$ and $V_d\left(q\right)$ representing the desired inertia matrix and the desired potential energy, respectively. It is assumed that there is an isolated minimum at the desired equilibrium point $q_{*}$, i.e.

$$q_{*}=argmin{H}_d\left(q\right)=argmin{V}_d\left(q\right).$$

(6)

The equation (6) holds if the following conditions are met.

Condition 1:${\nabla }_q{V}_d\left(q_{*}\right)=0$.

Condition 2:${\nabla }_q^2{V}_d\left(q_{*}\right)>0$.

It can easily verified that the first line of the equation (4) is satisfied. The second line of the equation (4) can be written as

$$-{\nabla }_{q}H+G\left(q\right)u_{es}=J_2\left(q,p\right){\nabla }_p{H}_d-M_d{M}^{-1}{\nabla }_q{H}_d,$$

which is equivalent

$$G\left(q\right)u_{es}={\nabla }_{q}H+J_2\left(q,p\right){\nabla }_p{H}_d-M_d{M}^{-1}{\nabla }_q{H}_d.$$

(7)

Let $G^{\perp }$ represent a full rank left annihilator of G, i.e. $G^{\perp }G=0$. As a result, multiplying (7) by $G^{\perp }$ from the left hand side gives

$$G^{\perp }\left\{{\nabla }_{q}H+J_2{\nabla }_p{H}_d-M_d{M}^{-1}{\nabla }_q{H}_d\right\}=0.$$

(8)

The PDE (8) can be equivalently deduced as the following two PDEs:

$$\begin{array}{ccc}\hfill G^{\perp }\left\{{\nabla }_q\left(p^T{M}^{-1}p\right)-M_d{M}^{-1}{\nabla }_q\left(p^T{M}_d^{-1}p\right)+2J_2{\nabla }_p{H}_d\right\}& =& 0,\hfill \end{array}$$

(9)

$$\begin{array}{ccc}\hfill G^{\perp }\left\{{\nabla }_{q}V-M_d{M}^{-1}{\nabla }_q{V}_d\right\}& =& 0.\hfill \end{array}$$

(10)

The energy shaping control law $u_{es}$ can be determined as

$$u_{es}={\left(G^{T}G\right)}^{-1}{G}^T\left({\nabla }_{q}H-M_d{M}^{-1}{\nabla }_q{H}_d+J_2{M}_d^{-1}p\right).$$

(11)

by solving (9) for $J_2$ and $M_d$, (10) for $V_d$.

2) Damping Injection:The object is to design a damping injection controller

$$u_{di}=-K_v{G}^T{\nabla }_p{H}_d,$$

(12)

where $K_v=K_v^T>0$ is a parameter matrix.

By using the controller

$$u=u_{es}+u_{di},$$

(13)

the given PCH system (3) is made to have the following expected PCH dynamics

$$\begin{array}{ccc}\hfill \left[\begin{array}{c}\dot{q}\\ \dot{p}\end{array}\right]& =& \left[J_d\left(q,p\right)-R_d\left(q,p\right)\right]\left[\begin{array}{c}{\nabla }_q{H}_d\\ {\nabla }_p{H}_d\end{array}\right]\hfill \end{array}$$

(14)

$$\begin{array}{ccc}\hfill y_d& =& G^T\left(q\right){\nabla }_p{H}_d,\hfill \end{array}$$

(15)

where $J_d$ and $R_d$ are the redistributed expected interconnection and damping matrices defined by

$$\begin{array}{ccc}\hfill J_d& =& -J_d^T=\left[\begin{array}{cc}0& M^{-1}{M}_d\\ -M_d{M}^{-1}& J_2\end{array}\right],\hfill \\ \hfill R_d& =& R_d^T=\left[\begin{array}{cc}0& 0\\ 0& G{K}_v{G}^T\end{array}\right]\ge 0.\hfill \end{array}$$

The desired Hamiltonian function (5) is considered as a candidate Lyapunov function. Its derivative is

$$\begin{array}{ccc}\hfill {\dot{H}}_d& =& {\left({\nabla }_p{H}_d\right)}^T\dot{p}+{\left({\nabla }_q{H}_d\right)}^T\dot{q}\hfill \\ & =& p^T{M}_d^{-1}\left(-{\nabla }_{q}H+G\left(u_{es}+u_{di}\right)\right)+{\left({\nabla }_q{H}_d\right)}^T{M}^{-1}p\hfill \\ & =& p^T{M}_d^{-1}{J}_2{M}_d^{-1}p-p^T{M}_d^{-1}G{K}_v{G}^T{\nabla }_p{H}_d\hfill \\ & =& -{\left({\nabla }_p{H}_d\right)}^{T}G{K}_v{G}^T{\nabla }_p{H}_d\le 0.\hfill \end{array}$$

Thus, $\left(q_{*},0\right)$ is a stable equilibrium point of the system (14) and (15). In addition, if the zero-state of the closed loop system (14) and (15) can be detectable from its output (15), then the equilibrium point $\left(q_{*},0\right)$ is asymptotically stable.

2.2. Possible uncertainties

It is usually assumed that all parameters of the system are known when using the IDA-PBC method. However, uncertainties exist inevitably in reality, which may lead to poor control performance. In PCH system (3), uncertainties might occur in the the Hamiltonian function H. In addation, frictions f exist in almost all mechanical systems.

Consider the following dynamic model

$$M\ddot{q}+\Delta M\ddot{q}+\left(M-\frac{1}{2}\frac{\partial {\dot{q}}^{T}M}{\partial q}\right)\dot{q}+\left(\Delta M-\frac{1}{2}\frac{\partial {\dot{q}}^T\Delta M}{\partial q}\right)\dot{q}+\left({\nabla }_{q}V+{\nabla }_q\Delta V\right)=G{u}^{*}-f,$$

(16)

where $\Delta M$ and $\Delta V$ denote unknown terms in M and V, respectively and $f=diag\left\{{\dot{q}}_1,{\dot{q}}_2,\cdots ,{\dot{q}}_n\right\}\theta$ represents frictions with $\theta ={\left[{\theta }_1,{\theta }_2,\cdots ,{\theta }_n\right]}^T$ being frictional coefficients. Define H as (2) and $p=M\dot{q}$. Then, it can be verified that (16) can be changed to the following PCH form

$$\left[\begin{array}{c}\dot{q}\\ \dot{p}\end{array}\right]=\left[\begin{array}{cc}0& I_n\\ -I_n& 0\end{array}\right]\left[\begin{array}{c}{\nabla }_{q}H\\ {\nabla }_{p}H\end{array}\right]+\left[\begin{array}{c}0\\ G\end{array}\right]u^{*}-\left[\begin{array}{c}0\\ {\Phi }^T\beta \end{array}\right],$$

(17)

where ${\Phi }^T\beta =\Delta M\ddot{q}+\left(\Delta M-\frac{1}{2}\frac{\partial {\dot{q}}^T\Delta M}{\partial q}\right)\dot{q}+{\nabla }_q\left(\Delta V\right)+f$ represents a parameterization of the uncertainties with $\beta ={\left[{\beta }_1,{\beta }_2,\cdots ,{\beta }_n\right]}^T$ being a vector of unknown constant parameters.

3. Controller design and stability analysis

In this section, an adaptive controller is designed to compensate for uncertainties, which was discussed in the previous section.

Define

$$H_d^{*}=\frac{1}{2}{p}^T{M}_d^{-1}p+V_d+\frac{1}{2\gamma }{\tilde{\beta }}^T\tilde{\beta },$$

(18)

where $\tilde{\beta }=\beta -\widehat{\beta }$, $\widehat{\beta }$ is the estimated value of unknown constant parameter and $\gamma =diag\left\{{\gamma }_1,{\gamma }_2,\cdots ,{\gamma }_n\right\}>0$ is the controller parameters.

An energy shaping controller $u_{es}^{*}$ should be constructed so that

$$\begin{array}{ccc}& & \left[\begin{array}{cc}0& I_n\\ -I_n& 0\end{array}\right]\left[\begin{array}{c}{\nabla }_{q}H\\ {\nabla }_{p}H\end{array}\right]+\left[\begin{array}{c}0\\ G\end{array}\right]u_{es}^{*}-\left[\begin{array}{c}0\\ {\Phi }^T\beta \end{array}\right]\hfill \\ & =& \left[\begin{array}{cc}0& M^{-1}{M}_d\\ -M_d{M}^{-1}& J_2\end{array}\right]\left[\begin{array}{c}{\nabla }_q{H}_d^{*}\\ {\nabla }_p{H}_d^{*}\end{array}\right]-\left[\begin{array}{c}0\\ {\Phi }^T\tilde{\beta }\end{array}\right].\hfill \end{array}$$

(19)

The first line of the equation is satisfied automatically, but the second line becomes

$$G{u}_{es}^{*}={\nabla }_{q}H-M_d{M}^{-1}{\nabla }_q{H}_d+J_2{\nabla }_p{H}_d^{*}.$$

(20)

By premultiplying $G^{\perp }$, it follows from the above equation that

$$G^{\perp }\left\{{\nabla }_{q}H+J_2{\nabla }_p{H}_d^{*}-M_d{M}^{-1}{\nabla }_q{H}_d\right\}=0.$$

(21)

(21) can be divided into the following two PDEs:

$$\begin{array}{ccc}\hfill G^{\perp }\left\{{\nabla }_q\left(p^T{M}^{-1}p\right)-M_d{M}^{-1}{\nabla }_q\left(p^T{M}_d^{-1}p\right)+2J_2{\nabla }_p{H}_d^{*}\right\}& =& 0,\hfill \end{array}$$

(22)

$$\begin{array}{ccc}\hfill G^{\perp }\left\{{\nabla }_{q}V-M_d{M}^{-1}{\nabla }_q{V}_d\right\}& =& 0.\hfill \end{array}$$

(23)

Premultiplying (20)by$G^T$ and solving for $u_{es}^{*}$ produce

$$u_{es}^{*}={\left(G^{T}G\right)}^{-1}{G}^T\left({\nabla }_{q}H-M_d{M}^{-1}{\nabla }_q{H}_d^{*}+J_2{\nabla }_p{H}_d^{*}\right).$$

(24)

The derivative of (18) along the trajectories of (30) and (31) is

$$\begin{array}{ccc}\hfill {\dot{H}}_d^{*}& =& p^T{M}_d^{-1}\dot{p}+{\left({\nabla }_q{V}_d\right)}^T\dot{q}-\frac{1}{\gamma }{\tilde{\beta }}^T\stackrel{·}{\widehat{\beta }}\hfill \\ & =& p^T{M}_d^{-1}\left(-{\nabla }_{q}H+G\left(u_{es}^{*}+u_{di}^{*}\right)-{\Phi }^T\beta \right)+{\left({\nabla }_q{V}_d\right)}^T{\nabla }_{p}H-\frac{1}{\gamma }{\tilde{\beta }}^T\stackrel{·}{\widehat{\beta }}\hfill \\ & =& p^T{M}_d^{-1}\left(-M_d{M}^{-1}{\nabla }_q{H}_d^{*}+J_2{\nabla }_p{H}_d^{*}+G{u}_{di}^{*}-{\Phi }^T\beta \right)+{\left({\nabla }_q{V}_d\right)}^T{M}^{-1}p-\frac{1}{\gamma }{\tilde{\beta }}^T\stackrel{·}{\widehat{\beta }}\hfill \\ & =& -{\left(p^T{M}^{-1}{\nabla }_q{H}_d^{*}\right)}^T+p^T{M}_d^{-1}{J}_2{\nabla }_p{H}_d^{*}-p^T{M}_d^{-1}{\Phi }^T\beta +p^T{M}_d^{-1}G{u}_{di}^{*}\hfill \\ & & +{\left({\nabla }_q{V}_d\right)}^T{M}^{-1}p-\frac{1}{\gamma }{\tilde{\beta }}^T\stackrel{·}{\widehat{\beta }}\hfill \\ & =& -\frac{1}{\gamma }{\tilde{\beta }}^T\stackrel{·}{\widehat{\beta }}+p^T{M}_d^{-1}G{u}_{di}^{*}-p^T{M}_d^{-1}{\Phi }^T\beta .\hfill \end{array}$$

(25)

Define an adaptive law

$$\stackrel{·}{\widehat{\beta }}=-\gamma \Phi {\nabla }_p{H}_d^{*},$$

(26)

and choose a damping injection controller

$$u_{di}^{*}={\left(G^{T}G\right)}^{-1}{G}^T\left(-G{K}_v{G}^T{\nabla }_p{H}_d+{\Phi }^T\widehat{\beta }\right),$$

(27)

By substituting (26) and (27) into (25), the following inequality is obtained

$${\dot{H}}_d^{*}=-p^T{M}_d^{-1}G{K}_v{G}^T{M}_d^{-1}p\le 0.$$

(28)

It can be proved that with the controller

$$u^{*}=u_{es}^{*}+u_{di}^{*},$$

(29)

and the adaptive law (26), the closed-loop system can be expressed as

$$\begin{array}{ccc}\hfill \left[\begin{array}{c}\dot{q}\\ \dot{p}\end{array}\right]& =& \left[\begin{array}{cc}0& M^{-1}{M}_d\\ -M_d{M}^{-1}& J_2\end{array}\right]\left[\begin{array}{c}{\nabla }_q{H}_d^{*}\\ {\nabla }_p{H}_d^{*}\end{array}\right]-\left[\begin{array}{c}0\\ {\Phi }^T\tilde{\beta }\end{array}\right]\hfill \end{array}$$

(30)

$$\begin{array}{ccc}\hfill y_d^{*}& =& G^T{\nabla }_p{H}_d^{*}.\hfill \end{array}$$

(31)

Proposition 1. 

Assume that the detectability condition of the output (31) is satisfied. Then, $\left(q_{*},0\right)=\left(0,0\right)$ is a locally asymptotically stable equilibrium point of the closed-loop system (30) and (31). However, $\tilde{\beta }$ is just a stable equilibrium point.

Proof of Proposition 1. 

It follows from (25) that the desired equilibrium point $\left(q_{*},0,\tilde{\beta }\right)$ is stable. Furthermore, since the output (31) is detectable, the local asymptotic stability of the state is guaranteed. □

4. Example: the ball and beam system

In this section, the well-known ball and beam system [28] is considered to have uncertainties including uncertainties of friction coefficients and uncertainties in the Hamiltonian function H.

4.1. System model

The dynamic behavior of the ball and beam system is described as

$$\begin{array}{ccc}\hfill {\ddot{q}}_1+gsin\left(q_2\right)-q_1{\dot{q}}_2^2+{\beta }_1{\dot{q}}_1& =& 0\hfill \\ \hfill \left(L^2+q_1^2\right){\ddot{q}}_2+2q_1{\dot{q}}_1{\dot{q}}_2+g{q}_{1}cos\left(q_2\right)+{\beta }_2{\dot{q}}_2& =& u,\hfill \end{array}$$

(32)

where $q_1$ and $q_2$ are the position of the ball and the angle of the bar, respectively, L is the length of the bar, ${\beta }_1,{\beta }_2$ are friction coefficients. Due to the ball being always maintained on the beam, the angle of the bar $q_2$ is assumed to be $q_2\in \left(-\pi ,\pi \right)$. Additionally, from (32), the inertia matrix M and potential energy function V are attained

$$\begin{array}{ccc}\hfill M\left(q\right)& =& \left[\begin{array}{cc}1& 0\\ 0& L^2+q_1^2\end{array}\right],\hfill \end{array}$$

(33)

$$\begin{array}{ccc}\hfill V\left(q\right)& =& g{q}_{1}sin\left(q_2\right),\hfill \end{array}$$

(34)

and $G={\left[\begin{array}{cc}0& 1\end{array}\right]}^T$.

The friction coefficients, the length of the bar L in the inertia matrix and the gravitational constant g in the potential energy function are assumed to be unknown, which includes both matched and unmatched. Therefore,

$$\begin{array}{ccc}\hfill {\Phi }^T\beta & =& \Delta M\ddot{q}+\left(\Delta M-\frac{1}{2}\frac{\partial {\dot{q}}^T\Delta M}{\partial q}\right)\dot{q}+{\nabla }_q\left(\Delta V\right)+f\hfill \\ & =& \left[\begin{array}{c}{\beta }_{3}sin\left(q_2\right)+{\beta }_1{\dot{q}}_1\\ {\beta }_4{\ddot{q}}_2+{\beta }_4{\dot{q}}_2+{\beta }_3{q}_{1}cos\left(q_2\right)+{\beta }_2{\dot{q}}_2\end{array}\right].\hfill \end{array}$$

(35)

where $\Delta M=\left[\begin{array}{cc}0& 0\\ 0& {\beta }_4\end{array}\right]$, $\Delta V={\beta }_3{q}_{1}sin\left(q_2\right)$,${\beta }_3$ and ${\beta }_4$are the unknown part of the gravitational constant $g$and the length of the bar L, respectively.

4.2. Controller design

Notice that (22) for solving $M_d$ is the same as the PDE given by [1], so $M_d$ is adopted as

$$M_d=\left(L^2+q_1^2\right)\left[\begin{array}{cc}\sqrt{2}{\left(L^2+q_1^2\right)}^{-\frac{1}{2}}& 1\\ 1& \sqrt{2\left(L^2+q_1^2\right)}\end{array}\right].$$

(36)

After substituting $M_d$ into (22), j is calculated as

$$j=q_1\left(p_1-\sqrt{\frac{2}{L^2+q_1^2}}{p}_2\right).$$

(37)

The potential energy is solved by substituting $M_d$ into PDE (23), which can be represented as

$$\sqrt{2\left(L^2+q_1^2\right)}\frac{\partial V_d}{\partial q_1}+\frac{\partial V_d}{\partial q_2}=gsin\left(q_2\right).$$

(38)

Solving (38) for $V_d$ results in

$$V_d\left(q\right)=g\left[1-cos\left(q_2\right)\right]+\frac{k_p}{2}{\left(q_2-\frac{1}{\sqrt{2}}arcsinh\left(\frac{q_1}{L}\right)\right)}^2,$$

where $k_p$ is the controller parameter.

By substituting $M_d$, j, and $H_d^{*}$ into (24), and after some straightforward calculations, the energy shaping term $u_{es}^{*}$ is expressed as

$$\begin{array}{ccc}\hfill u_{es}^{*}& =& \frac{q_1}{\sqrt{2}\left(L^2+q_1^2\right)}\times \left(-\sqrt{\left(L^2+q_1^2\right)}{p}_1^2+\sqrt{2}{p}_1{p}_2+\frac{p_2^2}{\sqrt{\left(L^2+q_1^2\right)}}\right)\hfill \\ & & +\phi \left(q\right),\hfill \end{array}$$

(39)

where

$$\begin{array}{ccc}\hfill \phi \left(q\right)& =& g{q}_{1}cos\left(q_2\right)-g\sqrt{2\left(L^2+q_1^2\right)}sin\left(q_2\right)\hfill \\ & & -k_p\sqrt{\frac{\left(L^2+q_1^2\right)}{2}}\left(q_2-\frac{1}{\sqrt{2}}arcsinh\left(\frac{q_1}{L}\right)\right).\hfill \end{array}$$

(40)

Furthermore, from (27), the damping injection term $u_{di}^{*}$ is determined as

$$\begin{array}{ccc}\hfill u_{di}^{*}& =& \frac{k_v}{\left(L^2+q_1^2\right)}\left(p_1-\sqrt{\frac{2}{\left(L^2+q_1^2\right)}}{p}_2\right)\hfill \\ & & +{\widehat{\beta }}_2{\dot{q}}_2+{\widehat{\beta }}_3{q}_{1}cos\left(q_2\right)+{\widehat{\beta }}_4{\dot{q}}_2+{\widehat{\beta }}_4{\ddot{q}}_2,\hfill \end{array}$$

(41)

where $k_v$ is a damping injection gain. And the adaptive law from (26) is constructed as

$$\begin{array}{ccc}\hfill \stackrel{·}{\widehat{\beta }}& =& -\gamma \Phi M_d^{-1}p\hfill \\ & =& -\left[\begin{array}{c}{\gamma }_1\left(\sqrt{2}\frac{p_1}{\sqrt{L^2+q_1^2}}-\frac{p_2}{L^2+q_1^2}\right){\dot{q}}_1\\ {\gamma }_2\left(\sqrt{2}\frac{p_2}{{\left(L^2+q_1^2\right)}^{\frac{3}{2}}}-\frac{p_1}{L^2+q_1^2}\right){\dot{q}}_2\\ {\gamma }_3\left(sin\left(q_2\right)\left(\sqrt{2}\frac{p_1}{\sqrt{L^2+q_1^2}}-\frac{p_2}{L^2+q_1^2}\right)+q_{1}cos\left(q_2\right)\left(\sqrt{2}\frac{p_2}{{\left(L^2+q_1^2\right)}^{\frac{3}{2}}}-\frac{p_1}{L^2+q_1^2}\right)\right)\\ {\gamma }_4\left(\sqrt{2}\frac{p_2}{{\left(L^2+q_1^2\right)}^{\frac{3}{2}}}-\frac{p_1}{L^2+q_1^2}\right)\left({\ddot{q}}_2+{\dot{q}}_2\right)\end{array}\right]\hfill \end{array}$$

(42)

4.3. Stability analysis

In this example, the LaSalle’s invariance principle is applied to prove the asymptotic stability of the closed-loop system. Under the control of (29) and (26), the state equations of the ball and beam system can be described as

$$\begin{array}{ccc}\hfill {\dot{q}}_1& =& p_1,\hfill \\ \hfill {\dot{p}}_1& =& -gsin\left(q_2\right)+q_1{\dot{q}}_2^2-{\beta }_1{\dot{q}}_1\hfill \\ & =& -gsin\left(q_2\right)+\frac{1}{{\left(L^2+q_1^2\right)}^2}{q}_1{p}_2^2-{\beta }_1{p}_1,\hfill \\ \hfill {\dot{q}}_2& =& \frac{1}{\left(L^2+q_1^2\right)}{p}_2,\hfill \\ \hfill {\dot{p}}_2& =& -g{q}_{1}cos\left(q_2\right)-{\beta }_2{\dot{q}}_2+u^{*}\hfill \\ & =& -\frac{{\beta }_2{p}_2}{\left(L^2+q_1^2\right)}+\frac{q_1}{\sqrt{2}\left(L^2+q_1^2\right)}\left(-\sqrt{\left(L^2+q_1^2\right)}{p}_1^2+\sqrt{2}{p}_1{p}_2+\frac{p_2^2}{\sqrt{\left(L^2+q_1^2\right)}}\right)\hfill \\ & & -g\sqrt{2\left(L^2+q_1^2\right)}sin\left(q_2\right)-\frac{\sqrt{2\left(L^2+q_1^2\right)}}{2}\left(q_2-\frac{1}{\sqrt{2}}arcsinh\left(\frac{q_1}{L}\right)\right)\hfill \\ & & +\frac{K_v}{\left(L^2+q_1^2\right)}\left(p_1-\sqrt{\frac{2}{\left(L^2+q_1^2\right)}}{p}_2\right)+{\widehat{\beta }}_2{\dot{q}}_2+{\widehat{\beta }}_3{q}_{1}cos\left(q_2\right)+{\widehat{\beta }}_4{\dot{q}}_2+{\widehat{\beta }}_4{\ddot{q}}_2.\hfill \end{array}$$

(43)

and the output equation is

$$\begin{array}{ccc}\hfill y_d^{*}& =& G^T{\nabla }_p{H}_d^{*}\hfill \\ & =& \frac{1}{\left(L^2+q_1^2\right)}\left(\sqrt{\frac{2}{\left(L^2+q_1^2\right)}}{p}_2-p_1\right).\hfill \end{array}$$

(44)

Restricted by manifold $y_d^{*}\equiv 0,$$\forall t$, the trajectories of the system (43) are analyzed as follows. It follows from $y_d^{*}=0$ that

$$p_2=\sqrt{\frac{L^2+q_1^2}{2}}{p}_1.$$

(45)

According to ${\dot{y}}_d^{*}=0$, one has

$$\begin{array}{ccc}& & gsin\left(q_2\right)+q_2-\frac{1}{\sqrt{2}}arcsinh\left(\frac{q_1}{L}\right)+\frac{q_1{p}_1^2}{\left(L^2+q_1^2\right)}\hfill \\ & & -{\beta }_1{p}_1+\frac{{\tilde{\beta }}_2{p}_1}{\left(L^2+q_1^2\right)}-\frac{\sqrt{2}{\widehat{\beta }}_3{q}_{1}cos\left(q_2\right)}{\sqrt{\left(L^2+q_1^2\right)}}\hfill \\ & & -\frac{{\widehat{\beta }}_4}{\left(L^2+q_1^2\right)}\left(p_1-\frac{q_1{p}_1^2}{2\left(L^2+q_1^2\right)}-gsin\left(q_2\right)-{\beta }_1{p}_1\right)\hfill \\ & =& 0.\hfill \end{array}$$

(46)

From equation (45), the ball and beam system can be reduced to the following system

$$\begin{array}{ccc}\hfill {\dot{q}}_1& =& p_1,\hfill \end{array}$$

(47)

$$\begin{array}{ccc}\hfill {\dot{p}}_1& =& \frac{1}{2}\frac{q_1{p}_1^2}{\left(L^2+q_1^2\right)}-gsin\left(q_2\right)-{\beta }_1{p}_1,\hfill \end{array}$$

(48)

$$\begin{array}{ccc}\hfill {\dot{q}}_2& =& \frac{1}{\sqrt{2\left(L^2+q_1^2\right)}}{p}_1.\hfill \end{array}$$

(49)

According to (47) and (49), it can be easily verified that

$$\begin{array}{ccc}& & \frac{d}{dt}\left(q_2-\frac{1}{\sqrt{2}}arcsinh\left(\frac{q_1}{L}\right)\right)\hfill \\ & =& \frac{1}{\sqrt{2\left(L^2+q_1^2\right)}}{p}_1-\frac{1}{\sqrt{2\left(L^2+q_1^2\right)}}{\dot{q}}_1\hfill \\ & =& 0,\hfill \end{array}$$

which means that

$$q_2\left(t\right)-\frac{1}{\sqrt{2}}arcsinh\left(\frac{q_1\left(t\right)}{L}\right)=2\delta ,\forall t,$$

(50)

It is known that the origin is an equilibrium point, which requires $\delta =0$. According to (50), $q_2\equiv \frac{1}{\sqrt{2}}arcsinh\left(\frac{q_1}{L}\right)$. After a series of simplifications, the final second-order system is described as:

$$\begin{array}{ccc}\hfill q_1& =& Lsinh\left(\sqrt{2}{q}_2\right),\hfill \end{array}$$

(51)

$$\begin{array}{ccc}\hfill {\dot{q}}_2& =& \frac{1}{\sqrt{2\left(L^2+q_1^2\right)}}{p}_1,\hfill \end{array}$$

(52)

$$\begin{array}{ccc}\hfill {\dot{p}}_1& =& \frac{1}{2}\frac{q_1{p}_1^2}{\left(L^2+q_1^2\right)}-gsin\left(q_2\right)-{\beta }_1{p}_1.\hfill \end{array}$$

(53)

Assumption A1. 

Assume that the following inequality is true.

$$-L^2{\beta }_1+{\tilde{\beta }}_2+{\widehat{\beta }}_4{\beta }_1-{\widehat{\beta }}_4\ne 0.$$

Lemma 1. 

With Assumption 1, if the trajectories of (51)-(53) are confined to (45) and (46), then $q_2\left(t\right)\ge 0$, $\forall t$ or $q_2\left(t\right)\le 0$, $\forall t$.

Proof of Lemma 1. 

Assume $q_2\left(t^{*}\right)=0$ at $t^{*}$. It can be obtained from (51) that $q_1\left(t^{*}\right)=0$. Substituting $\left(q_1,q_2\right)=\left(0,0\right)$ into (46) leads to

$$\left(-L^2{\beta }_1+{\tilde{\beta }}_2+{\widehat{\beta }}_4{\beta }_1-{\widehat{\beta }}_4\right)p_1=0.$$

(54)

According to Assumption 1, (54) has only one real solution $p_1=0$, which implies that $p_2=0$. Since $q_1\left(t^{*}\right)=0$ and $q_2\left(t^{*}\right)=0$, the fact that ${\dot{p}}_1\left(t^{*}\right)=0$ can be obtained from (51). Because $p_1\left(t^{*}\right)=0$ and ${\dot{p}}_1\left(t^{*}\right)=0$, $p_1\left(t\right)=0$ for $t>t^{*}$, which implies that $p_1\left(t\right)=0$ for $t>t^{*}$ due to (45). Since $q_2\left(t^{*}\right)=0$ and ${\dot{q}}_2\left(t\right)=0$ for $t>t^{*}$ due to $p_1\left(t\right)=0$ and (52), $q_2\left(t\right)=0$ for $t>t^{*}$. With (51), it follows from $q_2\left(t\right)=0$ that $q_1\left(t\right)=0$ for $t>t^{*}$. Finally, it can be concluded that if $q_2\left(t^{*}\right)=0$, the system will stay at the origin for $t>t^{*}$. □

Next, the final second-order system (52) and (53) is linearized at the origin to get

$$\begin{array}{ccc}\hfill {\dot{q}}_2& =& \frac{1}{\sqrt{2}L}{p}_1,\hfill \\ \hfill {\dot{p}}_1& =& -g{q}_2-{\beta }_1{p}_1.\hfill \end{array}$$

(55)

The system matrix A can be obtained from (55) as follows

$$A=\left[\begin{array}{cc}0& \frac{1}{\sqrt{2}L}\\ -g& -{\beta }_1\end{array}\right],$$

and eigenvalues of matrix A is calculated as

$$\begin{array}{ccc}\hfill {\lambda }_1& =& \frac{1}{2}\left(-{\beta }_1+\sqrt{{\beta }_1^2-\frac{2\sqrt{2}g}{L}}\right),\hfill \\ \hfill {\lambda }_2& =& \frac{1}{2}\left(-{\beta }_1-\sqrt{{\beta }_1^2-\frac{2\sqrt{2}g}{L}}\right).\hfill \end{array}$$

Whether ${\beta }_1^2-2\sqrt{2}\frac{g}{L}>0$ or ${\beta }_1^2-2\sqrt{2}\frac{g}{L}<0$ , the eigenvalues have a negative real part. So the linearized system (55) is asymptotically stable. Furthermore, the system (52) and (53) is locally asymptotically stable. Since $q_2\in \left(-\pi ,\pi \right)$, it can be easily verified that the origin is the only equilibrium point, so $q_2\to 0$. According to the proof of Lemma 1, $\left(q,p\right)=\left(0,0\right)$ can be deduced. Hence, the conclusion that $p=0$ is proved. The equilibrium point $\left(q_{*},0\right)=\left(0,0\right)$ is locally asymptotically stable.

4.4. Numerical simulation results

In this section, the proposed controller is simulated and compared with the article [1] under different initial conditions and controller parameters. The simulation results are shown in Figure 1, Figure 2, Figure 3, Figure 4, Figure 5, Figure 6, Figure 7, Figure 8, Figure 9, Figure 10, Figure 11 and Figure 12.

Case 1: The initial conditions are chosen as $\left[q,p\right]=\left[0,91.7^{\circ },0,0\right]$, $\left[{\widehat{\beta }}_1,{\widehat{\beta }}_2,{\widehat{\beta }}_3,{\widehat{\beta }}_4\right]=\left[0.1,0.2,0.5,0.25\right]$, the controller parameters are set to ${\gamma }_1=0.001$, ${\gamma }_2=0.1$, ${\gamma }_3=0.5$, ${\gamma }_4=0.01$, $k_p=1$, $k_v=20$ and the length of the bar is selected as $L=0.5m$. Figure 1, Figure 2, Figure 3 and Figure 4 show the results under Case 1. It can be obviously seen from Figure 1 that using the proposed adaptive controller (29) and (42), the ball reaches the expected position in about $20s$. But under the existing controller (13), the ball gradually stabilizes to the expected position after $70s$. The bar keeps swinging in the range of $\left(-41.6^{\circ },27.6^{\circ }\right)$ from $7s$ to $30s$, and the swing amplitude decreases significantly after $30s$ and the bar gradually reaches a stable state with the proposed adaptive controller in Figure 2. In addition, it can be observed from Figure 3 that the control signal of the proposed adaptive controller is much smaller and settles down faster than the existing controller. As depicted in Figure 4, the estimated values of the four parameters are all convergent, which proves the boundedness of the estimated parameters.

Case 2: The initial conditions are consistent with the selection of the first group, the controller parameters are given as ${\gamma }_1=0.001$, ${\gamma }_2=0.5$, ${\gamma }_3=1$, ${\gamma }_4=0.01$, $k_p=1$, $k_v=1$ and the length of the bar is also selected as $L=0.5m$. The results of this set of conditions are shown in Figure 5, Figure 6, Figure 7 and Figure 8. In terms of these conditions, the position of the ball attains balance in approximately $10s$ under the proposed adaptive controller in Figure 5. As shown in Figure 6, the angle of the bar reaches stability in about $20s$ under the proposed adaptive controller. However, with the existing controller, the angle of the bar undergos violent oscillation around $30s$ and reaches ${71}^{\circ }$ and the bar stabilizes after $50s$. Besides, the control effects are presented in Figure 7. It can be seen from Figure 8 that the estimated values of the four parameters converge after about $10s$.

Case 3: The initial conditions are selected as $\left[q,p\right]=\left[10,0^{\circ },0,0\right]$, $\left[{\widehat{\beta }}_1,{\widehat{\beta }}_2,{\widehat{\beta }}_3,{\widehat{\beta }}_4\right]=\left[0.1,0.2,0.05,0.25\right]$, the controller parameters are given as ${\gamma }_1=0.01$, ${\gamma }_2=0.2$, ${\gamma }_3=0.1$, ${\gamma }_4=0.01$, $k_p=1$, $k_v=50$ and the length of the bar is chosen as $L=10m$. At about $50s$, the ball oscillates and gradually settles down at the desired position with the proposed adaptive controller in Figure 9. At the same time, under the existing controller, the ball keeps swinging around the equilibrium position. Figure 10 shows that the swing amplitude of the bar gradually becomes smaller and smaller and the angle reaches steady state by the proposed adaptive controller. The vibration amplitude of the control signal significantly decreases under the application of the proposed adaptive controller in Figure 11. It can be seen from Figure 12 that the four estimated parameters are bounded.

5. Conclusions

In this paper, an adaptive control law has been designed for a class of underactuated mechanical systems with matched and unmatched uncertainties. With this controller, the locally asymptotic stability of the underactuated mechanical system is ensured under uncertainties. The estimate values of the unknown terms are placed in the damping injection controller $u_{di}$, which simplifies the design of the controller. The locally asymptotic stability of the ball and beam system is proved by using the LaSalle’s invariance principle and approximate linearization.