Robust Control of Collaborative Manipulators - Flexible Object System

Abstract

Abstract In many manufacturing and automobile industries, flexible components need to be positioned with the help of coordinated operations of manipulators. This paper deals with the robust design of a control system for two planar rigid manipulators moving a flexible object in the prescribed trajectory while suppressing the vibration of the flexible object. Dynamic equations of the flexible object are derived using the Hamiltonian principle, which is expressed as a partial differential equation (PDE) with appropriate boundary conditions. Then, a combined dynamics is formulated by combining the manipulators and object dynamics without any approximation. The resulting dynamics are thus described by the PDEs, having rigid as well as flexible parameters coupled together. This paper attempts to develop a robust control scheme without approximating the PDE in order to avoid measurements of flexible coordinates and their time derivatives. For this purpose, the two subsystems, namely slow and fast subsystems, are identified by using the singular perturbation technique. Specific robust controllers for both the subsystems are developed. In general, usage of the singular perturbation technique necessitates exponential stability of both subsystems, which is evaluated by satisfying Tikhnov's theorem. Hence, the exponential stability analysis is performed for both subsystems. Focusing on two three-link manipulators holding a flexible beam, simulations are performed and simulation results demonstrate the versatility of the proposed robust composite control scheme.

1. Introduction

The past four decades saw the tremendous growth of the use of robots in automobile, electronic, construction, manufacturing and medical applications. Robotic manipulator systems relieved humans of boring repetitive tasks, tasks in dangerous environments such as space and underwater, and tasks in high radiation environments. In some situations a single robot may not be able to grasp and move slender flexible objects along a prescribed trajectory in a safe and efficient way. Owing to the single arm structure, present day robots are called handicapped operators for performing complex tasks. Most tasks in assembly/disassembly, handling large or heavy objects are done efficiently with the help of two robot arms. Collaborative manipulators have the following advantages compared to single arm manipulators: increased load carrying capacity, reduced need for extra auxiliary equipment, efficient use of available workspace and safer handling of delicate and flexible objects. Considering these facts, two manipulator systems are employed in a wide range of tasks and the first master/slave teleoperated manipulator was used in the nuclear industry in the 1940s at Oak Ridge and Argonne National Laboratories [1].

In modern automobile body assembly there are more than 200 sheet metal parts which must be assembled in a precise way. Handling them needs special equipment and skilled operations. Two collaboratively operated robots can grasp flexible pieces of sheet metal and manipulate them together for assembly. Furthermore, in industry many deformable objects, such as rubber tubes, sheet metal, cords and leather products, are handled by special equipment or human operators. Satellite solar arrays and the aerospace industry have composite materials with allowable flexibility that must be handled with the help of at least two robotic manipulators. In order to effectively manipulate complex flexible objects and bend them in a desired manner, for example, in the insertion of a flexible beam into a hole and the assembling of sheet metal in the required place, at least two robot arms are needed. Robot assisted surgery is another area where collaborative robots in a vibration free environment will be useful.

Early studies on the collaboration of two arms handling objects mainly focused on rigid bodies. Some of them can be found in [2] and [3] and their references. There are some studies on the handling of flexible objects in the literature. Two manipulators handling flexible objects is a complex problem compared with manipulating rigid objects. In particular, the dynamics and control of flexible objects are more tedious than that for rigid objects. In 1991, Zheng et al. [4] studied inserting flexible cables or pegs into a small hole by using a single arm robot. Dellinger and Anderson [5] developed a mathematical model for the case of two manipulators handling a non-rigid object and they also analysed the interactive forces between the manipulators and the object. Zheng and Chen [6] extended the work of [4] on single arm manipulation to two arm manipulation for the alignment of flexible sheets in printed circuit boards. Yukawa et al. [7] proposed a position control algorithm to achieve total system stability and they also suppressed the vibration at the intermediate points of the flexible beam. Kosuge et al. [8] employed finite element model of the sheet metal and also implemented their control algorithm experimentally to reduce the vibration of the object. Sun and Liu [9] and [10] developed a new approach of modelling the beam using assumed modes and proposed a hybrid impedance controller which is used to stabilize the system while suppressing the vibrations and controlling internal forces. Azadi et al. [11] considered robust control of two 5 DOF robot manipualtors. Caccavale et al. [14] studied the impedance behaviour of dual arm cooperative manipulators. Gueaieb [12] et al. developed hybrid control and recently Moosavian et al. [16] and Yagiz [17] et al. considered nonlinear control for the cooperative manipualtor scenario. Al-Yahmali and Hsia [13] derived the sliding mode control algorithm to provide robustness against the model imperfections and uncertainty. The stability of the system was proved and the manipulators were able to move the object in the desired trajectory while suppressing the vibration. Ali et al. [15] used two time scale controllers to track the desired trajectory for the rigid body motion and also to suppress the vibration. They implemented the PD control scheme for the rigid subsystem and after linearizing the flexible dynamics, they used the pole placement technique and the designed linear observer was used to estimate the vibration parameters.

It is important to note that in all the above studies on the handling of flexible objects, the dynamics of the flexible object is derived using the finite element method or assumed mode method. The truncation of the original model with infinite degrees of freedom of a flexible beam to a finite dimensional model leads to the following problems needing resolution [18]:

Requirement of as many sensors as the locations of the measurement of vibration and the difficulty in implementation.

Presence of control and observation spill over due to the ignored high frequency dynamics.

It is often not clear how many modes must be considered to approximate the PDE-based model into the ordinary differential equation (ODE) model.

Destabilization of the system due to the negligence of the higher order modes.

Necessity of higher order controller.

The above mentioned disadvantages may be avoided by using the PDE-based system. As a matter of fact, there are some results available for single-link flexible manipulators where the PDE is adopted without the need of discretization or approximation, see, for example, [18] and [19]. In these schemes, the strain signal is used instead of the measurements of flexible coordinates and their time derivatives. Quite recently, [20] and [21] also considered the PDE in their controller designs. In the present study the dynamics for the flexible object moved by two planar rigid manipulators will be described by the PDE, and the controller will be designed without using the assumed mode method. For this purpose, the dynamic equations of the object are derived using the Hamiltonian principle with appropriate boundary conditions. Then, a combined dynamics is formulated by combining the manipulators and object dynamics. The resulting system dynamics are thus described by the PDE, having rigid as well as flexible parameters coupled together.

In order to control the manipulators-flexible object system without using any approximation methods, a suitable control approach is the two-time scale theory which considers the high frequency phenomenon of flexible motion in different time scales. The basic idea for the two-time scale theory is to identify the slow and fast subsystems in separate time scales. Then, one can design a control algorithm for each subsystem and they together form a composite control input to achieve the desired objective. Therefore, the two subsystems, namely the slow system describing the rigid body motion and the fast subsystems describing the transverse vibration, are identified by using the singular perturbation technique. Focusing on these two subsystems, specific controllers are developed. The use of the singular perturbation technique necessitates the exponential stability of each subsystem which satisfies Tikhnov's theorem for the closed loop stability of the system which are presented here. The proposed robust controller has the following features: 1) it does not need the information of the modes; 2) exact knowledge of the manipulators parameters or the beam are not required; and 3) it satisfies Tikhnov's theorem. It should be noted that the singular perturbation technique has also been applied to the singlelink flexible beams, say, for example, in [20] and [21], where a slow system and fast subsystems were identified as well. Comparing with the case of the manipulators handling a flexible beam, the major difference lies in the so-called slow system, where it couples with manipulator dynamics. In this case, the slow motion controller design is completely different, especially for the case of unknown parameters, similar to the case of manipulators handling a rigid body, though the structure and properties of the system are not altered compared to a single manipulator.

The content of the paper is organized as follows. Section II presents the kinematics and dynamics of the beam without approximating or discretizing the beam. Section III deals with the manipulators' dynamics and combined dynamics resulting from the amalgamation of beam and manipulator dynamics. A singularly perturbed complete system of the dynamic model has been developed in section IV. Then, by incorporating the typical steps of the singular perturbation approach, the slow and fast subsystems are derived in section V, so that, one can design a controller for each subsystem. In section VI, a composite control strategy and corresponding stability analysis are carried out. For the control of slow subsystem, a regressor-based sliding mode control algorithm is developed. The main advantages of the proposed robust control scheme is that it can handle rapidly varying parameters and alleviate the problem of choosing the upper bounds of the uncertainties in the robust approaches; it also does not need persistency of excitation and guarantees the exponential convergence of transient behaviour. In the case of control of a fast subsystem, a dissipator operator forms the closed loop system and it damps out the vibration of the flexible object. Exponential stability analysis of slow and fast subsystem results validate the composite control strategy with the help of Tikhnov's theorem. Simulations are performed by considering two three-link manipulators rigidly grasping and moving a flexible beam in section VII. Simulation results demonstrate that the proposed control approach can achieve satisfactory tracking performance while suppressing the vibration of the flexible object.

2. Kinematics and dynamics of the beam

In various applications such as turbine rotor blades, spacecrafts with flexible appendages, flexible robot arms and aerospace systems, flexible objects can be modelled as beams. In the following analysis, a flexible object will be modelled as an Euler-Bernoulli beam. Since the two manipulators are used to move the object in the desired position and orientation which necessitates the allowance of rotation at the two ends of the beam, simply supported end conditions are considered in deriving the dynamic equations of motion of the beam.

A. Kinematics of the Beam

The co-ordinate frames X₁Y₁ and X₂Y₂ shown in Fig. 1 are fixed frames, and the xy-frame is a moving coordinate frame which is attached to the beam at its midpoint. Consider a beam of length L and mass m = ρL, where ρ is mass per unit length. The mass centre position and orientation of the beam with respect to X₁Y₁-frame are given by{c₀} = {x₀,y₀,θ}^T. The Z and z directions are perpendicular to the plane. F_1x,F_1y,F_2x,F_2y are the forces applied by the manipulator at the two ends of the beam. The transverse displacement η(x,t) is measured with respect to xy-frame and deformation in the longitudinal direction is neglected. All kinematic relations are written with respect to the X₁Y₁ frame, as shown in Fig. 2. Any point on the beam can be written as,

Figure 1.

Two planar rigid manipulators grasping a flexible object

Figure 2.

Flexible beam with coordinate frames

$\mathrm{X}={\mathrm{x}}_0+\mathrm{x}cos\theta -\eta sin\theta$

(1)

$\mathrm{Y}={\mathrm{y}}_0+\mathrm{x}sin\theta +\eta cos\theta$

(2)

where x is the spatial co-ordinate ranging from $\frac{-\mathrm{L}}{2}$ to $\frac{+\mathrm{L}}{2}$ and η is the transverse displacement.

The beam has rigid body motion on which the flexible motion or vibration of the beam is superimposed. It is evident from Fig. 1 that the left and right end point of the beam shares the left and right end-effector grasping point of the two end-effectors. Assuming that the slope due to transverse deflection is negligible compared to the orientation of the beam, the slopes at the two ends of the beam are neglected.

The following relations are obtained from Fig. 2.

The left end pose (position and orientation) of the beam is given by

$\begin{array}{l}\left\{{\mathrm{e}}_1\right\}=\left\{{\mathrm{c}}_{\mathrm{o}}\right\}-{\left\{\begin{array}{ccc}\frac{\mathrm{L}}{2}cos\theta & \frac{\mathrm{L}}{2}sin\theta & 0\end{array}\right\}}^{\mathrm{T}}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+{\left\{\begin{array}{ccc}-\eta \sin \theta & \eta \cos \theta & 0\end{array}\right\}}^{\mathrm{T}}\end{array}$

(3)

and the right end pose of the beam is given by

$\begin{array}{l}\left\{{\mathrm{e}}_2\right\}=\left\{{\mathrm{c}}_{\mathrm{o}}\right\}-{\left\{\begin{array}{ccc}\frac{\mathrm{L}}{2}cos\theta & \frac{\mathrm{L}}{2}sin\theta & 0\end{array}\right\}}^{\mathrm{T}}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+{\left\{\begin{array}{ccc}-\eta \sin \theta & \eta \cos \theta & 0\end{array}\right\}}^{\mathrm{T}}\end{array}$

(4)

Differentiating (3) and (4) results in

$\begin{array}{l}\left\{{\stackrel{.}{e}}_1\right\}=\left\{\begin{array}{c}{\stackrel{.}{x}}_0+\frac{\mathrm{L}}{2}sin\theta \stackrel{.}{\theta }-\eta \cos \theta \stackrel{.}{\theta }-\stackrel{.}{\eta }\sin \theta \\ {\stackrel{.}{y}}_0+\frac{\mathrm{L}}{2}cos\theta \stackrel{.}{\theta }-\eta \sin \theta \stackrel{.}{\theta }-\stackrel{.}{\eta }\cos \theta \\ \stackrel{.}{\theta }\end{array}\right\}\\ \left\{{\stackrel{.}{e}}_2\right\}=\left\{\begin{array}{c}{\stackrel{.}{x}}_0+\frac{\mathrm{L}}{2}sin\theta \stackrel{.}{\theta }-\eta \cos \theta \stackrel{.}{\theta }-\stackrel{.}{\eta }\sin \theta \\ {\stackrel{.}{y}}_0+\frac{\mathrm{L}}{2}cos\theta \stackrel{.}{\theta }-\eta \sin \theta \stackrel{.}{\theta }-\stackrel{.}{\eta }\cos \theta \\ \stackrel{.}{\theta }\end{array}\right\}\end{array}$

Above relations can be written in compact form with respect to the Cartesian co-ordinates as

$\left\{\stackrel{.}{e}\right\}=\left[\mathrm{R}\right]\left\{{\stackrel{.}{X}}_{rf}\right\}$

(5)

where ${\stackrel{.}{X}}_{rf}\stackrel{△}{=}{\left\{\begin{array}{ccc}{\stackrel{.}{x}}_0& {\stackrel{.}{y}}_0& \stackrel{.}{\theta }\end{array}\right\}}^{\mathrm{T}}$ denotes the velocity of the mass centre of the object, and

$\mathrm{R}=\left[\begin{array}{ccc}1& 0& \frac{\mathrm{L}}{2}sin\theta -\eta \cos \theta \\ 0& 1& -\frac{\mathrm{L}}{2}cos\theta -\eta \sin \theta \\ 0& 0& 1\\ 1& 0& -\frac{\mathrm{L}}{2}sin\theta -\eta \cos \theta \\ 0& 1& \frac{\mathrm{L}}{2}cos\theta -\eta \sin \theta \\ 0& 0& 1\end{array}\right]$

Differentiating (5) gives the acceleration

$\left\{\stackrel{..}{e}\right\}=\left[\stackrel{.}{R}\right]\left\{{\stackrel{.}{X}}_{rf}\right\}+\left[\mathrm{R}\right]\left\{{\stackrel{..}{X}}_{rf}\right\}$

(6)

where ${\stackrel{..}{X}}_{rf}$ describes the acceleration of the mass centre of the object.

B. Dynamics of the Beam

In order to determine the dynamic equations of motion of the beam, the kinetic energy, potential energies due to elasticity of the beam and due to gravity must be obtained. Then, by applying the extended Hamilton principle [22], the dynamic equations of the beam will be derived.

Kinetic energy of the beam is defined as

$\mathrm{T}=\frac{1}{2}{\displaystyle {\int }_{\frac{-\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\rho \left({\stackrel{.}{X}}^2+{\stackrel{.}{Y}}^2\right)dx}$

(7)

Differentiating (1) and (2) gives

$\stackrel{.}{X}={\stackrel{.}{x}}_0-\left(\mathrm{x}sin\theta +\eta cos\theta \right)\stackrel{.}{\theta }-sin\theta \stackrel{.}{\eta }$

(8)

$\stackrel{.}{Y}={\stackrel{.}{y}}_0-\left(\mathrm{x}cos\theta +\eta sin\theta \right)\stackrel{.}{\theta }+cos\theta \stackrel{.}{\eta }$

(9)

Squaring (8) and (9) and substituting into (7) yields

$\begin{array}{l}\mathrm{T}=\frac{1}{2}{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\rho }[{\stackrel{.}{x}}_0^2+{\stackrel{.}{y}}_0^2+{\stackrel{.}{\theta }}^2{\eta }^2+{\left(\mathrm{x}\stackrel{.}{\theta }+\stackrel{.}{\eta }\right)}^2\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}-2\stackrel{.}{\theta }\eta \left({\stackrel{.}{x}}_{0}cos\theta +{\stackrel{.}{y}}_0\sin \theta \right)\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+2\left(\stackrel{.}{\theta }\mathrm{x}+\stackrel{.}{\eta }\right)\left({\stackrel{.}{y}}_{0}cos\theta -{\stackrel{.}{x}}_{0}sin\theta \right)]dx\end{array}$

(10)

Neglecting the shear deformation and considering the bending of the beam, the potential energy of the beam due to elasticity can be obtained as

${\mathrm{U}}_{\mathrm{e}}=\frac{1}{2}{\displaystyle {\int }_{\mathrm{V}}{\sigma }_{xx}{\epsilon }_{xx}dV}$

(11)

where,

${\sigma }_{xx}=\mathrm{E}{\epsilon }_{xx},\hspace{0.17em}{\epsilon }_{xx}=-{\mathrm{y}}_{\mathrm{d}}{\eta }^{″},\hspace{0.17em}dV=dxdydz$

(12)

Further, σ_xx,ε_xx, dV and E denote the stress, strain component, infinitesimal volume and Young's Modulus of a beam element, respectively. Bending strain is measured at a distance y_d from the neutral axis of the beam.

Substituting (12) into (11) gives

$\begin{array}{l}{\mathrm{U}}_{\mathrm{e}}=\frac{1}{2}{\displaystyle {\int }_{\mathrm{V}}\mathrm{E}{\epsilon }_{xx}^{2}dV}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}=\frac{1}{2}{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\left[{{\eta }^{″}}^2\mathrm{E}{\displaystyle {\int }_0^{\mathrm{A}}{\mathrm{y}}_{\mathrm{d}}^{2}dA}\right]}dx\\ {\mathrm{U}}_{\mathrm{e}}=\frac{1}{2}{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\left[EI{{\eta }^{″}}^2\right]dx}\end{array}$

(13)

where, $\mathrm{I}={\displaystyle {\int }_0^{\mathrm{A}}{\mathrm{y}}_{\mathrm{d}}^{2}dA}$

Potential energy due to the gravitational force can be obtained as

${\mathrm{U}}_{\mathrm{g}}=\rho g{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\left({\mathrm{y}}_0+\mathrm{x}sin\theta \right)}dx=mg{y}_0$

(14)

Total potential energy can be calculated by the following

$\mathrm{U}={\mathrm{U}}_g+{\mathrm{U}}_{\mathrm{e}}$

(15)

Work done due to the external forces (Fig. 2) are formulated as

$\begin{array}{l}\mathrm{W}={\mathrm{F}}_{1\mathrm{x}}\left({\mathrm{x}}_0-\frac{\mathrm{L}}{2}\cos \theta \right)+{\mathrm{F}}_{1\mathrm{y}}\left({\mathrm{y}}_0-\frac{\mathrm{L}}{2}\sin \theta \right)\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+{\mathrm{F}}_{2\mathrm{x}}\left({\mathrm{x}}_0+\frac{\mathrm{L}}{2}\cos \theta \right)+{\mathrm{F}}_{2\mathrm{y}}\left({\mathrm{y}}_0+\frac{\mathrm{L}}{2}\sin \theta \right)\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left({\mathrm{F}}_{2\mathrm{y}}+{\mathrm{F}}_{1\mathrm{y}}\right)\eta \cos \theta -\left({\mathrm{F}}_{1\mathrm{x}}+{\mathrm{F}}_{2\mathrm{x}}\right)\eta \sin \theta \end{array}$

(16)

By applying the extended Hamilton principle,

${\displaystyle {\int }_{{\mathrm{t}}_1}^{{\mathrm{t}}_2}\left(\delta \mathrm{U}-\delta \mathrm{T}-\delta \mathrm{W}\right)}dt=0$

(17)

the following equations of rigid motion of the beam along X, Y and Z directions are obtained.

The equation of motion for translation along the X - direction is derived as

$\begin{array}{l}\mathrm{m}{\stackrel{..}{x}}_0-\left(\rho \cos \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\eta dx}\right)\stackrel{..}{\theta }-\rho \sin \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\stackrel{..}{\eta }dx}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}-2\rho \stackrel{.}{\theta }\cos \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\stackrel{.}{\eta }dx}+{\stackrel{.}{\theta }}^2\rho \sin \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\eta }dx\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}={\mathrm{F}}_{1\mathrm{x}}+{\mathrm{F}}_{2\mathrm{x}}\end{array}$

(18)

The translation in the Y-direction is described by the equation

$\begin{array}{l}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathrm{m}{\stackrel{..}{y}}_0-\left(\rho \sin \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\eta dx}\right)\stackrel{..}{\theta }-\rho \cos \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\stackrel{..}{\eta }dx}\\ -2\rho \stackrel{.}{\theta }\sin \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\stackrel{.}{\eta }dx}+{\stackrel{.}{\theta }}^2\rho \cos \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\eta }dx+mg\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}={\mathrm{F}}_{1\mathrm{y}}+{\mathrm{F}}_{2\mathrm{y}}\end{array}$

(19)

Rotation about the Z axis is described by

$\begin{array}{l}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\left(-\rho \cos \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\eta dx}\right){\stackrel{..}{x}}_0-\left(\rho \sin \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\eta dx}\right){\stackrel{..}{y}}_0\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\left(\frac{\mathrm{m}{\mathrm{L}}^2}{12}+\rho {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}{\eta }^{2}dx}\right)\stackrel{..}{\theta }+\rho {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\mathrm{x}\stackrel{..}{\eta }dx}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+2\rho \stackrel{.}{\theta }{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\eta \stackrel{.}{\eta }dx}={\mathrm{F}}_{1\mathrm{x}}\left(\frac{\mathrm{L}}{2}\sin \theta -\eta \cos \theta \right)\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+{\mathrm{F}}_{1\mathrm{y}}\left(-\frac{\mathrm{L}}{2}\cos \theta -\eta \sin \theta \right)\\ +{\mathrm{F}}_{2\mathrm{x}}\left(-\frac{\mathrm{L}}{2}\sin \theta -\eta \cos \theta \right)+{\mathrm{F}}_{2\mathrm{y}}\left(\frac{\mathrm{L}}{2}\cos \theta -\eta \sin \theta \right)\end{array}$

(20)

The differential equation of motion of transverse vibration of the beam is expressed as

$\begin{array}{l}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}-\sin \theta {\stackrel{..}{x}}_0+\cos \theta {\stackrel{..}{y}}_0+\mathrm{x}\stackrel{..}{\theta }+\stackrel{..}{\eta }-\eta {\stackrel{.}{\theta }}^2+\frac{EI}{\rho }{\eta }^{iv}\\ =\frac{1}{\mathrm{m}}\left[-{\mathrm{F}}_{1\mathrm{x}}\sin \theta +{\mathrm{F}}_{1\mathrm{y}}\cos \theta -{\mathrm{F}}_{2\mathrm{x}}\sin \theta +{\mathrm{F}}_{1\mathrm{y}}\cos \theta \right]\end{array}$

(21)#

Hence, from (18), (19) and (20) the beam dynamics of the rigid motion with respect to Cartesian coordinates {x₀ y₀ θ}^T can be written in a compact form as

${\mathrm{M}}_{rf}{\stackrel{..}{X}}_{rf}+{\mathrm{C}}_{rf}+{\eta }_{rf}+{\mathrm{G}}_{rf}={\mathrm{F}}_{rf}\left(-\mathrm{f}\right)$

(22)

where

$\begin{array}{c}{\mathrm{M}}_{rf}=\left[\begin{array}{ccc}\mathrm{m}& 0& {\mathrm{M}}_{rf1}\\ 0& \mathrm{m}& {\mathrm{M}}_{rf2}\\ {\mathrm{M}}_{rf1}& {\mathrm{M}}_{rf2}& {\mathrm{M}}_{rf3}\end{array}\right]\\ {\mathrm{C}}_{rf}=\left\{\begin{array}{c}{\stackrel{.}{\theta }}^2\left(\rho \sin \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\eta dx}\right)\\ -{\stackrel{.}{\theta }}^2\left(\rho \cos \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\eta dx}\right)\\ 0\end{array}\right\}\\ {\eta }_{rf}=\left\{\begin{array}{c}-\rho \sin \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\stackrel{..}{\eta }dx}-2\rho \stackrel{.}{\theta }\cos \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\stackrel{.}{\eta }dx}\\ \rho \cos \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\stackrel{..}{\eta }dx}-2\rho \stackrel{.}{\theta }\sin \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\stackrel{.}{\eta }dx}\\ \rho {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\mathrm{x}\stackrel{..}{\eta }dx}+2\rho \stackrel{.}{\theta }{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\eta \stackrel{.}{\eta }dx}\end{array}\right\}\\ {\mathrm{G}}_{rf}=\left\{\begin{array}{c}0\\ mg\\ 0\end{array}\right\};\hspace{0.17em}\hspace{0.17em}{\mathrm{F}}_{rf}=\left[\begin{array}{cccccc}1& 0& 0& 1& 0& 0\\ 0& 1& 0& 0& 1& 0\\ {\mathrm{F}}_1& {\mathrm{F}}_2& 0& {\mathrm{F}}_3& {\mathrm{F}}_4& 0\end{array}\right]\\ {\stackrel{..}{X}}_{rf}=\left\{\begin{array}{c}{\stackrel{..}{x}}_0\\ {\stackrel{..}{y}}_0\\ \stackrel{..}{\theta }\end{array}\right\};\hspace{0.17em}\mathrm{f}=\left\{\begin{array}{c}{\mathrm{F}}_{1\mathrm{x}}\\ {\mathrm{F}}_{1\mathrm{y}}\\ {\mathrm{M}}_{{\mathrm{O}}_1}\\ {\mathrm{F}}_{2\mathrm{x}}\\ {\mathrm{F}}_{2\mathrm{y}}\\ {\mathrm{M}}_{{\mathrm{O}}_2}\end{array}\right\}\end{array}$

where,

$\begin{array}{l}{\mathrm{M}}_{rf1}=-\rho \cos \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\eta dx}\\ {\mathrm{M}}_{rf2}=-\rho \sin \theta {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\eta dx}\\ {\mathrm{M}}_{rf3}=\frac{\mathrm{m}{\mathrm{L}}^2}{12}+\rho {\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}{\eta }^{2}dx}\end{array}$

and also the force terms are,

$\begin{array}{l}{\mathrm{F}}_1=\frac{\mathrm{L}}{2}\sin \theta -\eta \cos \theta \\ {\mathrm{F}}_2=-\frac{\mathrm{L}}{2}\cos \theta -\eta \sin \theta \\ {\mathrm{F}}_3=-\frac{\mathrm{L}}{2}\sin \theta -\eta \cos \theta \\ {\mathrm{F}}_4=\frac{\mathrm{L}}{2}\cos \theta -\eta \sin \theta \end{array}$

The flexible motion is the transverse vibration of beam dynamics (21) which can be rewritten as

$-\sin \theta {\stackrel{..}{x}}_0+\cos \theta {\stackrel{..}{y}}_0+\mathrm{x}\stackrel{..}{\theta }+\stackrel{..}{\eta }-\eta {\stackrel{.}{\theta }}^2+\frac{EI}{\rho }{\eta }^{iv}={\mathrm{F}}_{ff}\left(\mathrm{f}\right)$

(23)

The dynamic equations of the beam have been expressed in (22) and (23). In the following development, the dynamics of two manipulators rigidly grasping the beam will be combined with the beam dynamics (22) and (23) to form complete system dynamics of the manipulators-beam system.

3. Manipulator dynamics and combined dynamics

A. Manipulator Dynamics

The general manipulator dynamic equation can be expressed as [23]

${\mathrm{M}}_{\mathrm{i}}\left({\mathrm{q}}_{\mathrm{i}}\right){\stackrel{..}{q}}_{\mathrm{i}}+{\mathrm{C}}_{\mathrm{i}}\left({\mathrm{q}}_{\mathrm{i}},{\stackrel{.}{q}}_{\mathrm{i}}\right){\stackrel{.}{q}}_{\mathrm{i}}+{\mathrm{G}}_{\mathrm{i}}\left({\mathrm{q}}_{\mathrm{i}}\right)={\tau }_{\mathrm{i}}+{\mathrm{J}}_{\mathrm{i}}^{\mathrm{T}}{\mathrm{f}}_{\mathrm{i}}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}wherei=1,2$

(24)

M_i(q_i) represents symmetric positive definite inertia matrix, C_i(q_i,q̇_i)q̇_i is the vector due to Coriolis and centrifugal components, G_i (q_i) represents the vector of gravitational components, τ_i is the vector of input torque applied at each joints of the manipulator, f_i is the interaction force between the manipulator and the flexible beam, J_i is Jacobian matrix of manipulator, and q_i is the vector of joint angles. The two manipulator equations are assembled and can be written in joint space as

${\mathrm{M}}_{\mathrm{r}}\stackrel{..}{q}+{\mathrm{C}}_{\mathrm{r}}\stackrel{.}{q}+{\mathrm{G}}_{\mathrm{r}}=\tau +{\mathrm{J}}^{\mathrm{T}}\mathrm{f}$

(25)

where,

$\begin{array}{l}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{\mathrm{M}}_{\mathrm{r}}=\left[\begin{array}{cc}{\mathrm{M}}_1& 0\\ 0& {\mathrm{M}}_2\end{array}\right];\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{\mathrm{C}}_{\mathrm{r}}=\left[\begin{array}{cc}{\mathrm{C}}_1& 0\\ 0& {\mathrm{C}}_2\end{array}\right];\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{\mathrm{G}}_{\mathrm{r}}=\left\{\begin{array}{c}{\mathrm{G}}_1\\ {\mathrm{G}}_2\end{array}\right\}\\ \tau =\left\{\begin{array}{c}{\tau }_1\\ {\tau }_2\end{array}\right\};\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathrm{J}=\left[\begin{array}{cc}{\mathrm{J}}_1& 0\\ 0& {\mathrm{J}}_2\end{array}\right];\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathrm{f}=\left\{\begin{array}{c}{\mathrm{f}}_1\\ {\mathrm{f}}_2\end{array}\right\};\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathrm{q}=\left\{\begin{array}{c}{\mathrm{q}}_1\\ {\mathrm{q}}_2\end{array}\right\}\end{array}$

It can be seen from (22) that the beam dynamics is represented in Cartesian space and the manipulator dynamics (25) is given in joint space. In order to formulate a complete system of dynamic equations in Cartesian space, the manipulator dynamics will be converted into Cartesian space and the result will be combined with the beam dynamics to form the combined dynamics.

B. Combined Dynamics

The following procedure illustrates the procedure to convert the manipulator dynamics into the Cartesian space and finally to obtain combined dynamics. A well known kinematic relation between the end-effector velocity and joint velocity gives [24]

$\left\{{\stackrel{.}{e}}_{\mathrm{i}}\right\}=\left[{\mathrm{J}}_{\mathrm{i}}\right][\left\{{\stackrel{.}{q}}_{\mathrm{i}}\right\}$

(26)

For the two manipulators

$\left\{{\stackrel{.}{e}}_1\right\}=\left[{\mathrm{J}}_1\right]\left\{{\stackrel{.}{q}}_1\right\},\hspace{0.17em}\left\{{\stackrel{.}{e}}_2\right\}=\left[{\mathrm{J}}_2\right]\left\{{\stackrel{.}{q}}_2\right\}$

(27)

Equation (27) in an assembled form can be written as

$\left\{\stackrel{.}{e}\right\}=\left[\mathrm{J}\right]\left\{\stackrel{.}{q}\right\}$

(28)

where

$\stackrel{.}{e}=\left\{\begin{array}{c}{\stackrel{.}{e}}_1\\ {\stackrel{.}{e}}_2\end{array}\right\};\hspace{0.17em}\hspace{0.17em}\mathrm{J}=\left[\begin{array}{cc}{\mathrm{J}}_1& 0\\ 0& {\mathrm{J}}_2\end{array}\right]$

Substituting (5) into (28) yields

$\left\{\stackrel{.}{q}\right\}=\left[{\mathrm{J}}^{-1}\right]\left[\mathrm{R}\right]\left\{{\stackrel{.}{X}}_{rf}\right\}$

(29)

Differentiating (29) gives

$\stackrel{..}{q}={\mathrm{J}}^{-1}\stackrel{.}{R}{\stackrel{.}{X}}_{rf}+{\mathrm{J}}^{-1}\mathrm{R}{\stackrel{..}{X}}_{rf}+{\stackrel{.}{J}}^{-1}\mathrm{R}{\stackrel{.}{X}}_{rf}$

(30)

For simplicity, parentheses for vectors and square brackets for matrices are omitted in the following.

Substituting (29) and (30) into (25), we obtain the manipulator dynamics in Cartesian space

${\mathrm{M}}_{\mathrm{r}}{\mathrm{J}}^{-1}\mathrm{R}{\stackrel{..}{X}}_{rf}+\left({\mathrm{M}}_{\mathrm{r}}{\stackrel{.}{J}}^{-1}\mathrm{R}+{\mathrm{M}}_{\mathrm{r}}{\mathrm{J}}^{-1}\stackrel{.}{R}+{\mathrm{C}}_{\mathrm{r}}{\mathrm{J}}^{-1}\mathrm{R}\right){\stackrel{.}{X}}_{rf}+{\mathrm{G}}_{\mathrm{r}}\hspace{0.17em}=\tau +{\mathrm{J}}^{\mathrm{T}}\mathrm{f}$

(31)

Premultiplying (31) by R^TJ^–T gives

$\begin{array}{r}{\mathrm{R}}^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}{\mathrm{M}}_{\mathrm{r}}{\mathrm{J}}^{-1}\mathrm{R}{\stackrel{..}{X}}_{rf}\\ +{\mathrm{R}}^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}\left({\mathrm{M}}_{\mathrm{r}}{\stackrel{.}{J}}^{-1}\mathrm{R}+{\mathrm{M}}_{\mathrm{r}}{\mathrm{J}}^{-1}\stackrel{.}{R}+{\mathrm{C}}_{\mathrm{r}}{\mathrm{J}}^{-1}\mathrm{R}\right){\stackrel{.}{X}}_{rf}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+{\mathrm{R}}^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}{\mathrm{G}}_{\mathrm{r}}={\mathrm{R}}^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}\tau +{\mathrm{R}}^{\mathrm{T}}\mathrm{f}\end{array}$

(32)

In view of the assumption of simply supported beam boundary conditions, the moments at the two ends are zero. However, in reality manipulators experience forces as well as moments at the two ends of the beam [2]. Utilizing this fact, the moments at the two ends of the beam should be included in the matrix F_rf of (22), which results in

${\mathrm{F}}_{rf}=\left[\begin{array}{cccccc}1& 0& 0& 1& 0& 0\\ 0& 1& 0& 0& 1& 0\\ {\mathrm{F}}_1& {\mathrm{F}}_2& 1& {\mathrm{F}}_3& {\mathrm{F}}_4& 1\end{array}\right]$

In this case, F_rf is equivalent to R^T. Combining the rigid motion of the beam dynamics (22) with the manipulator dynamics in Cartesian space (32) yields the combined rigid dynamics as

$\begin{array}{r}{\mathrm{R}}^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}{\mathrm{M}}_{\mathrm{r}}{\mathrm{J}}^{-1}\mathrm{R}{\stackrel{..}{X}}_{rf}\\ +{\mathrm{R}}^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}\left({\mathrm{M}}_{\mathrm{r}}{\stackrel{.}{J}}^{-1}\mathrm{R}+{\mathrm{M}}_{\mathrm{r}}{\mathrm{J}}^{-1}\stackrel{.}{R}+{\mathrm{C}}_{\mathrm{r}}{\mathrm{J}}^{-1}\mathrm{R}\right){\stackrel{.}{X}}_{rf}+{\mathrm{R}}^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}{\mathrm{G}}_{\mathrm{r}}\\ ={\mathrm{R}}^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}\tau -\left({\mathrm{M}}_{rf}{\stackrel{..}{X}}_{rf}+{\mathrm{C}}_{rf}+{\eta }_{rf}+{\mathrm{G}}_{rf}\right)\end{array}$

(33)

The above combined rigid dynamics can be rewritten as

${\mathrm{M}}_{orf}{\stackrel{..}{X}}_{rf}+{\mathrm{C}}_{orf}{\stackrel{.}{X}}_{rf}+{\mathrm{G}}_{orf}+{\eta }_{orf}={\mathrm{u}}_{orf}$

(34)

where,

$\begin{array}{l}{\mathrm{M}}_{orf}={\mathrm{R}}^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}{\mathrm{M}}_{\mathrm{r}}{\mathrm{J}}^{-1}\mathrm{R}+{\mathrm{M}}_{rf}\\ {\mathrm{C}}_{orf}={\mathrm{R}}^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}\left({\mathrm{M}}_{\mathrm{r}}{\stackrel{.}{J}}^{-1}\mathrm{R}+{\mathrm{M}}_{\mathrm{r}}{\mathrm{J}}^{-1}\stackrel{.}{R}+{\mathrm{C}}_{\mathrm{r}}{\mathrm{J}}^{-1}\mathrm{R}\right)+{\mathrm{C}}_{rf}\\ {\mathrm{G}}_{orf}={\mathrm{R}}^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}{\mathrm{G}}_{\mathrm{r}}+{\mathrm{G}}_{rf}\\ {\eta }_{orf}={\eta }_{rf}\\ {\mathrm{u}}_{orf}={\mathrm{R}}^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}\tau \end{array}$

Taking into account the transverse vibration of beam dynamics, the complete manipulator-beam system dynamics is represented as

${\mathrm{M}}_{orf}{\stackrel{..}{X}}_{rf}+{\mathrm{C}}_{orf}{\stackrel{.}{X}}_{rf}+{\mathrm{G}}_{orf}+{\eta }_{orf}={\mathrm{u}}_{orf}$

(35)

$-\sin \theta {\stackrel{..}{x}}_0+\cos \theta {\stackrel{..}{y}}_0+\mathrm{x}\stackrel{..}{\theta }+\stackrel{..}{\eta }-\eta {\stackrel{.}{\theta }}^2+\frac{EI}{\rho }{\eta }^{iv}={\mathrm{F}}_{ff}\left(\mathrm{f}\right)$

(36)

4. Singular perturbation model

The control task is stated as follows: for any given desired bounded trajectory of the mass centre of the beam X_rfd and X_rfd, with some or all of the manipulator and beam parameters unknown, derive a controller for the manipulators τ such that the beam centre X_rf tracks X_rfd, while suppressing the vibration of the flexible object, to zero.

The system dynamics derived without using any approximation methods given in (35) and (36) have rigid as well as flexible parameters that are coupled together. To achieve the above control objective for such a complicated nonlinear system, a suitable control approach is the two-time scale theory which considers the high frequency phenomenon of flexible motion in different time scales. The basic idea for the two-time scale theory is to identify the slow and fast subsystems in separate time scales by employing the singular perturbation approach [25]. Then, control algorithms are designed for each of the subsystems, which will be combined to yield the composite control strategy for the original system. However, the challenge is to design the subcontrollers to satisfy the so-called Tikhnov theorem in order to guarantee that the composite controller can be applied to the original system, especially when the parameters of the system are unknown.

Since the manipulators are assumed to be rigid, the corresponding inertia matrix M_r, Coriolis and centrifugal matrix C_r, gravitational vectors G_r and Jacobian matrix J do not have any flexible parameters. The beam has rigid as well as flexible parameters and are coupled together in M_rf, C_rf, η_rf and also in R. These flexible parameters have to be uncoupled from the above matrices and vectors by using the singular perturbation technique. This technique also accounts for the neglected high frequency characteristics when the beam undergoes vibration [20]. Using a perturbation parameter, say ε², the order of the system dynamics can be changed and this small parameter depends upon the system variable. Keeping that in mind, the term EI / ρ in (36), which has large magnitude compared to other coefficients [20], can be redefined as,

$\frac{EI}{\rho }\hspace{0.17em}\triangleq \mathrm{a}.\mathrm{K}$

(37)

where K is a dimensionless parameter which has a large value for the different materials [20] and [25], and its order is equal to EI / ρ and also the variable “a” satisfies the equalities. For example, an aluminium rod with a diameter of 0.05m, E=71 GPa and ρ = 2700kg / m³ has the value of the co-efficient EI / ρ ≅ 4.1 × 10³ and therefore a = 4.1 and K = 103.

The beam has rigid motion with respect to the state variables X_rf = {x₀ y₀ θ}^T and also the transverse vibration η with respect to the state variable occurs in different time scales. Then, one needs to introduce a new variable w(x,t) in the same order as the state variable by using the following

$\eta \left(\mathrm{x},\mathrm{t}\right)\hspace{0.17em}\triangleq \hspace{0.17em}{\epsilon }^2.\mathrm{w}\left(\mathrm{x},\mathrm{t}\right)$

(38)

where ε² = 1 / K is the so-called perturbed parameter.

Using (38) one can rewrite the rigid motion dynamics of the beam as follows:

The equation of motion (18) in the X-direction can be written in terms of perturbed parameter as

$\begin{array}{l}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathrm{m}{\stackrel{..}{x}}_0-\left(\rho \cos \theta {\epsilon }^2{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}wdx}\right)\stackrel{..}{\theta }-\rho \sin \theta {\epsilon }^2{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\stackrel{..}{w}dx}\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}-2\rho \stackrel{.}{\theta }\cos \theta {\epsilon }^2{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\stackrel{.}{w}dx}+{\stackrel{.}{\theta }}^2\rho \sin \theta {\epsilon }^2{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\mathrm{w}}dx\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}={\mathrm{F}}_{2\mathrm{x}}+{\mathrm{F}}_{1\mathrm{x}}\end{array}$

(39)

The equation of motion (19) in the Y-direction can be written in terms of perturbed parameters as

$\begin{array}{l}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\mathrm{m}{\stackrel{..}{y}}_0-\left(\rho \sin \theta {\epsilon }^2{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}wdx}\right)\stackrel{..}{\theta }+\rho \cos \theta {\epsilon }^2{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\stackrel{..}{w}dx}\\ -2\rho \stackrel{.}{\theta }\sin \theta {\epsilon }^2{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\stackrel{.}{w}dx}+{\stackrel{.}{\theta }}^2\rho \cos \theta {\epsilon }^2{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}\mathrm{w}}dx+mg\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}={\mathrm{F}}_{2\mathrm{y}}+{\mathrm{F}}_{1\mathrm{y}}\end{array}$

(40)

Rotation about the Z axis described in the equation of motion (20) can be rewritten in terms of perturbed parameter as

$\begin{array}{r}\left(-\rho \cos \theta {\epsilon }^2{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}wdx}\right){\stackrel{..}{x}}_0-\left(-\rho \sin \theta {\epsilon }^2{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}wdx}\right){\stackrel{..}{y}}_0\\ \hspace{0.17em}\left(\frac{\mathrm{m}{\mathrm{L}}^2}{12}+\rho {\epsilon }^4{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}{\mathrm{w}}^{2}dx}\right)\stackrel{..}{\theta }+\rho {\epsilon }^2{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}x\stackrel{..}{w}dx}\\ +2\rho \stackrel{.}{\theta }{\epsilon }^4{\displaystyle {\int }_{-\frac{\mathrm{L}}{2}}^{\frac{\mathrm{L}}{2}}w\stackrel{.}{w}dx\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}}\\ ={\mathrm{F}}_{1\mathrm{x}}\left(\frac{\mathrm{L}}{2}\sin \theta -{\epsilon }^2\mathrm{w}\cos \theta \right)+{\mathrm{F}}_{1\mathrm{y}}\left(-\frac{\mathrm{L}}{2}\cos \theta -{\epsilon }^2\mathrm{w}\sin \theta \right)\\ +{\mathrm{F}}_{2\mathrm{x}}\left(-\frac{\mathrm{L}}{2}\sin \theta -{\epsilon }^2\mathrm{w}\cos \theta \right)+{\mathrm{F}}_{2\mathrm{y}}\left(\frac{\mathrm{L}}{2}\cos \theta -{\epsilon }^2\mathrm{w}\sin \theta \right)\end{array}$

(41)

Accordingly, utilizing (37) and (38), the equation of motion for transverse vibration of the beam can be rewritten as

$-\sin \theta {\stackrel{..}{x}}_0+\cos \theta {\stackrel{..}{y}}_0+\mathrm{x}\stackrel{..}{\theta }+{\epsilon }^2\stackrel{..}{w}-{\epsilon }^2\mathrm{w}{\stackrel{.}{\theta }}^2+\mathrm{a}{\mathrm{w}}^{iv}={\mathrm{F}}_{ff}\left(\mathrm{f}\right)$

(42)

The equations (39), (40), (41) and (42) represent the singularly perturbed form which will be incorporated into the system of dynamic equations (35) and (36) to form the singularly perturbed model of the complete system.

5. Slow and fast subsystems

In this section, the two subsystems, namely slow and fast, are obtained by following the typical steps of the singular perturbation approach. Due to the presence of perturbed parameter ε², the complete system has two motions in the different time scales. Initially, the singularly perturbed model of rigid motion dynamics of the beam will be converted into the rigid dynamic model without involving any flexible parameter, and finally, it will be incorporated into (35) to form the slow subsystem. Similarly, after following the usual procedure of the singular perturbation approach, the fast subsystem will be obtained from (42).

A. Slow Subsystem

When the perturbation parameter ε approaches zero, the equivalent quasi-steady-state system [26] represents the slow subsystem. By setting ε = 0, (39), (40) and (41) we form the rigid dynamic model of the beam which is given in compact form as

${\mathrm{M}}_{rd}{\stackrel{..}{X}}_{rf}+{\mathrm{C}}_{rd}+{\eta }_{rd}+{\mathrm{G}}_{rd}={\mathrm{F}}_{rd}\left(-\mathrm{f}\right)$

(43)

$\begin{array}{c}{\mathrm{M}}_{rd}=\left[\begin{array}{ccc}\mathrm{m}& 0& 0\\ 0& \mathrm{m}& 0\\ 0& 0& \frac{\mathrm{m}{\mathrm{L}}^2}{12}\end{array}\right];\hspace{0.17em}{\stackrel{..}{X}}_{rf}=\left\{\begin{array}{c}{\stackrel{..}{x}}_0\\ {\stackrel{..}{y}}_0\\ \stackrel{..}{\theta }\end{array}\right\};\hspace{0.17em}\hspace{0.17em}{\mathrm{C}}_{rd}={\left\{\begin{array}{ccc}0& 0& 0\end{array}\right\}}^{\mathrm{T}};\\ {\eta }_{rd}={\left\{\begin{array}{ccc}0& 0& 0\end{array}\right\}}^{\mathrm{T}};\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{\mathrm{G}}_{rd}=\left\{\begin{array}{c}0\\ mg\\ 0\end{array}\right\}\\ {\mathrm{F}}_{rd}=\left[\begin{array}{cccccc}1& 0& 0& 1& 0& 0\\ 0& 1& 0& 0& 1& 0\\ \frac{\mathrm{L}}{2}\sin \theta & -\frac{\mathrm{L}}{2}\cos \theta & 1& -\frac{\mathrm{L}}{2}\sin \theta & \frac{\mathrm{L}}{2}\cos \theta & 1\end{array}\right]\end{array}$

The equation of motion for transverse vibration of the beam (42) becomes

${\left[-\sin \theta {\stackrel{..}{x}}_0+\cos \theta {\stackrel{..}{y}}_0+\mathrm{x}\stackrel{..}{\theta }+\mathrm{a}{\mathrm{w}}^{iv}\right]}_{|\mathrm{s}}={\mathrm{F}}_{ff}\left({\mathrm{f}}_{\mathrm{s}}\right)$

(44)

where f_s corresponds to f when ε = 0.

Again, substituting (38) into (5) and also setting ε = 0,R becomes R₁ which is given by

${\mathrm{R}}_1=\left[\begin{array}{ccc}1& 0& \frac{\mathrm{L}}{2}\sin \theta \\ 0& 1& -\frac{\mathrm{L}}{2}\cos \theta \\ 0& 0& 1\\ 1& 0& -\frac{\mathrm{L}}{2}\sin \theta \\ 0& 1& \frac{\mathrm{L}}{2}\cos \theta \\ 0& 0& 1\end{array}\right]$

Based upon the above results, the combined dynamic Equation (35) becomes

$\left({\mathrm{M}}_{\mathrm{o}}+{\mathrm{M}}_{rd}\right){\stackrel{..}{X}}_{rf}+{\mathrm{C}}_{\mathrm{o}}{\stackrel{.}{X}}_{rf}+{\mathrm{G}}_{\mathrm{o}}+{\mathrm{G}}_{rd}={\mathrm{u}}_{\mathrm{o}}$

(45)

where

$\begin{array}{l}{\mathrm{M}}_{\mathrm{o}}={\mathrm{R}}_1^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}{\mathrm{M}}_{\mathrm{r}}{\mathrm{J}}^{-1}{\mathrm{R}}_1\\ {\mathrm{C}}_{\mathrm{o}}={\mathrm{R}}_1^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}\left({\mathrm{M}}_{\mathrm{r}}{\stackrel{.}{J}}^{-1}{\mathrm{R}}_1+{\mathrm{M}}_{\mathrm{r}}{\mathrm{J}}^{-1}{\stackrel{.}{R}}_1+{\mathrm{C}}_{\mathrm{r}}{\mathrm{J}}^{-1}{\mathrm{R}}_1\right)\\ {\mathrm{G}}_{\mathrm{o}}={\mathrm{R}}_1^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}{\mathrm{G}}_{\mathrm{r}}\\ {\mathrm{u}}_{\mathrm{o}}={\mathrm{R}}_1^{\mathrm{T}}{\mathrm{J}}^{-\mathrm{T}}\tau \end{array}$

The above equation can be rewritten as

${\mathrm{M}}_{or}{\stackrel{..}{X}}_{rf}+{\mathrm{C}}_{or}{\stackrel{.}{X}}_{rf}+{\mathrm{G}}_{or}={\mathrm{u}}_{or}$

(46)

where

$\begin{array}{l}{\mathrm{M}}_{or}={\mathrm{M}}_{\mathrm{o}}+{\mathrm{M}}_{rd}\\ {\mathrm{C}}_{or}={\mathrm{C}}_{\mathrm{o}}\\ {\mathrm{G}}_{or}={\mathrm{G}}_{\mathrm{o}}+{\mathrm{G}}_{rd}\\ {\mathrm{u}}_{or}={\mathrm{u}}_{\mathrm{o}}\end{array}$

The slow subsystem given in (46) represents the rigid body motion without involving any flexible parameters. It will be used further to design a control algorithm to track a desired trajectory of the object.

The slow subsystem has the following properties which are important to the stability analysis and these properties can be proved.

Property 1: M_or is a symmetric positive definite matrix.

Property 2: The matrix M_or and C_or in (46) has to satisfy

${\mathrm{X}}^{\mathrm{T}}\left({\stackrel{.}{M}}_{or}-2{\mathrm{C}}_{or}\right)\mathrm{X}=0,\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\forall \mathrm{X}\ne 0$

(47)

where X is any arbitrary vector. Hence, (Ṁ_or − 2C_or) is a skew-symmetric matrix.

Property 3: A vector α_or ∈ R^v×1 exists which solely depends on manipulators and beam dynamic parameters (link lengths, masses and moment of inertia etc.) such that

${\mathrm{M}}_{or}{\stackrel{..}{X}}_{rf}+{\mathrm{C}}_{or}{\stackrel{.}{X}}_{rf}+{\mathrm{G}}_{or}={\mathrm{Y}}_{or}\left({\stackrel{..}{X}}_{rf},{\stackrel{.}{X}}_{rf},{\mathrm{X}}_{rf},\stackrel{.}{q},\mathrm{q}\right){\alpha }_{or}$

(48)

where Y_or ∈ R^u×v is called a regressor matrix of the manipulators-beam system in Cartesian space.

B. Fast Subsystem

Equation (42) represents the perturbed flexible model obtained with perturbation parameter ε which is very small and solely depends upon the E,I and ρ.

In order to study the dynamic behaviour of the fast system, the so-called boundary layer phenomenon [25] and [26] must be obtained. This can be identified by ensuring that the slow variables are kept constant in the fast time scale $\mathrm{v}=\frac{\mathrm{t}-{\mathrm{t}}_0}{\epsilon }.$ . From the typical steps of singular perturbation [25] one can define the fast variable w_f

${\mathrm{w}}_{\mathrm{f}}=\mathrm{w}-{\mathrm{w}}_{\mathrm{s}}$

(49)

Differentiating the fast time scale

$\begin{array}{l}\mathrm{d}\nu =\frac{dt}{\epsilon }\\ \frac{\mathrm{d}\nu }{dt}=\frac{1}{\epsilon }\end{array}$

(50)

Differentiating (49) gives

$\stackrel{.}{w}={\stackrel{.}{w}}_{\mathrm{s}}+\frac{\mathrm{d}}{dt}{\mathrm{w}}_{\mathrm{f}}={\stackrel{.}{w}}_{\mathrm{s}}+\frac{\mathrm{d}\nu }{dt}\frac{\mathrm{d}}{\mathrm{d}\nu }{\mathrm{w}}_{\mathrm{f}}={\stackrel{.}{w}}_{\mathrm{s}}+\frac{1}{\epsilon }{\widehat{w}}_{\mathrm{f}}$

(51)

where ŵ_f denotes differentiating fast variable with respect to fast time scale.

Differentiating (51) again yields

$\stackrel{..}{w}={\stackrel{..}{w}}_{\mathrm{s}}+\frac{\mathrm{d}\nu }{dt}\frac{\mathrm{d}}{\mathrm{d}\nu }{\widehat{w}}_{\mathrm{f}}={\stackrel{..}{w}}_{\mathrm{s}}+\frac{1}{{\epsilon }^2}{\widehat{\widehat{w}}}_{\mathrm{f}}$

(52)

Utilizing (49), (42) can be rewritten as

$\begin{array}{l}-\sin \theta {\stackrel{..}{x}}_0+\cos \theta {\stackrel{..}{y}}_0+\mathrm{x}\stackrel{..}{\theta }+{\epsilon }^2\stackrel{..}{w}-{\epsilon }^2\left({\mathrm{w}}_{\mathrm{s}}+{\mathrm{w}}_{\mathrm{f}}\right){\stackrel{.}{\theta }}^2\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}+\mathrm{a}\left({\mathrm{w}}_{\mathrm{s}}^{iv}+{\mathrm{w}}_{\mathrm{f}}^{iv}\right)={\mathrm{F}}_{ff}\left(\mathrm{f}\right)\end{array}$

(53)

Utilizing (44) and (52), (53) becomes

${\mathrm{F}}_{ff}\left({\mathrm{f}}_{\mathrm{s}}\right)+\widehat{\widehat{w}}+{\epsilon }^2{\stackrel{..}{w}}_{\mathrm{s}}-{\epsilon }^2\left({\mathrm{w}}_{\mathrm{s}}+{\mathrm{w}}_{\mathrm{f}}\right){\stackrel{.}{\theta }}^2+\mathrm{a}{\mathrm{w}}_{\mathrm{f}}^{iv}={\mathrm{F}}_{ff}\left(\mathrm{f}\right)$

(54)

By defining F_{f f}(f_f) = F_{f f}(f) – F_{f f}(f_s), the above equation can be rewritten as

${\widehat{\widehat{w}}}_{\mathrm{f}}+{\epsilon }^2{\stackrel{..}{w}}_{\mathrm{s}}-{\epsilon }^2\left({\mathrm{w}}_{\mathrm{s}}+{\mathrm{w}}_{\mathrm{f}}\right){\stackrel{.}{\theta }}^2+\mathrm{a}{\mathrm{w}}_{\mathrm{f}}^{iv}={\mathrm{F}}_{ff}\left({\mathrm{f}}_{\mathrm{f}}\right)$

(55)

However, in the boundary layer system, the slow variable w_s is constant which implies ẅ = 0 and also ε = 0 [26]. Then, the fast dynamics can be represented as

${\widehat{\widehat{w}}}_{\mathrm{f}}+\mathrm{a}{\mathrm{w}}_{\mathrm{f}}^{iv}={\mathrm{F}}_{ff}\left({\mathrm{f}}_{\mathrm{f}}\right)$

(56)

Following [19] and [27], an operator A is defined as

$\begin{array}{l}\mathrm{D}\left(\mathrm{A}\right)=\{{\mathrm{w}}_{\mathrm{f}}|{\mathrm{w}}_{\mathrm{f}}^{iv}\in \mathrm{H},{\mathrm{w}}_{\mathrm{f}}\left(0\right)={w^{″}}_{\mathrm{f}}\left(0\right)\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}={\mathrm{w}}_{\mathrm{f}}\left(\mathrm{L}\right)={w^{″}}_{\mathrm{f}}\left(\mathrm{L}\right)=0\}\end{array}$

(57)

$\mathrm{A}{\mathrm{w}}_{\mathrm{f}}=\mathrm{a}{\mathrm{w}}_{\mathrm{f}}^{iv},\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\forall {\mathrm{w}}_{\mathrm{f}}\in \mathrm{D}\left(\mathrm{A}\right)$

(58)

where, D(A) denotes the domain of the operator A and H denotes the Hilbert space. The operator A has the eigenvalues λ_i and the corresponding eigenfunctions ϕ_i satisfying the following conditions [27].

0 < λ₁ < λ₂ ……, limλ_i = ∞

Aϕ_i = λ_iϕ_i, i = 1, 2….∞

The set of the eigenfunctions forms a complete orthonormal system in the Hilbert space.

Utilizing (58), the partial differential Equation (56) can be rewritten as an abstract differential equation on H as

$\begin{array}{l}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}{\widehat{\widehat{w}}}_{\mathrm{f}}(\nu )+\mathrm{A}{\mathrm{w}}_{\mathrm{f}}(\nu )={\mathrm{F}}_{ff}\left({\mathrm{f}}_{\mathrm{f}}\right)\\ {\mathrm{w}}_{\mathrm{f}}\left(0\right)={\mathrm{w}}_{{\mathrm{f}}_0},\hspace{0.17em}\hspace{0.17em}{\stackrel{.}{w}}_{\mathrm{f}}\left(0\right)={\mathrm{w}}_{{\mathrm{f}}_1}\end{array}$

(59)

Equation (59) represents the fast subsystem which will be used for designing fast feedback control.

6. Composite control

A composite control for the two manipulators and beam system is given by

$\mathrm{u}={\mathrm{u}}_{\mathrm{s}}\left({\stackrel{.}{X}}_{rf},{\mathrm{X}}_{rf},\mathrm{t}\right)+{\mathrm{u}}_{\mathrm{f}}\left({\widehat{w}}_{\mathrm{f}},{\mathrm{w}}_{\mathrm{f}}\right)$

where u_s is designed based on slow subsystem (46) and u_f is designed to stabilize the fast subsystem (59).

A. Robust Controller Design for Slow Subsystem

The key issue is that with the unknown manipulator and beam parameters, the design of u_s(Ẋ_rf,X_rf,t) cannot be arbitrary. It has to guarantee the exponential tracking of the desired trajectories so that Tikhnov's theorem can be satisfied, which will be clear in the later development. For that purpose, a sliding mode control approach will be adopted.

Define the tracking error as

${\mathrm{e}}_{\mathrm{r}}={\mathrm{X}}_{rf}-{\mathrm{X}}_{rfd}$

(60)

where X_rfd is the desired trajectory and the auxiliary trajectory can also be defined as

${\stackrel{.}{X}}_{\mathrm{r}}={\stackrel{.}{X}}_{rfd}-\lambda {\mathrm{e}}_{\mathrm{r}}$

(61)

where λ is a positive definite matrix whose eigenvalues are strictly in the right half of the complex plane.

The sliding surface can be chosen as

$\mathrm{S}={\stackrel{.}{X}}_{rf}-{\stackrel{.}{X}}_{\mathrm{r}}={\stackrel{.}{e}}_{\mathrm{r}}+\lambda {\mathrm{e}}_{\mathrm{r}}$

(62)

The sliding mode controller can be given as

${\mathrm{u}}_{\mathrm{s}}={\mathrm{u}}_{or}={\mathrm{Y}}_{or}\psi -{\mathrm{K}}_{\mathrm{D}}\mathrm{S}$

(63)

where K_D is a positive definite gain matrix, ${\mathrm{Y}}_{or}({\stackrel{..}{X}}_{\mathrm{r}},{\stackrel{.}{X}}_{\mathrm{r}},\stackrel{.}{q},\mathrm{q})$ is the regressor matrix defined in (48), α_or is the parameter vector and ψ = [ψ₁….ψ_m]^T are the switching functions which are given by

$\psi =-\Gamma \frac{{\mathrm{Y}}_{or}^{\mathrm{T}}\mathrm{S}}{\Vert {\mathrm{Y}}_{or}^{\mathrm{T}}\mathrm{S}\Vert }$

(64)

where Γ ≥ ‖α_or‖ fall is the upperbound of α_or which is known though it could be conservatively selected.

B. Control Design for Fast Subsystem

The objective of the controller is to suppress the vibration of the flexible object by incorporating the following feedback control law

${\mathrm{u}}_{\mathrm{f}}=\left({\mathrm{f}}_{\mathrm{f}}\right)=-\Pi {\mathrm{F}}_{ff}^{†}{\widehat{w}}_{\mathrm{f}}(\nu )$

(65)

where F_ff^† can be found using pseudo inverse. The operator ∏# is neither self-adjoint nor positive definite, and it is also shown in [19] and [28] that it is A-symmetric and A- positive semi-definite. This operator was formulated in [28] as Π = kQA, where Q is a bounded and positive definite operator. In addition, the velocity signal ŵ_f (v) can be measured using a velocity sensor. It is to be noted that there are some established results available using the velocity feedback, for example in [29] and [30]. However, they have considered the assumed modes where, as in the present control algorithm, no form of approximation is used. Additionally, the recent boundary control methods available in the literature [20] and [21], utilizes velocity and slope of fast state variables as feedback which depends upon the boundary conditions of the beam. However, the present control algorithm uses only velocity as feedback which is irrespective of boundary conditions. This controller is simple to implement in real time and needs only a limited number of sensors.

Substituting (65) into (59) gives a closed loop system defined by

$\widehat{\widehat{w}}(\mathrm{v})+\mathrm{\Pi }{\widehat{w}}_{\mathrm{f}}(\mathrm{v})+\mathrm{A}{\mathrm{w}}_{\mathrm{f}}(\mathrm{v})=0$

(66)

Utilizing the operators Π and Q, (66) can be rewritten as

$\begin{array}{c}\hfill {\widehat{\widehat{w}}}_{\mathrm{f}}(\mathrm{v})+kQA{\widehat{w}}_{\mathrm{f}}(\mathrm{v})+\mathrm{A}{\mathrm{w}}_{\mathrm{f}}(\mathrm{v})=0\hfill \\ \hfill {\mathrm{w}}_{\mathrm{f}}(0)={\mathrm{w}}_{{\mathrm{f}}_0},\hspace{0.17em}\hspace{0.17em}{\stackrel{.}{w}}_{\mathrm{f}}(0)={\mathrm{w}}_{{\mathrm{f}}_1}\hfill \end{array}$

(67)

where k is the positive gain and the term QAŵ_f(v)is the special damping term [28]. This damping has been studied by various researchers especially [31]-[34]. The two operators Q and A are related by Q = A ^β and β varies between $\left[\frac{-1}{2},0\right].$ . It is shown analytically by Huang [32] that when $\beta =\frac{-1}{2}$ , the damping term QAŵ_f(v) becomes ${\mathrm{A}}^{\frac{1}{2}}{\widehat{w}}_{\mathrm{f}}(\mathrm{v}).$ . This corresponds to structural damping which can be also seen in [34]. If β = 0, then the damping term exhibits strong damping or overdamping characteristics [32] and [33].

C. Stability Analysis

Tikhnov's theorem [35] is used to justify a composite control strategy that uses a fast feedback control to stabilize the fast subsystem together with a slow control based on a slow subsystem. Tikhnov's theorem basically states the following with regard to the complete manipulators-flexible beam system under the composite control.

If the slow subsystem in (46) is exponentially stable defined on an interval [t₀, ∞] and the boundary-layer system or fast subsystem in (59) tends to be also exponentially stable, uniformly in [t, slow variables], then there exists ε^* such that for all ε < ε^*

$\begin{array}{l}\mathrm{w}(\mathrm{t})={\mathrm{w}}_{\mathrm{s}}(\mathrm{t})+{\mathrm{w}}_{\mathrm{f}}(\mathrm{v})+\mathrm{O}(\epsilon )\\ {\mathrm{X}}_{rf}(\mathrm{t})={\mathrm{X}}_{rf}(\mathrm{t})+\mathrm{O}(\epsilon )\end{array}$

holds uniformly for t ∈ [t₀, ∞], where O(ε) ≃ w(t,ε) – w_s(t), represents the order of the error which arises due to the difference of the original variable w starting at t₀ from a prescribed initial conditions δ(ε) and the quasi-steady state variable w_s is not free to start from a prescribed initial value. There can be a discrepancy between its initial value and already prescribed initial value which causes this order of error. It will be hold in the interval excluding t₀ that is, for t ∈ [t_b, t₁] where t_b > t₀.

1) Exponential Stability of Slow Subsystem: Tikhnov's theorem requires the slow subsystem to be exponentially stable. The following derivation will illustrate the exponential stability of the slow subsystem for all t > 0. . Differentiating the sliding surface (62) with respect to time. Gives

$\stackrel{.}{S}={\stackrel{..}{X}}_{rf}-{\stackrel{..}{X}}_{\mathrm{r}}$

(68)

multiplying both sides of (68) by M_or and using (46), (68) can be rewritten as

${\mathrm{M}}_{or}\stackrel{.}{S}={\mathrm{u}}_{or}-{\mathrm{C}}_{or}{\stackrel{.}{X}}_{rf}-{\mathrm{G}}_{or}-{\mathrm{M}}_{or}{\stackrel{..}{X}}_{\mathrm{r}}$

(69)

Adding and subtracting C_orX̄_r in (69) results in

${\mathrm{M}}_{\text{or}}\stackrel{.}{\mathrm{S}}={\mathrm{u}}_{\text{or}}-({\mathrm{M}}_{\text{or}}{\stackrel{..}{\mathrm{X}}}_{\mathrm{r}}+{\mathrm{C}}_{\text{or}}{\stackrel{.}{\mathrm{X}}}_{\mathrm{r}}+{\mathrm{G}}_{\text{or}})+{\mathrm{C}}_{\text{or}}{\stackrel{.}{\mathrm{X}}}_{\mathrm{r}}-{\mathrm{C}}_{\text{or}}{\stackrel{.}{\mathrm{X}}}_{\text{rf}}$

(70)

Using (62), (70) can be rewritten as

$({\mathrm{M}}_{or}{\stackrel{..}{X}}_{\mathrm{r}}+{\mathrm{C}}_{or}{\stackrel{.}{X}}_{\mathrm{r}}+{\mathrm{G}}_{or})={\mathrm{Y}}_{or}({\stackrel{..}{X}}_{\mathrm{r}},{\stackrel{.}{X}}_{\mathrm{r}},\stackrel{.}{q},\mathrm{q}){\alpha }_{or}$

(71)

where

$({\mathrm{M}}_{or}{\stackrel{..}{X}}_{\mathrm{r}}+{\mathrm{C}}_{or}{\stackrel{.}{X}}_{\mathrm{r}}+{\mathrm{G}}_{or})={\mathrm{Y}}_{or}({\stackrel{..}{X}}_{\mathrm{r}},{\stackrel{.}{X}}_{\mathrm{r}},\stackrel{.}{q},\mathrm{q}){\alpha }_{or}$

(72)

Consider a Lyapunov function candidate as

${\mathrm{V}}_1(\mathrm{t},\mathrm{S})=\frac{1}{2}{\mathrm{S}}^{\mathrm{T}}{\mathrm{M}}_{or}\mathrm{S}$

(73)

Differentiating (73) with respect to time gives

${\stackrel{.}{V}}_1(\mathrm{t},\mathrm{S})={\mathrm{S}}^{\mathrm{T}}{\mathrm{M}}_{or}\stackrel{.}{S}+\frac{1}{2}{\mathrm{S}}^{\mathrm{T}}{\stackrel{.}{M}}_{or}\mathrm{S}$

(74)

Substituting (71) into (74) and also using property 2 given in (47), the above equation yields

${\stackrel{.}{V}}_1(\mathrm{t},\mathrm{S})={\mathrm{S}}^{\mathrm{T}}[{\mathrm{u}}_{or}-{\mathrm{Y}}_{or}{\alpha }_{or}]$

(75)

Substituting the control law given in (63) and (64) into (75) results in

$\le -{\mathrm{S}}^{\mathrm{T}}{\mathrm{K}}_{\mathrm{D}}\mathrm{S}-\mathrm{\Gamma }\Vert {\mathrm{Y}}_{or}^{\mathrm{T}}\mathrm{S}\Vert +\Vert {\mathrm{S}}^{\mathrm{T}}{\mathrm{Y}}_{or}\Vert \Vert {\alpha }_{or}\Vert$

(76)

Taking transpose of ‖S^TY_or‖ and also Γ ≥ ‖α_or‖ gives

${\stackrel{.}{V}}_1(\mathrm{t},\mathrm{S})\le -{\mathrm{S}}^{\mathrm{T}}{\mathrm{K}}_{\mathrm{D}}\mathrm{S}$

(77)

It is known that [36] K_D = M_0rκ, where κ can be considered as a least eigenvalue. Hence, (77) can be rewritten as

$\frac{\mathrm{d}{\mathrm{V}}_1(\mathrm{t},\mathrm{S})}{dt}\le -{\mathrm{S}}^{\mathrm{T}}{\mathrm{M}}_{0\mathrm{r}}\kappa \mathrm{S}$

(78)

Using (73), (78) can be rewritten as

$\frac{\mathrm{d}{\mathrm{V}}_1(\mathrm{t},\mathrm{S})}{dt}\le -2\kappa {\mathrm{V}}_1(\mathrm{t},\mathrm{S})$

(79)

The solution of above equation is

${\mathrm{V}}_1(\mathrm{t},\mathrm{S})\le {\mathrm{V}}_1(0,\mathrm{S}(0)){\mathrm{e}}^{-2\kappa \mathrm{t}}$

(80)

It is evident from the above equation that the sliding surface will converge exponentially to zero. Thus the sliding surface is related to the tracking error e_r in (62) which also converges exponentially to zero which satisfies Tikhnov's theorem.

2) Exponential Stability of Fast Subsystem: Tikhnov's theorem requires that the fast subsystem must also be exponentially stable. The energy multiplier method used in [28] is followed to prove the exponential stability under the following theorem [37], [theorem 4.1] which guarantees exponential stability.

Theorem: Let A be the infinitesimal generator of a C_o semigroup T(t). If for some p, 1 ≤ p ≤ ∞

${\displaystyle {\int }_0^{\infty }{\Vert \mathrm{T}(\mathrm{t})\Vert }^{\mathrm{p}}dt<\infty }$

(81)

then there are constants M ≥ 1 and μ > 0 such that ‖T(t)‖‖ ≤ Me^–μt

Note: It is also shown in [28] and [38] that the property of L^p stability and exponential stability for a strongly continuous semigroup must satisfy the following

${\displaystyle {\int }_0^{\infty }{\Vert \mathrm{E}(\mathrm{v})\Vert }^{2}dv<\infty }$

(82)

Proof: Let the energy function for (67) be of the form

$\mathrm{E}(\mathrm{v})=\frac{1}{2}{\Vert \mathrm{A}{\mathrm{w}}_{\mathrm{f}}(\mathrm{v})\Vert }^2+\frac{1}{2}{\Vert {\mathrm{A}}^{\frac{1}{2}}{\widehat{w}}_{\mathrm{f}}(\mathrm{v})\Vert }^2$

(83)

where E(v) is weakly monotonically decreasing function with respect to fast time scale v [38]. The fast time scale derivative of E(v) from the above equation will be

$\widehat{E}(\mathrm{v})=-\mathrm{k}〈QA{\widehat{w}}_{\mathrm{f}}(\mathrm{v}).\mathrm{A}{\widehat{w}}_{\mathrm{f}}(\mathrm{v})〉\le 0$

(84)

Let us choose 0 < ε < 1 and the Lyapunov function candidate is given by

${\mathrm{V}}_2(\mathrm{v})=2(1-\epsilon )vE(\mathrm{v})+〈{\widehat{w}}_{\mathrm{f}}(\mathrm{v}).\mathrm{A}{\widehat{w}}_{\mathrm{f}}(\mathrm{v})〉$

(85)

We have the following relation

$〈{\widehat{w}}_{\mathrm{f}}(\mathrm{v}).\mathrm{A}{\widehat{w}}_{\mathrm{f}}(\mathrm{v})〉\le \frac{1}{2}\left({\Vert \mathrm{A}{\mathrm{w}}_{\mathrm{f}}(\mathrm{v})\Vert }^2+{\Vert {\mathrm{A}}^{\frac{-1}{2}}\Vert }^2{\Vert {\mathrm{A}}^{\frac{1}{2}}{\widehat{w}}_{\mathrm{f}}(\mathrm{v})\Vert }^2\right)$

There exists a constant C₁ such that

$\left[2(1-\epsilon )\mathrm{v}-{\mathrm{C}}_1\right]\mathrm{E}(\mathrm{v})\le {\mathrm{V}}_2(\mathrm{v})\le \left[2(1-\epsilon )\mathrm{v}+{\mathrm{C}}_1\right]\mathrm{E}(\mathrm{v})$

(86)

For v > v₁ the Lyapunov function is positive and v₁ is from

$2(1-\epsilon ){\mathrm{v}}_1-{\mathrm{c}}_1=0$

(87)

The derivative of V₂(v) with respect to the fast time scale is given by

$\begin{array}{r}{\widehat{V}}_2(\mathrm{v})=(2-\epsilon ){\Vert {\mathrm{A}}^{\frac{1}{2}}{\widehat{w}}_{\mathrm{f}}(\mathrm{v})\Vert }^2-\epsilon {\Vert \mathrm{A}{\mathrm{w}}_{\mathrm{f}}(\mathrm{v})\Vert }^2\\ \hspace{0.17em}\hspace{0.17em}-2\mathrm{k}(1-\epsilon )\mathrm{v}〈QA{\widehat{w}}_{\mathrm{f}}(\mathrm{v}).\mathrm{A}{\widehat{w}}_{\mathrm{f}}(\mathrm{v})\\ \hspace{0.17em}-\mathrm{k}〈QA{\widehat{w}}_{\mathrm{f}}(\mathrm{v}).\mathrm{A}{\mathrm{w}}_{\mathrm{f}}(\mathrm{v})〉\end{array}$

(88)

For any arbitrary constant, say, C₂ > 0 we have

$\begin{array}{l}-〈QA{\widehat{w}}_{\mathrm{f}}(\mathrm{v}).\mathrm{A}{\mathrm{w}}_{\mathrm{f}}(\mathrm{v})〉\le \frac{1}{2}{\lambda }_{\max }(\mathrm{Q})[{\mathrm{C}}_2^2({\Vert \mathrm{A}{\widehat{w}}_{\mathrm{f}}(\mathrm{v})\Vert }^2\\ +\frac{1}{{\mathrm{C}}_2{}^2}{\Vert \mathrm{A}{\mathrm{w}}_{\mathrm{f}}(\mathrm{v})\Vert }^2\end{array}$

(89)

Using (89), (88) can be rewritten as

$\begin{array}{r}{\widehat{V}}_2(\mathrm{v})=[(2-\epsilon ){\Vert {\mathrm{A}}^{-\frac{1}{2}}\Vert }^2+\frac{{\mathrm{C}}_2^2\mathrm{k}}{2}{\lambda }_{\max }(\mathrm{Q})\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}-2kv(1-\epsilon ){\lambda }_{\min }(\mathrm{Q})]{\Vert \mathrm{A}{\widehat{w}}_{\mathrm{f}}(\mathrm{v})\Vert }^2\\ \hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}\hspace{0.17em}-(\epsilon -\frac{\mathrm{k}}{2{\mathrm{C}}_2^2}{\lambda }_{\max }(\mathrm{Q})){\Vert \mathrm{A}{\mathrm{w}}_{\mathrm{f}}(\mathrm{v})\Vert }^2\end{array}$

(90)

where ${\lambda }_{\max }(\mathrm{Q})={\max }_{{\mathrm{w}}_{\mathrm{f}}\in \mathrm{H}}〈\mathrm{Q}{\mathrm{w}}_{\mathrm{f}}.{\mathrm{w}}_{\mathrm{f}}〉$ and ${\lambda }_{\min }(\mathrm{Q})={\min }_{{\mathrm{w}}_{\mathrm{f}}\in \mathrm{H}}〈\mathrm{Q}{\mathrm{w}}_{\mathrm{f}}.{\mathrm{w}}_{\mathrm{f}}〉.$ .

If C₂ can be chosen as large then $(\epsilon -\frac{\mathrm{k}}{2{\mathrm{C}}_2^2}{\lambda }_{\max }(\mathrm{Q}))>0$

${\widehat{V}}_2(\mathrm{v})\le 0\hspace{0.17em}\hspace{0.17em}\forall \mathrm{v}>{\mathrm{v}}_2$

(91)

where v₂ can be

$(2-\epsilon ){\Vert {\mathrm{A}}^{-\frac{1}{2}}\Vert }^2+\frac{{\mathrm{C}}_2^2\mathrm{k}}{2}{\lambda }_{\max }(\mathrm{Q})-2\mathrm{k}{\mathrm{v}}_2(1-\epsilon ){\lambda }_{\min }(\mathrm{Q})=0$

The above result in (91) shows that derivative of Lyapunov function has a decreasing nature for v > v₂, and it is also evident from (84) that the energy will also be dissipating for v > 0. Using these facts, for v > T ≔ max{v₁,v₂} and also from (86) E(v) can be estimated as

$\mathrm{E}(\mathrm{v})\le \frac{\mathrm{V}(\mathrm{T})}{2(1-\epsilon )\mathrm{v}-{\mathrm{C}}_1}\le \frac{[2(1-\epsilon )\mathrm{T}+{\mathrm{C}}_1]\mathrm{E}(0)}{2(1-\epsilon )\mathrm{v}-{\mathrm{C}}_1}$

(92)

Then, (92) can be rewritten as

${\displaystyle {\int }_{\mathrm{T}}^{\infty }\mathrm{E}{(\mathrm{v})}^{2}dv}\le {\displaystyle {\int }_{\mathrm{T}}^{\infty }\frac{{[2(1-\epsilon )\mathrm{T}+{\mathrm{C}}_1]}^2\mathrm{E}{(0)}^2}{2(1-\epsilon )\mathrm{v}-{\mathrm{C}}_1}}<\infty$

(93)

which confirms the exponential stability as given in (81) and hence it is proved.

By choosing the appropriate Lyapunov function candidate for the slow and fast subsystems, stability analysis has been carried out. It is evident from the above analysis that both the systems converge to zero exponentially which satisfies Tikhnov's theorem. Therefore, by using this composite control strategy, the flexible object will track the desired trajectory and simultaneously suppress the vibration.

Since the control law (63) and (64) is discontinuous across the sliding surface, such a control law leads to chattering. Chattering is undesirable in practice because it involves high control activity. To remedy this drawback, it usually uses $\frac{\mathrm{x}}{\left|\mathrm{x}\right|+\delta }$ to replace sgn(x) in the control law (63), where δ is x a constant. Let the switching function ψ in the control law (63) be replaced by

$\psi =-\mathrm{\Gamma }\frac{{\mathrm{Y}}_{or}^{\mathrm{T}}\mathrm{S}}{\Vert {\mathrm{Y}}_{or}^{\mathrm{T}}\mathrm{S}\Vert +\delta }$

(94)

where Γ becomes $\mathrm{\Gamma }\ge (1+\frac{\delta }{\epsilon })\Vert {\alpha }_{or}\Vert$ and ε is a design constant. Following the analysis given in [39], it can prove that e_r will exponentially converge to a small bound, which depends on selection of δ. Based on Tikhnov's theorem, the above stability analysis is still valid. However, in such a case, tracking error e_r will not tend to zero but is uniformly bounded. The admissible amplitude of the tracking error can be achieved by choosing a suitable δ.

7. Simulation studies

In the simulation, two identical planar rigid manipulators, each with three-links/revolute joints, are used to move the flexible beam along a desired trajectory specified by,

${\mathrm{X}}_{rfd}=\left[\begin{array}{c}\sin (\mathrm{t})\\ \cos (\mathrm{t})\\ 0\end{array}\right]$

The parameters of the identical manipulators [40] are given in Table I and the flexible beam parameters are given in Table II.

Table 1.

Parameters of the manipulator

Link	Length (m)	Mass (kg)	Moment of inertia (kgm2)
1	0.6	1.5	0.50
2	0.6	1.5	0.50
3	0.2	1.5	0.25

Table 2.

Parameters of the flexible beam

Parameter	Value
Mass (m)	1.5 kg
Length (L)	1 m
Density	2700kg/m3
Young's modulus (E)	71 GPa
Moment of Inertia (I)	0.13 kgm2

The beam initial position and orientation are X_rf = {0.51m,0.36m,0.1rad}^T and initial velocity and acceleration are considered to be zero. Initial joint angles of manipulator are q₁₁ = 0.2974 rad, q₁₂ = 1.6974 rad, q₁₃ = −1.6948 rad, q₂₁ = 0.2149 rad, q₂₂ = 1.4886 rad and q₂₃ = −1.4306 rad, respectively. The initial joint velocities of all the joints of manipulators are 0.001 rad/sec and joint accelerations are assumed to be zero. The simulation was carried out with a sampling period of 0.01 sec. The control parameters are given in Table III. The value of Γ was chosen based on the L_∞ norm of the time independent parameters of the regressor and in this case, vector α_or.

Table 3.

Control parameters

Parameter	Value
KD	diag(424.2)
λ	diag(7.9)
Γ	0.004
k	1

The regressor for the slow subsystem is very lengthy and due to space limitations it is not provided but can be made available upon request. Figs. 3 and 4 show the tracking of the planar motion of centre of the object along X, Y-directions. The orientation of the beam is shown in Fig. 5 and also the tracking error is shown in Fig. 6. It can be observed that the tracking of position and orientation is achieved within 1 sec, which shows the effectiveness of the controller. In order to achieve the desired trajectory of the object, both joints of the manipulators are moved in a smooth manner which are shown in Figs. 7 to 12. The sliding variables are shown in Fig. 13.

Figure 3.

X-position tracking

Figure 4.

Y-position tracking

Figure 5.

Orientation of the beam

Figure 6.

Tracking error

Figure 7.

Angular motion of joint 1 of manipulator 1

Figure 8.

Angular motion of joint 2 of manipulator 1

Figure 9.

Angular motion of joint 3 of manipulator 1

Figure 10.

Angular motion of joint 1 of manipulator 2

Figure 11.

Angular motion of joint 2 of manipulator 2

Figure 12.

Angular motion of joint 3 of manipulator 2

Figure 13.

Sliding variables

In the case of fast subsystem, the initial disturbance of 5 mm with zero initial velocity is considered for the simulation. Although, the beam dynamic equation is not assumed with the number of modes, simulation has to be performed with the assumption of modes. Normally, only the first few modes are dominant, yielding a higher vibration response and hence the first four modes are considered here. For the case of structural damping characteristics when β = −0.5, in Fig. 14, the vibration yields oscillatory motion that is completely suppressed at around 1 sec. It is also shown in Fig. 15 that exponential decay occurs at β = 0 which corresponds to the overdamping behaviour of the fast subsystem which is same for all modes. Simulations are also carried out for different damping ratios of 0.1 and 0.4. The results in Figs. 16 and 17 show that with an increasing damping ratio, the amplitude of vibration is significantly suppressed. Furthermore, the transverse displacement under the simply supported end condition of the beam with a damping ratio of 0.1 is evaluated at various locations of the beam using the modal summation method. Due to the symmetry of the boundary conditions, the deflection at 0.1 m from the left end and at the middle of the beam are considered. It can be seen from the Figs. 18 and 19 that the centre of the beam yields more deflection than any other point of the beam.

Figure 14.

Structural damping characteristics β=−0.5

Figure 15.

Strong damping characteristics β= 0

Figure 16.

Deflection of the beam for damping ratio=0.1

Figure 17.

Deflection of the beam for damping ratio=0.4

Figure 18.

Deflection at 0.1 m from the left end of the beam

Figure 19.

Deflection at mid point of the beam

Simulations were carried out with different masses of the beam and the results are shown in Figs. 20 and 21. It can be seen that the response is higher for lower masses. Results are also presented for different lengths of the beam in Figs. 22 and 23. It is seen that the responses are higher for longer beams.

Figure 20.

Mass of the beam 2 kg

Figure 21.

Mass of the beam 5 kg

Figure 22.

Length of Beam 1 m

Figure 23.

Length of Beam 2 m

8. Conclusion

In this article, the problem of two planar rigid manipulators used to move a flexible object along the desired trajectory while suppressing the vibration of the flexible object is considered. The kinematic and dynamic relations of the beam are derived without either approximating or discretizing the beam.

The complete system is then obtained by combining both the dynamics of the beam and the manipulators. To control this highly nonlinear system, the singular perturbation method is adopted, where the slow and fast subsystems are separately obtained. Based on Tikhnov's theorem, a regressor-based sliding mode control algorithm is designed for the slow subsystem, and a simple feedback control law is selected for the fast subsystem to form a composite controller. The simulation results show that the proposed composite controller is an efficient choice in tracking the desired trajectory while suppressing the vibration of the flexible object.