A model-based velocity controller for chaotization of flexible joint robot manipulators

Abstract

This article presents a model-based velocity controller able to induce a chaotic motion on n-degrees of freedom flexible joint robot manipulators. The proposed controller allows the velocity link vector of a robot manipulator to track an arbitrary, chaotic reference vector field. A rigorous theoretical analysis based on Lyapunov’s theory is used to prove the asymptotic stability of the tracking error signals when using the proposed controller, which implies that a chaotic motion is induced to the robotic system. Experimental results are provided using a flexible joint robot manipulator of two degrees of freedom. Finally, by using Poincaré maps and Lyapunov exponents, it is shown that the behavior exhibited by the robot joint positions is chaotic.

Introduction

Overview

Typically, rigid-link robot manipulators do not reflect a real scenario when a high performance is required. If the robot moves at high velocities or interacts with the environment, some deflections may occur in the robot components, inducing vibration to the robot. Undesired vibrations may occur by the joint flexibility, which is caused by elasticity of motion transmission elements, such as harmonic drives, belts, or long shafts. ¹

Some applications of robot manipulators where joint flexibility occurs include compliant manipulation in contact with a human environment, requiring the robot to have low reflected inertia and high accuracy in force control, ² which may be accomplished by increasing the joint flexibility of the robot manipulator. ³ Besides, the mass of flexible joint robot (FJR) manipulators is smaller than that of rigid robots. This feature allows them to move faster, while consuming a less amount of energy.

The term chaos is used to describe a nonperiodic and outwardly random behavior that happens in some nonlinear dynamic deterministic systems. Such a behavior exists because these nonlinear systems have ^{4 –6} high sensitivity to initial conditions (the butterfly effect ⁷ ), that is, neighboring states diverge as time proceeds; topological transitivity, that is, every non-empty open subset of the phase space is traveled by the trajectory originating from an arbitrary open subset; and denseness of periodic orbits, which yields a broad frequency spectrum.

There exists a large variety of applications in science and engineering where chaos can be beneficial. Some of them include applications in medicine and electronics. Other important applications are mixing of fluids and granular flows by inducing a chaotic motion on direct current (DC) motors and rigid-link robot manipulators, path planning for mobile robots, vibration control of FJR manipulators, parameter identification, motion control, manufacturing human-like robots, bipedal walking, among others. ^{8 –16} Therefore, the problem of inducing a chaotic motion to an FJR manipulator is an important issue that deserves to be studied.

Although the number of applications where the existence of chaos is required is increasing, there are a few works presenting methodologies for inducing a chaotic behavior in robot manipulators. In addition, the existing literature only considers applications on rigid-link robot manipulators or FJR manipulators of 1-degree of freedom (DOF). It is also worth noting that for a general n-DOF FJR manipulator the control task increases its complexity because the number of DOF is enlarged twice with respect to the DOF of a rigid-link robot manipulator. Hence, in order to contribute to this field in a more realistic way, this work presents a model-based velocity controller that allows inducing a chaotic motion on FJR manipulators of n-DOF. Furthermore, the proposed controller is validated through experimental results using a 2-DOF FJR manipulator.

Literature review

The mathematical model describing the dynamics of an FJR manipulator is important for the control synthesis. It is well known that the general dynamic model of an FJR manipulator may not satisfy neither the necessary conditions for input–output decoupling nor those for full linearization by static state feedback. ¹⁷ Instead, dynamic state feedback controllers can be employed to linearize the system, because FJRs are input–output invertible systems with stable zero dynamics. ¹⁸ On the other hand, a simplified model, which assumes that the angular part of the kinetic energy of each rotor is due only to its own rotation, can be used. ¹⁹ The advantage of using a simplified model is that it is globally feedback linearizable by static state feedback. Finally, note that many robust control techniques such as sliding modes can be used to control rigid n-link robot manipulators. ²⁰ However, the performance of such methodologies worsens when the joint flexibility is not negligible.

Methodologies for controlling FJR manipulators include state feedback controllers, ²¹ passivity-based impedance control, ^22,23 adaptive techniques, ^24,25 sliding mode control (SMC), ^26,27 adaptive SMC, ²⁸ terminal SMC, ²⁹ among others. Besides, an interesting survey on the control of FJR manipulators is presented in Ozgoli and Taghirad. ³⁰

Chaotic systems can be present and used in many areas including biological, electrical and mechanical sciences, control and synchronization, parameter identification and applications involving robot manipulators. ^{5,6,8 –16} Thus, it evolved a field of study known as chaotization, chaos synthesis or chaos anti-control, which aims to establish a desired chaotic behavior in an otherwise non-chaotic system, or to enhance the chaos of a chaotic system.

It is important to remark that the concepts of tracking, synchronization, and chaotization seem to be similar. Hence, it is useful to distinguish each from the others. To this end, let y(t) and $y_d\left(t\right)$ be the output and the reference of a given control system, respectively. The tracking problem consists in designing a control input u(t) which drives the output y(t) to the desired reference $y_d\left(t\right)$ , usually in asymptotic form, that is, ${{\displaystyle \underset{t\to \infty }{\text{lim}}}}_{\hspace{0.17em}}\left[y\left(t\right)-y_d\left(t\right)\right]=0$ . In a synchronization problem, the states of two or more systems change in a coordinated way or some scalar characteristics of the systems (e.g. amplitudes, frequencies, energies, powers, fractal dimensions) change concurrently. ^31,32 Finally, a chaotization problem can be stated as a tracking problem, but the goal is to meet some formal criteria of chaos. ³³

Chaotization for mechanical systems can be accomplished by means of several methods including proportional derivative control, ³⁴ saturated proportional control, ³⁵ time-delay feedback control applied to DC motors for industrial mixing, ³⁶ chaotization of brushless DC motor systems, ^{37 –39} adaptive techniques for rigid robot manipulators, ⁴⁰ feedback linearization for mobile robots, ⁴¹ proportional integral derivative control, ⁴² neural networks, ⁴³ path planning generators for autonomous robots, ^{44 –46} mechanism synthesis, ⁴⁷ computing robot kinematics, ⁴⁸ parameter identification, ⁴⁹ passive dynamic walking, ⁵⁰ among others. There are also some studies showing the importance of chaotization in FJR manipulators. In a study by Kandroodi et al., ⁵¹ chaotization was used to achieve trajectory tracking and vibration control of an FJR manipulator. Similar procedures for vibration suppression in FJR manipulators have been studied using a variable structure control, ⁵² a Lyapunov ruled-based fuzzy approach, ⁵³ feedback linearization, ⁵⁴ a Lyapunov-based approach, ⁵⁵ and neural networks. ⁵⁶ Finally, Singh et al. ⁵⁷ studied a 1-DOF FJR to show the occurrence of various chaotic behaviors.

The previous literature review shows the importance of inducing a chaotic behavior in non-chaotic systems. However, the aforementioned works only consider the chaotification of rigid-link robotic manipulators, or present methodologies that can be applied to 1-DOF FJR manipulators.

Contribution of the article

This work proposes a model-based velocity controller applied to an n-DOF FJR manipulator. An extended system is created by combining the dynamic model of the FJR manipulator and an arbitrary chaotic system. Then, it is theoretically shown that a chaotic behavior is generated on the robotic system when using the proposed controller.

A rigorous stability analysis based on a strict Lyapunov function demonstrates the asymptotic stability of all the tracking error signals. Finally, the proposed chaotification technique is validated in an experimental setup, which consists of a 2-DOF FJR manipulator from Quanser (Quanser Consulting, Markham, ON, Canada).

From the previous literature review, we concluded that the existing work in chaotification of robot manipulators has been focused on rigid-link manipulators or 1-DOF FJR manipulators. Hence, this article contributes by extending the existing results to a more realistic scenario, where the chaotification problem is solved for a general n-DOF robot manipulator where the manipulator joint flexibility is taken into account.

The main contributions of this article are:

developing a chaotification methodology that can be applied to an n-DOF FJR manipulator; and

presenting experimental results corroborating the validity of the proposed scheme.

Organization

The rest of the article is structured as follows. The system description and problem formulation are stated in the second section. The model-based velocity controller design, which accomplishes the chaotization goal, and its Lyapunov stability proof are presented in the third section. The fourth section is devoted to the experimental validation of the proposed chaotization methodology, using a Quanser 2-DOF FJR manipulator. Finally, some concluding remarks are given in the fifth section.

Notation

Let ℝ denote the set of real numbers, ℕ represents the natural numbers, ‖⋅‖ refers to the Euclidean norm of a given vector or matrix, and |⋅| denotes the absolute value of its argument. In general, vectors are given by bold lowercase letters and matrices are represented by uppercase letters; $0_{n\times n}$ is the zero matrix of order n × n and $I_{n\times n}$ the identity matrix of order n × n. Given a time-varying element v(t), the expressions $\dot{v}\left(t\right),\ddot{v}\left(t\right),\hspace{0.17em}\text{and}\hspace{0.17em}{v}^{\left(n\right)}\left(t\right)$ represent the first, second, and n-th time derivative of v(t), respectively. Matrix $A^{-1}$ denotes the inverse of matrix A, when it exists, and $B=\text{diag}\left\{b_1,\dots ,b_n\right\}\in {ℝ}^{n\times n}$ is a diagonal matrix with elements $b_1,\dots ,b_n$ on the main diagonal. For a given symmetric matrix, ${\lambda }_{\text{min}}(\cdot )\hspace{0.17em}\text{and}\hspace{0.17em}{\lambda }_{\text{max}}(\cdot )$ represent the minimum and maximum eigenvalue of the matrix, respectively. When required, the square brackets [⋅] are used to indicate the units of a given quantity. Finally, if no confusion is present, the arguments of functions will be omitted.

System description and problem formulation

The proposed chaotization procedure is applied to n-DOF FJR manipulators, that is, robots with rigid links but with elastic transmissions, where the flexibility is concentrated at the joints. In the following, the dynamic model of an n-DOF FJR manipulator is given. Then, the problem formulation for its chaotization is stated.

Dynamic model of the FJR

Let us consider an FJR manipulator of n-DOF, with rotational joints actuated by DC motors. It is assumed that linear torsional springs, with stiffness constants $k_{si}$ , $i=1,\dots ,n$ , represent the elasticity of the i-th joint. This implies that the rotor of each actuator can be modeled as an additional rigid body in the chain with its own inertia. It is also assumed that the elastic forces are limited to the linear range of operation. ¹⁹

A schematic representation of an n-DOF FJR manipulator is depicted in Figure 1. In this figure, each rotor is coupled to the link by means of a torsional spring. Variable q_i, $i=1,\dots ,n$ , corresponds to the position of link i, which is measured around the axis z_i, according to the Denavit–Hartenberg convention. ^58,59 The variable ϕ_i, $i=1,\dots ,n$ , stands for the position of rotor i after the gearbox. Rotors or links positive displacements are measured in the counterclockwise direction.

The dynamic equations of an FJR manipulator are given by ^19,59,60

$M_1\left(\mathit{q}\right)\ddot{\mathit{q}}+N\left(\mathit{q},\dot{\mathit{q}}\right)+K\left[\mathit{q}-\varphi \right]+G\left(\mathit{q}\right)=0$

1

$M_2\ddot{\varphi }+F_2\left(\dot{\varphi }\right)+K\left[\varphi -\mathit{q}\right]=\mathit{\tau }$

2

with

$\mathit{q}={\left[q_1,\dots ,q_n\right]}^{\text{T}}\in {ℝ}^n$

$\varphi ={\left[{\varphi }_1,\dots ,{\varphi }_n\right]}^{\text{T}}\in {ℝ}^n$

$\mathit{\tau }={\left[{\tau }_1,\dots ,{\tau }_n\right]}^{\text{T}}\in {ℝ}^n$

$N\left(\mathit{q},\dot{\mathit{q}}\right)=C_A\left(\mathit{q},\dot{\mathit{q}}\right)\dot{\mathit{q}}+F_1\left(\dot{\mathit{q}}\right)$

$F_1\left(\dot{\mathit{q}}\right)=F_{v1}\dot{\mathit{q}}+{\mathit{f}}_{c1}\left(\dot{\mathit{q}}\right)\in {ℝ}^n$

$F_2\left(\dot{\varphi }\right)=F_{v2}\dot{\varphi }+{\mathit{f}}_{c2}\left(\dot{\varphi }\right)\in {ℝ}^n$

${\mathit{f}}_{c1}\left(\dot{\mathit{q}}\right)={\left[f_{c11}\text{tanh}\left(\beta {\dot{q}}_1\right),\dots ,f_{c1n}\text{tanh}\left(\beta {\dot{q}}_n\right)\right]}^{\text{T}}$

${\mathit{f}}_{c2}\left(\dot{\varphi }\right)={\left[f_{c21}\text{tanh}\left(\beta {\dot{\varphi }}_1\right),\dots ,f_{c2n}\text{tanh}\left(\beta {\dot{\varphi }}_n\right)\right]}^{\text{T}}$

3

Figure 1.

Schematic representation of an n-DOF FJR manipulator. DOF: degree of freedom; FJR: flexible joint robot.

In equations (1) and (2), $M_1\left(\mathit{q}\right)\in {ℝ}^{n\times n}$ is the mass matrix of the rigid robot, $M_2\in \hspace{0.17em}{ℝ}^{n\times n}$ is the matrix corresponding to the rotor inertias, $C_A\left(\mathit{q},\dot{\mathit{q}}\right)\in {ℝ}^{n\times n}$ is the matrix of centrifugal and Coriolis torques, $K\in {ℝ}^{n\times n}$ is the matrix associated with the spring stiffness, $F_{v1},\hspace{0.17em}{F}_{v2}\in {ℝ}^{n\times n}$ are diagonal matrices containing the viscous friction coefficients, ${\mathit{f}}_{c1}\left(\dot{\mathit{q}}\right),{\mathit{f}}_{c2}\left(\dot{\varphi }\right)\in {ℝ}^n$ are the vectors corresponding to a continuous version of the Coulomb friction with β > 0 large enough, $G\left(\mathit{q}\right)\in {ℝ}^n$ is the vector of gravitational torques, and $\tau \in {ℝ}^n$ is the torque input vector. Coordinates q , $\varphi \in {ℝ}^n$ refer to the angular position of the links and the rotors after the gearbox, respectively. Finally, the mass matrix $M_1\left(\mathit{q}\right)$ is a symmetric, positive definite matrix, and the matrix $C_A\left(\mathit{q},\dot{\mathit{q}}\right)$ , defined via the Christoffel symbols, possesses the skewed-symmetric matrix property ${\mathit{x}}^{\text{T}}\left[{\dot{M}}_1-2C_A\right]\mathit{x}=0$ for any $\mathit{x}\in {ℝ}^n$ . ⁵⁹

Remark 1

The discontinuous function $\text{sign} (\cdot )$ defined as

$\text{sign}\left(x\right)=\{\begin{array}{c}-1\\ \in \left[-1,1\right]\\ 1\end{array}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\begin{array}{c}, x<0\\ , x=0\\ , x>0\end{array}$

4

in ${\mathit{f}}_{ci}(\cdot )$ ( $i=1,2$ ) provides a better description of the Coulomb friction. ⁶¹ If the sign function was considered, the dynamic equations (1) and (2) can be properly represented as a multivalued differential inclusion $\dot{\mathit{z}}\left(t\right)\in F\left(t,\mathit{z},\mathit{u}\right)$ , where $\mathit{z}\left(t\right)$ are the states, $\mathit{u}\left(t\right)$ is the control input, and $F\left(t,\mathit{z},\mathit{u}\right)$ is a piecewise continuous function. The solution of the differential equation can be defined in the sense of Filippov, Aizerman–Pyatnitskiy, or Gelig–Leonov–Yakubovich. Any solution to such differential inclusion has different behavior when compared to that obtained using a continuous approximation of the sign(⋅) function. Some examples demonstrating this change of behavior in the sliding regimes and the rest segment can be found in Leonov et al. ⁶² However, the appropriate analysis required to prove the convergence of the solutions is different from our approach. Hence, in this work, the continuous approximation of the Coulomb friction, given by the tanh(⋅) function, will be used in our control scheme for the model compensation. It should be pointed out that a better approximation to the discontinuous function can be obtained by increasing the magnitude of the parameter β. ⁶³

In the next development, it is assumed that position and velocity measurements of robot joints as well as links are available for feedback.

Chaotic velocity field and extended system

Let us consider a chaotic system described by the next autonomous system

$\frac{\text{d}}{\text{d}t}\mathit{x}=\mathit{f}\left(\mathit{x}\right)$

5

with $\mathit{x}\in {ℝ}^N$ . An example of a chaotic system expressed in equation (5) is the Lorenz system ⁶⁴

$\frac{\text{d}}{\text{d}t}\left[\begin{array}{c}{x}_1\\ x_2\\ x_3\end{array}\right]=\left[\begin{array}{c}\sigma \left[x_2-x_1\right]\\ r{x}_1-x_2-x_1{x}_3\\ x_1{x}_2-b{x}_3\end{array}\right]$

6

where the parameters σ and r are the Prandtl and the Rayleigh numbers, respectively. A chaotic behaviour is obtained by using $\sigma =10$ , $r=28$ , and $b=8/3$ . This chaotic motion is depicted in Figure 2.

Let us define the extended vector $\mathit{r}\left(t\right)$ as

$\mathit{r}\left(t\right)=\left[\begin{array}{c}\mathit{q}\left(t\right)\\ \mathit{p}\left(t\right)\end{array}\right]\in {ℝ}^N$

7

with $\mathit{p}\left(t\right)\in {ℝ}^{N-n}$ . Now, $\mathit{r}\left(t\right)$ and the structure of the chaotic system (equation (5)) are used to define the next extended system

$\frac{\text{d}}{\text{d}t}\dot{\mathit{r}}=\mathit{f}\left(\mathit{r}\right)=\left[\begin{array}{c}{\mathit{f}}_1\left(\mathit{r}\right)\\ {\mathit{f}}_2\left(\mathit{r}\right)\end{array}\right]$

8

with ${\mathit{f}}_1:{ℝ}^N\mapsto {ℝ}^n$ and ${\mathit{f}}_2\left(\mathit{r}\right):{ℝ}^N\mapsto {ℝ}^{N-n}$ .

Figure 2.

Chaotic attractor of the Lorenz system (equation (6)) under the parameters $\sigma =10$ , $\mathit{r}=28$ , and $b=8/3$ .

By considering the chaotic system (equation (5)) and the extended system (equation (8)), the state $\mathit{p}\left(t\right)$ can be computed by solving the next differential equation

$\dot{\mathit{p}}={\mathit{f}}_2\left(\mathit{r}\right),\mathit{p}\left(0\right)\in {ℝ}^{N-n}$

9

To clarify the previous development, let us consider a 2-DOF FJR manipulator. Then, the link position vector is $\mathit{q}={\left[q_1,q_2\right]}^{\text{T}}$ . By using the Lorenz system (equation (6)) and the extended vector $\mathit{r}={\left[q_1,q_2,p_1\right]}^{\text{T}}$ , we can write the vector field $\mathit{f}\left(\mathit{r}\right)$ from equation (8) as follows

$\mathit{f}\left(\mathit{r}\right)=\left[\begin{array}{c}{\mathit{f}}_1\left(\mathit{r}\right)\\ f_2\left(\mathit{r}\right)\end{array}\right]=\left[\begin{array}{c}\sigma \left[q_2-q_1\right]\\ r{q}_1-q_2-q_1{p}_1\\ q_1{q}_2-b{p}_1\end{array}\right]$

with

$\begin{array}{c}{\mathit{f}}_1\left(\mathit{r}\right)=\left[\begin{array}{c}\sigma \left[q_2-q_1\right]\\ r{q}_1-q_2-q_1{p}_1\end{array}\right]\\ f_2\left(\mathit{r}\right)=q_1{q}_2-b{p}_1\end{array}$

Note that the dimension of the extended vector $\mathit{r}\left(t\right)$ and the link position vector $\mathit{q}\left(t\right)$ are not necessarily distinct. If $n=N$ , then $\mathit{f}\left(\mathit{r}\right)={\mathit{f}}_1\left(\mathit{r}\right)$ and ${\mathit{f}}_2\left(\mathit{r}\right)=0$ . Besides, if $N>n$ , we are free to compose any ${\mathit{f}}_1\left(\mathit{r}\right)$ from $\mathit{f}\left(\mathit{r}\right)$ .

Problem formulation

The aim of this work is to show that by means of a model-based velocity controller a chaotic motion is induced to the link position vector of an n-DOF FJR manipulator. Then, this chaotification problem can be formulated as follows.

Control problem

Find a control law $\mathit{\tau }\left(t\right)$ for the system (equations (1) and (2)) such that

$\underset{t\to \infty }{{\displaystyle \text{lim}}}\hspace{0.17em}{\tilde{\mathit{f}}}_1\left(t\right)=0$

10

where

${\tilde{\mathit{f}}}_1\left(t\right)={\mathit{f}}_1\left(\mathit{r}\right)-\dot{\mathit{q}}\left(t\right)$

11

is the chaotization error.

Chaotification of the FJR

This section shows how to design a control law $\mathit{\tau }\left(t\right)$ such that the chaotification goal (equation (10)) is achieved. The block diagram of the proposed velocity controller, which has a nested feedback structure, is illustrated in Figure 3. Stage A shows an outer loop consisting of a link model-based velocity controller ${\tau }_{rc}\left(t\right)$ which guarantees that the difference between the link velocity $\dot{\mathit{q}}\left(t\right)$ and the desired chaotic reference field ${\mathit{f}}_1\left(t\right)$ converges to zero. Using this controller, a desired rotor position signal ${\varphi }_d\left(t\right)$ , together with its first and second time derivatives, is generated. Stage B contains an inner loop that, using the desired rotor position ${\varphi }_d\left(t\right)$ and its first two time derivatives, decouples the rotor dynamics (equation (2)) of the FJR manipulator from equation (1), and makes it globally exponentially stable. In the following development, the mathematical description and stability analysis of the proposed controller are presented.

Figure 3.

Schematic representation of the proposed chaotization scheme for an n-DOF FJR manipulator. DOF: degree of freedom; FJR: flexible joint robot.

Controller design

In order to decouple the rotor dynamics (equation (2)) from the system (equation (1)), let us consider the controller

$\mathit{\tau }=K\left[\varphi -\mathit{q}\right]+F_2\left(\dot{\varphi }\right)+K_{pf}\tilde{\varphi }+K_{df}\dot{\tilde{\varphi }}+M_2{\ddot{\varphi }}_d$

12

where $K_{pf},K_{df}\in {ℝ}^{n\times n}$ are diagonal positive definitive matrices, $\tilde{\varphi }={\varphi }_d-\varphi$ is the rotor position error, and ${\varphi }_d$ is the desired rotor position to be defined later.

By substituting equation (12) into equation (2), the next closed-loop equation is obtained

$M_2\ddot{\tilde{\varphi }}+K_{df}\dot{\varphi }+K_{pf}\tilde{\varphi }=0$

13

which can be written in state-space form as follows

$\begin{array}{l}\frac{\text{d}}{\text{d}t}\tilde{\varphi }=\dot{\tilde{\varphi }}\\ \frac{\text{d}}{\text{d}t}\dot{\tilde{\varphi }}=-M_2^{-1}\left[K_{df}\dot{\tilde{\varphi }}+K_{pf}\tilde{\varphi }\right]\end{array}$

14

Notice that the system (equation (14)) is linear with respect to variables $\tilde{\varphi }$ and $\dot{\tilde{\varphi }}$ . Then, by a proper selection of matrices $K_{pf}$ and $K_{df}$ , the state-space origin of system (equation (14)) can be rendered exponentially stable, which implies that $\varphi \left(t\right)$ exponentially converges to ${\varphi }_d\left(t\right)$ .

The next step in the control design consists of proving that, when $\varphi ={\varphi }_d$ , the link dynamics (equation (1)) has solutions that converge asymptotically to zero. To this end, let us define the desired rotor position ${\varphi }_d\left(t\right)$ as follows

${\varphi }_d=K^{-1}{\mathit{\tau }}_{rc}+\mathit{q}$

15

The controller ${\mathit{\tau }}_{rc}\left(t\right)$ in equation (15) is the model-based velocity controller that accomplishes the chaotification objective (equation (10)). The mathematical description of ${\mathit{\tau }}_{rc}\left(t\right)$ is given by

$\begin{array}{l}{\tau }_{rc}=M_1\left(\mathit{q}\right){\ddot{\mathit{q}}}_{\mathit{r}}+C_A\left(\mathit{q},\dot{\mathit{q}}\right){\dot{\mathit{q}}}_{\mathit{r}}+F_1\left(\dot{\mathit{q}}\right)+G\left(\mathit{q}\right)\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+K_{pr}\xi +K_{dr}{\tilde{\mathit{f}}}_1\end{array}$

16

where

$\begin{array}{l}{\dot{\mathit{q}}}_r={\mathit{f}}_1\left(\mathit{r}\right)+\Lambda \xi \\ \dot{\xi }={\tilde{\mathit{f}}}_1\end{array}$

17

and $K_{pr},K_{dr},\Lambda \in {ℝ}^{n\times n}$ are diagonal positive definite matrices.

The controller $\mathit{\tau }\left(t\right)$ in equation (12) requires to compute the signals ${\dot{\varphi }}_d\left(t\right)$ and ${\ddot{\varphi }}_d\left(t\right)$ . At first sight, it seems that link acceleration $\ddot{\mathit{q}}\left(t\right)$ and jerk ${\mathit{q}}^{\left(3\right)}\left(t\right)$ would be necessary. However, the structure of the rigid model (equation (1)) allows obtaining ${\dot{\varphi }}_d\left(t\right)$ and ${\ddot{\varphi }}_d\left(t\right)$ by writing the link acceleration and jerk as a function of the link position and velocity.

Note that

${\dot{\mathit{q}}}_r-\dot{\mathit{q}}={\mathit{f}}_1-\dot{\mathit{q}}+\Lambda \xi ={\tilde{\mathit{f}}}_1+\Lambda \xi$

Then, by using equations (11), (15), and (17), the rigid robot dynamics (equation (1)) can be written as follows

$\begin{array}{l}{M}_1\left[{\dot{\tilde{\mathit{f}}}}_1+\Lambda {\tilde{\mathit{f}}}_1\right]+C_A\left[{\tilde{\mathit{f}}}_1+\Lambda \xi \right]+K_{pr}\xi \\ +K_{dr}{\tilde{\mathit{f}}}_1=K\tilde{\varphi }\end{array}$

which yields

$\begin{array}{l}{\dot{\tilde{\mathit{f}}}}_1=-\Lambda {\tilde{\mathit{f}}}_1-M_1^{-1}[C_A\left(\mathit{q},\dot{\mathit{q}}\right)\left[{\tilde{\mathit{f}}}_1+\Lambda \xi \right]+K_{pr}\xi \hspace{0.17em}\\ +\hspace{0.17em}{K}_{dr}{\tilde{\mathit{f}}}_1]+M_1^{-1}K\tilde{\varphi }\end{array}$

18

Using equations (14), (17), and (18), the overall closed-loop system can be written as the next first-order system

$\frac{\text{d}}{\text{d}t}\left[\begin{array}{c}\xi \\ {\tilde{\mathit{f}}}_1\\ \tilde{\varphi }\\ \dot{\tilde{\varphi }}\end{array}\right]=\left[\begin{array}{c}{\tilde{\mathit{f}}}_1\\ {\dot{\tilde{\mathit{f}}}}_1\\ \dot{\tilde{\varphi }}\\ -M_2^{-1}\left[K_{df}\dot{\tilde{\varphi }}+K_{pf}\tilde{\varphi }\right]\end{array}\right]$

19

Note that ${\left[{\xi }^{\text{T}},{\tilde{\mathit{f}}}_1^{\text{T}},{\tilde{\varphi }}^{\text{T}},{\dot{\tilde{\varphi }}}^{\text{T}}\right]}^{\text{T}}=\mathbf{0}\in {ℝ}^{4n}$ is an equilibrium of the closed-loop system (equation (19)). Besides, by using equations (7), (8), and (11), the time derivative of the extended vector $\mathit{r}\left(t\right)$ can be written as

$\frac{\text{d}}{\text{d}t}\mathit{r}=\mathit{f}\left(\mathit{r}\right)-\left[\begin{array}{c}{\tilde{\mathit{f}}}_1\\ 0_{N-n}\end{array}\right]$

20

Thus, the asymptotic stability of the state-space origin of system (equation (19)) implies that the chaotification objective given in equation (10) is accomplished.

Stability analysis

To prove the stability of the origin of equation (19), let us consider the Lyapunov candidate function $V\left(\xi ,{\tilde{\mathit{f}}}_1,\tilde{\varphi },\dot{\tilde{\varphi }}\right)$ defined as

$\begin{array}{l}V=\frac{1}{2}{\left[{\tilde{\mathit{f}}}_1+\Lambda \xi \right]}^{\text{T}}{M}_1\left(\mathit{q}\right)\left[{\tilde{\mathit{f}}}_1+\Lambda \xi \right]\\ \hspace{0.17em}\hspace{0.17em}+\frac{1}{2}{\xi }^{\text{T}}\left[K_{pr}+\Lambda K_{dr}\right]\xi +\frac{{\alpha }_3}{2}{\mathit{z}}^{\text{T}}P\mathit{z}\end{array}$

21

where ${\alpha }_3$ is a strictly positive constant whose range of values will be obtained later

$\mathit{z}={\left[{\tilde{\varphi }}^{\text{T}},{\dot{\tilde{\varphi }}}^{\text{T}}\right]}^{\text{T}}$

22

and P is a $4\hspace{0.17em}\times \hspace{0.17em}4$ symmetric positive definite matrix and the unique solution of the Lyapunov equation

$\frac{1}{2}\left[A^{\text{T}}P+PA\right]=-Q_1$

23

with $Q_1\in {ℝ}^{4\times 4}$ a symmetric positive definite matrix, and

$A=\left[\begin{array}{cc}{0}_{2\times 2}& I_{2\times 2}\\ -M_2^{-1}{K}_{pf}& -M_2^{-1}{K}_{df}\end{array}\right]$

24

By virtue that $M_1\left(\mathit{q}\right)$ , $K_{pr}$ , $K_{dr}$ , $\Lambda$ , and P are symmetric and positive definite matrices, then V in equation (21) is a globally positive definite and radially unbounded function.

By considering that $Q_1=I_{2\times 2}$ , the time derivative of equation (21) along the trajectories of equation (19) becomes

$\begin{array}{l}\dot{V}={\left[{\tilde{\mathit{f}}}_1+\Lambda \xi \right]}^{\text{T}}{M}_1\left[{\dot{\tilde{\mathit{f}}}}_1+\Lambda \dot{\xi }\right]+{\xi }^{\text{T}}\left[K_{p\mathit{r}}+\Lambda K_{d\mathit{r}}\right]{\tilde{\mathit{f}}}_1\\ \hspace{0.17em}\hspace{0.17em}+\frac{1}{2}{\left[{\tilde{\mathit{f}}}_1+\Lambda \xi \right]}^{\text{T}}{\dot{M}}_1\left[{\tilde{\mathit{f}}}_1+\Lambda \xi \right]\\ \hspace{0.17em}\hspace{0.17em}+\frac{{\alpha }_3}{2}\left[{\dot{\mathit{z}}}^{\text{T}}P\mathit{z}+{\mathit{z}}^{\text{T}}P\dot{\mathit{z}}\right]\\ ={\left[{\tilde{\mathit{f}}}_1+\Lambda \xi \right]}^{\text{T}}[K\tilde{\varphi }-C_A\left[{\tilde{\mathit{f}}}_1+\Lambda \xi \right]-K_{p\mathit{r}}\xi \\ \hspace{0.17em}\hspace{0.17em}-K_{d\mathit{r}}{\tilde{\mathit{f}}}_1]+\frac{1}{2}{\left[{\tilde{\mathit{f}}}_1+\Lambda \xi \right]}^{\text{T}}{\dot{M}}_1\left[{\tilde{\mathit{f}}}_1+\Lambda \xi \right]\\ \hspace{0.17em}\hspace{0.17em}+{\xi }^{\text{T}}\left[K_{p\mathit{r}}+\Lambda K_{d\mathit{r}}\right]{\tilde{\mathit{f}}}_1+\frac{{\alpha }_3}{2}\left[{\dot{\mathit{z}}}^{\text{T}}P\mathit{z}+{\mathit{z}}^{\text{T}}P\dot{\mathit{z}}\right]\\ ={\tilde{\mathit{f}}}_1^{\text{T}}K\tilde{\varphi }-{\tilde{\mathit{f}}}_1^{\text{T}}{K}_{p\mathit{r}}\xi -{\tilde{\mathit{f}}}_1^{\text{T}}{K}_{d\mathit{r}}{\tilde{\mathit{f}}}_1+{\xi }^{\text{T}}\Lambda K\tilde{\varphi }\\ \hspace{0.17em}\hspace{0.17em}-{\xi }^{\text{T}}\Lambda K_{p\mathit{r}}\xi -{\xi }^{\text{T}}\Lambda K_{d\mathit{r}}{\tilde{\mathit{f}}}_1+{\tilde{\mathit{f}}}_1^{\text{T}}{K}_{p\mathit{r}}\xi \\ \hspace{0.17em}\hspace{0.17em}+{\xi }^{\text{T}}\Lambda K_{d\mathit{r}}{\tilde{\mathit{f}}}_1+\frac{{\alpha }_3}{2}\left[{\dot{\mathit{z}}}^{\text{T}}P\mathit{z}+{\mathit{z}}^{\text{T}}P\dot{\mathit{z}}\right]\\ =-{\tilde{\mathit{f}}}_1^{\text{T}}{K}_{d\mathit{r}}{\tilde{\mathit{f}}}_1-{\xi }^{\text{T}}\Lambda K_{p\mathit{r}}\xi +{\xi }^{\text{T}}\Lambda K\tilde{\varphi }+{\tilde{\mathit{f}}}_1^{\text{T}}K\tilde{\varphi }\\ \hspace{0.17em}\hspace{0.17em}+\frac{{\alpha }_3}{2}\left[{\dot{\mathit{z}}}^{\text{T}}P\mathit{z}+{\mathit{z}}^{\text{T}}P\dot{\mathit{z}}\right]\end{array}$

25

where the skew symmetry property that relates ${\dot{M}}_1$ and $C_A$ was used.

By using equations (19) and (22), the system

$\begin{array}{l}\dot{\mathit{z}}=\left[\begin{array}{c}\dot{\tilde{\varphi }}\\ \ddot{\tilde{\varphi }}\end{array}\right]=\left[\begin{array}{cc}{0}_{2\times 2}& I_{2\times 2}\\ -M_2^{-1}{K}_{p\mathit{f}}& -M_2^{-1}{K}_{d\mathit{f}}\end{array}\right]\left[\begin{array}{c}\tilde{\varphi }\\ \dot{\tilde{\varphi }}\end{array}\right]=A\mathit{z}\end{array}$

26

is obtained. Thus, substitution of equation (26) into equation (25) yields

$\begin{array}{l}\dot{V}=-{\tilde{\mathit{f}}}_1^{\text{T}}{K}_{d\mathit{r}}{\tilde{\mathit{f}}}_1-{\xi }^{\text{T}}\Lambda K_{p\mathit{r}}\xi +{\xi }^{\text{T}}\Lambda K\tilde{\varphi }+{\tilde{\mathit{f}}}_1^{\text{T}}K\tilde{\varphi }\\ +\frac{{\alpha }_3}{2}\left[{\mathit{z}}^{\text{T}}{A}^{\text{T}}P\mathit{z}+{\mathit{z}}^{\text{T}}PA\mathit{z}\right]\\ =-{\tilde{\mathit{f}}}_1^{\text{T}}{K}_{d\mathit{r}}{\tilde{\mathit{f}}}_1-{\xi }^{\text{T}}\Lambda K_{p\mathit{r}}\xi +{\xi }^{\text{T}}\Lambda K\tilde{\varphi }+{\tilde{\mathit{f}}}_1^{\text{T}}K\tilde{\varphi }\\ -{\alpha }_3{\mathit{z}}^{\text{T}}\mathit{z}\\ \le -{\lambda }_{\text{min}}\left(\Lambda K_{p\mathit{r}}\right){\Vert \xi \Vert }^2-{\lambda }_{\text{min}}\left(K_{d\mathit{r}}\right){\Vert {\tilde{\mathit{f}}}_1\Vert }^2\\ -{\alpha }_3\left[{\Vert \tilde{\varphi }\Vert }^2+{\Vert \dot{\tilde{\varphi }}\Vert }^2\right]\\ +{\lambda }_{\text{max}}\left(\Lambda K\right)\Vert \xi \Vert \Vert \tilde{\varphi }\Vert +{\lambda }_{\text{max}}\left(K\right)\Vert {\tilde{\mathit{f}}}_1\Vert \Vert \tilde{\varphi }\Vert \\ =-{\mathit{\eta }}^{\text{T}}{Q}_2\mathit{\eta }\end{array}$

27

where

$\begin{array}{l}\eta ={\left[\begin{array}{cccc}\Vert \xi \Vert & \Vert {\tilde{\mathit{f}}}_1\Vert & \Vert \tilde{\varphi }\Vert & \Vert \dot{\tilde{\varphi }}\Vert \end{array}\right]}^{\text{T}}\\ Q_2=\left[\begin{array}{cccc}{Q}_{11}& 0& Q_{13}& 0\\ 0& Q_{22}& Q_{23}& 0\\ Q_{31}& Q_{32}& {\alpha }_3& 0\\ 0& 0& 0& {\alpha }_3\end{array}\right]\end{array}$

with

$\begin{array}{l}{Q}_{11}={\lambda }_{\text{min}}\left(\Lambda K_{p\mathit{r}}\right)\\ Q_{13}=-\frac{{\lambda }_{\text{max}}\left(\Lambda K\right)}{2}\\ Q_{22}={\lambda }_{\text{min}}\left(K_{d\mathit{r}}\right)\\ Q_{23}=-\frac{{\lambda }_{\text{max}}\left(K\right)}{2}\\ Q_{31}=-\frac{{\lambda }_{\text{max}}\left(\Lambda K\right)}{2}\\ Q_{32}=-\frac{{\lambda }_{\text{max}}\left(K\right)}{2}\end{array}$

Let us define

$\begin{array}{l}\Delta ={\lambda }_{\text{min}}\left(\Lambda K_{p\mathit{r}}\right){\left[{\lambda }_{\text{max}}\left(K\right)\right]}^2\\ \hspace{0.17em}\hspace{0.17em}+{\lambda }_{\text{max}}\left(\Lambda K\right){\lambda }_{\text{min}}\left(K_{d\mathit{r}}\right){\lambda }_{\text{max}}\left(\Lambda K\right)\end{array}$

28

Then, a sufficient condition for the matrix Q₂ to be positive definite is that the inequality

$0<{\alpha }_3<\frac{\Delta }{{\lambda }_{\text{min}}\left(\Lambda K_{p\mathit{r}}\right){\lambda }_{\text{min}}\left(K_{d\mathit{r}}\right)}$

29

holds, which is always possible because α₃ can be chosen arbitrarily. Then, by assuming that condition (equation (29)) is fulfilled, it can be concluded that $\dot{V}$ is globally negative definite. Therefore, the limit

$\underset{t\to \infty }{{\displaystyle \text{lim}}}\left[\xi \left(t\right),{\tilde{\mathit{f}}}_1\left(t\right),\tilde{\varphi }\left(t\right),\dot{\tilde{\varphi }}\left(t\right)\right]=\mathbf{0}$

30

is satisfied for any initial condition.

The next proposition summarizes the previous development and provides the equations of the proposed model-based velocity controller for the chaotization of a general n-DOF FJR manipulator.

Proposition 1

Consider the system (equations (1) and (2)) in closed loop with the decoupling controller

$\mathit{\tau }=K\left[\varphi -\mathit{q}\right]+F_2\left(\dot{\varphi }\right)+K_{p\mathit{f}}\tilde{\varphi }+K_{d\mathit{f}}\dot{\tilde{\varphi }}+M_2{\ddot{\varphi }}_d$

31

where

${\varphi }_d=K^{-1}{\mathit{\tau }}_{\mathit{r}c}+\mathit{q}$

is the desired rotor position, and

$\begin{array}{l}{\mathit{\tau }}_{\mathit{r}c}=M_1\left(\mathit{q}\right){\ddot{\mathit{q}}}_{\mathit{r}}+C_A\left(\mathit{q},\dot{\mathit{q}}\right){\dot{\mathit{q}}}_{\mathit{r}}+F_1\left(\dot{\mathit{q}}\right)+G\left(\mathit{q}\right)\\ +K_{p\mathit{r}}\xi +K_{d\mathit{r}}{\tilde{\mathit{f}}}_1\end{array}$

32

describes the model-based velocity controller with

$\begin{array}{l}{\dot{\mathit{q}}}_{\mathit{r}}={\mathit{f}}_1\left(\mathit{r}\right)+\Lambda \xi \\ \dot{\xi }={\tilde{\mathit{f}}}_1\end{array}$

Then, the signals $\xi ,{\tilde{\mathit{f}}}_1,\tilde{\varphi },\dot{\tilde{\varphi }}$ converge globally asymptotically to zero, and the chaotification objective (equation (10)) is guaranteed.

Remark 2

Note from equation (14) that the gain matrices $K_{p\mathit{f}}$ and $K_{d\mathit{f}}$ are directly related to the aysmptotic convergence of the rotor position error $\tilde{\varphi }\left(t\right)$ . Hence, large values of $K_{p\mathit{f}}$ ensure a fast convergence of $\varphi \left(t\right)$ to ${\varphi }_d\left(t\right)$ . Besides, the oscillations can be regulated through the values of $K_{d\mathit{f}}$ . On the other hand, the gain matrices $K_{p\mathit{r}}$ and $K_{d\mathit{r}}$ are related to the model-based velocity controller ${\mathit{\tau }}_{\mathit{r}c}$ as indicated by equation (16) which in turn ensures that ${\tilde{\mathit{f}}}_1\left(t\right)$ and $\xi \left(t\right)$ converge to zero. Therefore, by increasing the values of $K_{p\mathit{r}}$ helps to reduce the magnitude of the chaotization error ${\tilde{\mathit{f}}}_1\left(t\right)$ , while reducing the oscillations through a proper selection of the values of $K_{d\mathit{r}}$ . Finally, let us define the variable $\mathit{s}={\tilde{\mathit{f}}}_1+\Lambda \xi$ . From the “Stability analysis” section, the asymptotic convergence of ${\tilde{\mathit{f}}}_1\left(t\right)$ and $\xi \left(t\right)$ was proved. Hence, by using equation (17), we conclude that the system solutions converge to the manifold $\mathbf{s}=0$ , which yields the dynamics ${\dot{\tilde{\mathit{f}}}}_1=-\Lambda {\tilde{\mathit{f}}}_1$ . Therefore, we are allowed to reduce the magnitude of the chaotization error by increasing the values of the elements of matrix Λ.

Remark 3

Before implementing the controller (equation (31)), it is worth noting that the complexity of its implementation may increase due to the presence of the terms ${\dot{\varphi }}_d\left(t\right)$ and ${\ddot{\varphi }}_d\left(t\right)$ . This stage can be complicated and time-consuming and becomes even more difficult when the number of DOF increases. In this work, the experimental results were obtained by computing analytically every signal. However, for more complex systems, a possibility to overcome such a complication may be by resorting to an arbitrary-order exact differentiator. ⁶⁵

Experimental results

This section illustrates the performance of the previous chaotization procedure using a 2-DOF FJR manipulator. First, the experimental setup which consists of a Quanser 2-DOF FJR manipulator is described. Then, the chaotic reference field used in the experiments is presented. A methodology for computing the rotor and link velocities is also described. Finally, the analysis of the experimental results is provided.

Platform description

The experiments were carried out using Matlab-Simulink R2012b and the Quanser 2-DOF FJR manipulator depicted in Figure 4. Simulink models were built using a fixed step size of 1 ms and a solver ode3 Bogacki–Shampine. In the robotic platform, the joints exhibit visible harmonics during accelerations, while the links are rigid in comparison. The system has two DC motors, the Maxon (Maxon motor, Fall River, Massachusetts, USA) 273759 and 118752 for joints 1 and 2, respectively. Each motor drives harmonic gearboxes with null backlash in a two-bar serial linkage. Both links are rigid, the primary is coupled to the first drive by means of a flexible joint and it carries at its end the second harmonic drive, which is coupled to the second rigid link via another flexible joint. Both motors and flexible joints are instrumented with quadrature optical encoders of 4096 pulses per revolution. Note that the robot moves on the horizontal plane, which implies a null effect of the gravity on it.

Figure 4.

Experimental platform of a 2-DOF FJR used in the experiments. DOF: degree of freedom; FJR: flexible joint robot.

The DC motors delivering torque at the joints are driven by servo amplifiers configured in current mode. ^15,66 Hence, the torques applied at the joints are related with the motor current as $\mathit{\tau }\left(t\right)=K_m{i}_m\left(t\right)$ , where $i_m\in {ℝ}^2$ is the motor current vector and $K_m\in {ℝ}^{2\times 2}$ is a positive definite diagonal matrix containing the motor torque constants. Then, the actual motor current becomes $i_m\left(t\right)=i_d\left(t\right)=K_{sa}\mathit{u}\left(t\right)$ , where $i_d\left(t\right)\in {ℝ}^2$ is the desired motor current vector, $K_{sa}\in {ℝ}^{2\times 2}$ is a diagonal matrix with the servo amplifier gains and $\mathit{u}\left(t\right)\in {ℝ}^2$ is the servo amplifier input voltage. Finally, the DC motor output torques can be written as follows

$\mathit{\tau }\left(t\right)=K_m{K}_{sa}\mathit{u}\left(t\right)$

33

indicating that the delivered torque is proportional to the input servo amplifier voltage.

In the experimental platform, the motor torque constant matrix K_m and servo amplifier gain matrix $K_{sa}$ have the values

$\begin{array}{l}{K}_m=\left[\begin{array}{cc}0.119& 0\\ 0& 0.0234\end{array}\right]\hspace{0.17em}\hspace{0.17em}[\text{Nm/A}]\\ K_{sa}=\left[\begin{array}{cc}0.5& 0\\ 0& 0.5\end{array}\right]\hspace{0.17em}\hspace{0.17em}[\text{A/V}]\end{array}$

respectively. Besides, the maximum torques delivered by motors 1 and 2 are $59.5 \hspace{0.17em}[\text{Nm}]$ and $5.8 \hspace{0.17em}[\text{Nm}]$ , respectively.

Regarding the dynamic model equations (1) and (2), the corresponding matrices for the experimental platform are ^59,60

$M_1\left(\mathit{q}\right)=\left[\begin{array}{cc}{\theta }_1+2{\theta }_2\text{cos}\left(q_2\right)& {\theta }_2\text{cos}\left(q_2\right)+{\theta }_3\\ {\theta }_2\text{cos}\left(q_2\right)+{\theta }_3& {\theta }_3\end{array}\right]$

$M_2=\left[\begin{array}{cc}{\theta }_4& 0\\ 0& {\theta }_5\end{array}\right]$

$C_A\left(\mathit{q},\dot{\mathit{q}}\right)=\left[\begin{array}{cc}-{\theta }_2{\dot{q}}_2\text{sin}\left(q_2\right)& -{\theta }_2\left({\dot{q}}_1+{\dot{q}}_2\right)\text{sin}\left(q_2\right)\\ {\theta }_2{\dot{q}}_1\text{sin}\left(q_2\right)& 0\end{array}\right]$

$\begin{array}{l}{F}_1\left(\dot{\mathit{q}}\right)=F_{v1}\dot{\mathit{q}}+{\mathit{f}}_{c1}\left(\dot{\mathit{q}}\right),F_2\left(\dot{\varphi }\right)\\ =F_{v2}\dot{\varphi }+{\mathit{f}}_{c2}\left(\dot{\varphi }\right)\end{array}$

$F_{v1}=\left[\begin{array}{cc}{\theta }_6& 0\\ 0& {\theta }_7\end{array}\right],F_{v2}=\left[\begin{array}{cc}{\theta }_8& 0\\ 0& {\theta }_9\end{array}\right]$

$K=\left[\begin{array}{cc}{\theta }_{10}& 0\\ 0& {\theta }_{11}\end{array}\right]$

${\mathit{f}}_{c1}\left(\dot{\mathit{q}}\right)=\left[\begin{array}{c}{\theta }_{12}\text{tanh}\left(\beta {\dot{q}}_1\right)\\ {\theta }_{13}\text{tanh}\left(\beta {\dot{q}}_2\right)\end{array}\right]$

${\mathit{f}}_{c2}\left(\dot{\varphi }\right)=\left[\begin{array}{c}{\theta }_{14}\text{tanh}\left(\beta {\dot{\varphi }}_1\right)\\ {\theta }_{15}\text{tanh}\left(\beta {\dot{\varphi }}_2\right)\end{array}\right]$

$\mathit{\tau }={\left[\begin{array}{cc}{\tau }_1& {\tau }_2\end{array}\right]}^{\text{T}},\mathit{\theta }={\left[\begin{array}{ccc}{\theta }_1& \dots & {\theta }_{15}\end{array}\right]}^{\text{T}}$

The definitions of parameters θ_i are given in Table 1 using the following notation

m_li     mass of link i,

m _r1    mass of rotor i,

L ₁    length of link 1,

l_i     distance to the center of gravity of link i,

I_lzzi     principal inertia moment at z axis of link i,

I_rzzi     principal inertia moment at z axis of rotor i,

r_i     gear ratio of actuator i,

k _s1    first FJR torsional stiffness,

k _s2    second FJR torsional stiffness,

f ₁    viscous friction coefficient of link 1,

f ₂    viscous friction coefficient of link 2,

f ₃    viscous friction coefficient of joint 1,

f ₄    viscous friction coefficient of joint 2,

f _c1    Coulomb friction coefficient of link 1,

f _c2    Coulomb friction coefficient of link 2.

f _c3    Coulomb friction coefficient of joint 1,

f _c4    Coulomb friction coefficient of joint 2,

with $i=1,2$ .

Table 1.

Parameters ${\theta }_i,i=1,\dots \mathrm{,15}$ , for the 2-DOF FJR.

Parameter	Definition	Estimated value	Unit
${\theta }_1$	$\left(m_{l1}+m_{\mathit{r}1}\right)l_1^2+I_{l\mathit{z}\mathit{z}1}+\left(m_{l2}+m_{\mathit{r}2}\right)\left(L_1^2+l_2^2\right)+I_{l\mathit{z}\mathit{z}2}$	0.207184	Kg m2
θ 2	$\left(m_{l2}+m_{\mathit{r}2}\right)L_1{l}_2$	0.017580	Kg m2
θ 3	$\left(m_{l2}+m_{\mathit{r}2}\right)l_2^2+I_{l\mathit{z}\mathit{z}2}$	0.013163	Kg m2
θ 4	$r_1^2{I}_{rzz1}$	0.216776	Kg m2
θ 5	$r_2^2{I}_{rzz2}$	0.006842	Kg m2
θ 6	f 1	0.037675	Nms/rad
θ 7	f 2	0.002959	Nms/rad
θ 8	f 3	8.435283	Nms/rad
θ 9	f 4	0.135564	Nms/rad
θ 10	$k_{s1}$	9.358730	Nm/rad
θ 11	$k_{s2}$	4.212811	Nm/rad
θ 12	$f_{c1}$	0.015771	Nms
θ 13	$f_{c2}$	0.005018	Nms
θ 14	$f_{c3}$	1.360210	Nms
θ 15	$f_{c4}$	0.162529	Nms

DOF: degree of freedom; FJR: flexible joint robot.

Note that not all the system parameters θ_i are provided by the robot manufacturer. ⁶⁰ Furthermore, during robot operation, some of these values may vary, which makes it necessary to include a parameter identification stage. ⁶⁷ In this work, a dynamic filtering of the robot model, together with the application of a least squares parameter identification method, was used to estimate the values of the parameter vector θ . A detailed procedure description of this parameter identification methodology is given in the literature. ⁶⁰ The numerical values of the estimated parameters for the 2-DOF FJR manipulator are shown in Table 1.

Chaotic reference field

The velocity reference field $\mathit{f}\left(\mathit{r}\right)$ given in equation (8) was implemented using the Genesio–Tesi system described by ⁶⁸

$\frac{\text{d}}{\text{d}t}\left[\begin{array}{c}{x}_1\\ x_2\\ x_3\end{array}\right]=\left[\begin{array}{c}{x}_2\\ x_3\\ x_1^2+a_1{x}_1+a_2{x}_2+a_3{x}_3\end{array}\right]$

34

where $a_1,\hspace{0.17em}{a}_2,\hspace{0.17em}\text{and}\hspace{0.17em}{a}_3$ are positive constants. Due to the mechanical restrictions in amplitude of our experimental platform, we considered as an output of the system (equation (34)) the function

$y={\left[a_4{x}_1,\hspace{0.17em}{a}_5{x}_2,\hspace{0.17em}{a}_6{x}_3\right]}^{\text{T}}$

35

with $a_4,\hspace{0.17em}{a}_5,\hspace{0.17em}\text{and}\hspace{0.17em}{a}_6$ being positive constants that allow reducing the amplitude of the corresponding signals $x_1,\hspace{0.17em}{x}_2,\hspace{0.17em}\text{and}\hspace{0.17em}{x}_3$ . Then, by using the extended vector $\mathit{r}={\left[{\mathit{q}}^{\text{T}},p_1\right]}^{\text{T}}$ and equations (34) and (35), we can write the velocity reference field $\mathit{f}\left(\mathit{r}\right)$ as follows

$\begin{array}{l}\mathit{f}\left(\mathit{r}\right)=\left[\begin{array}{c}{\mathit{f}}_1\left(\mathit{r}\right)\\ f_2\left(\mathit{r}\right)\end{array}\right]\\ =a_0\left[\begin{array}{c}{a}_4{q}_2\\ a_5{p}_1\\ a_6\left[q_1^2+a_1{q}_1+a_2{q}_2+a_3{p}_1\right]\end{array}\right]\end{array}$

36

where ${\mathit{f}}_1\left(\mathit{r}\right)\in {ℝ}^2$ , $f_2\left(\mathit{r}\right)\in ℝ$ , $a_0>0$ is a constant used to speed up the velocity, $\mathit{q}\in {ℝ}^2$ are the FJR link positions.

By using equation (34), the coordinate p₁ can be obtained by solving the next differential equation

$\frac{\text{d}}{\text{d}t}{p}_1=q_1^2+a_1{q}_1+a_2{q}_2+a_3{p}_1$

37

Velocity estimation

For a given real robotic system, it is usual to have only encoders for position measurement. Then, velocity values must be estimated using the encoder measurements. This is the case for our experimental robot. The velocity values $\dot{\mathit{q}}\left(t\right)\hspace{0.17em}\text{and}\hspace{0.17em}\dot{\varphi }\left(t\right)$ were computed as follows. Consider the polynomials ^60,69

$\begin{array}{l}{r}_1\left(z^{-1}\right)=1+z^{-1}+z^{-2}+\dots +z^{-{\iota }_1}\\ r_2\left(z^{-1}\right)=1+z^{-1}+z^{-2}+\dots +z^{-{\iota }_2}\end{array}$

where ${\iota }_1,{\iota }_2\in \mathbb{N}$ . The average $\overline{\mathit{w}}\left(k{T}_s\right)$ of the last ${\iota }_3+1$ position samples at time instant $k{T}_s$ is given by

$\overline{\mathit{w}}{\left(kT\right)}_s=\frac{r_1\left(z^{-1}\right)}{{\iota }_3+1}\mathit{w}\left(k{T}_s\right)$

where k represents the time index, T_s is the sampling period, and w represents any of the signals $\mathit{q}\left(t\right)\hspace{0.17em}\text{and}\hspace{0.17em}\varphi \left(t\right)$ . Then, the discrete dirty derivative of $\overline{\mathit{w}}{\left(k{T}_s\right)}_{\hspace{0.17em}}$ can be computed as follows

$\dot{\overline{\mathit{w}}}\left(k{T}_s\right)=\frac{1-z^{-1}}{T_s}\overline{\mathit{w}}\left(k{T}_s\right)$

Finally, the velocity values were estimated by averaging the last ${\iota }_2+1$ samples of $\dot{\overline{\mathit{w}}}\left(k{T}_s\right)$ , that is

$\dot{\mathit{w}}\left(k{T}_s\right)\approx \frac{r_2\left(z^{-1}\right)}{{\iota }_2+1}\dot{\overline{\mathit{w}}}\left(k{T}_s\right)$

38

In this document, the velocity estimates were obtained using equation (38) with ${\iota }_1={\iota }_2=4$ .

Results and discussion

The experiments with the 2-DOF FJR were achieved by using

$\begin{array}{l}{a}_0=1, \hspace{0.17em}{a}_1=-1, \hspace{0.17em}{a}_2=-1.1, \hspace{0.17em}{a}_3=-0.45,\\ a_4=0.7, \hspace{0.17em}{a}_5=0.5, \hspace{0.17em}{a}_6=1\end{array}$

in the velocity reference field $\mathit{f}\left(r\right)$ in equation (36). The matrices of the controller (equation (31)) were set to

$\begin{array}{l}{K}_{pr}=\left[\begin{array}{cc}0.1& 0\\ 0& 5\end{array}\right],\hspace{0.17em}{K}_{dr}=\left[\begin{array}{cc}1& 0\\ 0& 1\end{array}\right],\\ K_{p\mathit{f}}=\left[\begin{array}{cc}1500& 0\\ 0& 5\end{array}\right],\hspace{0.17em}{K}_{d\mathit{f}}=\left[\begin{array}{cc}50& 0\\ 0& 0.1\end{array}\right],\\ \Lambda =\left[\begin{array}{cc}0.1& 0\\ 0& 0.1\end{array}\right]\end{array}$

In Figure 5, the time evolution of the torques applied at joints 1 and 2 is depicted. The control signals corresponding to each joint have high-frequency components. This behavior can be attributed to the fluctuation due to the pulse width modulation switching of the servo amplifier, quantized measurements of the angular position of the joints and the links, and the velocity estimation algorithm used to compute $\dot{\varphi }\left(t\right)$ and $\dot{\mathit{q}}\left(t\right)$ , introducing noise directly to the applied control action.

Figure 5.

Applied torques ${\tau }_1\left(t\right)$ and ${\tau }_2\left(t\right)$ for the chaotization experiment.

Using the measurements of the link positions $q_1\left(t\right)\hspace{0.17em}\text{and}\hspace{0.17em}{q}_2\left(t\right)$ , the plot of the chaotic attractor was obtained, which is depicted in the left-hand side of Figure 6. Similarly, using the values corresponding to p₁ obtained from equation (37) and the link positions allows obtaining the chaotic attractor shown in the right-hand side of Figure 6.

Figure 6.

Chaotic attractor: view using the coordinates $q_1\left(t\right)\hspace{0.17em}\text{and}\hspace{0.17em}{q}_2\left(t\right)$ (left); and view using the coordinates $p_1\left(t\right),\hspace{0.17em}{q}_1\left(t\right),\hspace{0.17em}\text{and}\hspace{0.17em}{q}_2\left(t\right)$ (right).

From the experimental data corresponding to the time evolution of q(t) and $p_1\left(t\right)$ , the Poincaré maps for the planes $q_1-p_1$ at $q_2=0$ , $q_2-p_1$ at $q_1=0$ , and $q_1-q_2$ at $p_1=0$ were obtained as shown in Figure 7.

Figure 7.

Poincaré maps generated for the planes: $q_1-p_1$ at $q_2=0$ (top), $q_2-p_1$ at $q_1=0$ (middle), and $q_1-q_2$ at $p_1=0$ (bottom).

The plot of the Poincaré maps is included as a numerical test, indicating that the robot motion is chaotic. Distinct set of points in the Poincaré map indicate the existence of chaos, ^70,71 which is the case for the maps given in Figure 7. In addition, a way to determine if a system is chaotic consists of computing the Lyapunov exponents. ⁷² If the largest Lyapunov exponent is positive, then the associated system dynamics has a chaotic behavior. From the experimental data, the largest Lyapunov exponent was ${\lambda }_e=0.252$ . Hence, we can conclude that the robot link motion corresponding to the FJR manipulator is chaotic. Finally, the signals corresponding to the components of the chaotization error vector ${\tilde{\mathit{f}}}_1\left(t\right)$ are depicted in Figure 8. We observe the existence of a small error, which may be due to parametric error, velocity estimation error or disturbances. Although this error can be reduced by increasing the controller gains (see remark 2), it does not affect the chaotization procedure as corroborated by the Poincaré maps and the value of the Lyapunov exponents.

Figure 8.

Behavior of the components of the chaotization error vector ${\tilde{\mathit{f}}}_1\left(t\right)$ in the chaotization experiment. From equation (36), the solid line corresponds to the term $a_0{a}_4{q}_2-{\dot{q}}_1$ and the dotted line represents the element $a_0{a}_5{p}_1-{\dot{q}}_2$ .

The experimental results indicate that the decoupling controller $\mathit{\tau }\left(t\right)$ described by equation (31) together with the model-based velocity controller ${\mathit{\tau }}_{rc}\left(t\right)$ defined in equation (32) accomplish the chaotification objective (equation (10)), inducing a chaotic behavior to the links of the FJR manipulator. Furthermore, this behavior is achieved despite the disturbances affecting the system such as position quantization error, error in velocity estimation, and parametric uncertainty, which demonstrates the robustness of the proposed controller and chaotification technique.

Many future studies can be made in the chaotization field. For instance, for parameter identification of robot manipulators, it is essential to consider the spectral richness (SR) of the reference signal and the trajectories followed by the system. If the SR is high enough, it improves the accuracy of the parameter identification technique. This methodology may be approached by considering that the robot links follow a chaotic reference. Some interesting works dealing with parameter identification are covered by Boubaker and Iriarte ⁷³ and Lassoued and Boubaker. ⁷⁴ In addition, the improvements in the chaotization of FJR manipulators when using robust controllers can be verified. ⁷³ Finally, other studies related to real-world applications of chaotic systems can be found in Boubaker and Jafary. ⁷⁵

Conclusion

This article proposed a model-based velocity controller that allows achieving the chaotization of an n-DOF FJR manipulator. The methodology employed a chaotic reference field, namely a Genesio–Tesi system. A Lyapunov-based analysis demonstrated that the robotic system in closed loop with the proposed controller is asymptotically stable, guaranteeing the chaotization objective. The proposed controller was experimentally tested on a Quanser 2-DOF FJR manipulator. The experimental data were used to compute the Poincaré maps and the largest Lyapunov exponent. Because these maps presented distinct set of points and the largest Lyapunov exponent was positive, it was concluded that the robot links exhibited a chaotic motion, validating the effectiveness and robustness of the proposed controller. As a future work, robust control techniques will be used to induce a chaotic behavior to the FJR manipulator affected by matched and unmatched disturbances. Besides, we will make use of the chaotization procedure to verify improvements in parameter identification.