Nonlinear dynamic analysis and state space representation of a manipulator under viscoelastic material conditions

 

 

H. Esfandiar^{I, 1, *}; S. Daneshmand^{II, 2}

^IDepartment of Mechanical Engineering, Firuzkooh Branch Islamic Azad University, Firuzkooh, Iran. habib_esfandiar_2006@yahoo.com
^IIDepartment of Mechanical Engineering, Majlesi Branch Islamic Azad University, Isfahan, Iran. S.daneshmand@iau maj lesi.ac.ir

 

 

---

ABSTRACT

In this paper, based on the Voigt-Kelvin constitutive model, nonlinear dynamic modelling and state space representation of a viscoelastic beam acting as a flexible robotic manipulator is investigated. Complete nonlinear dynamic modelling of a viscoelastic beam without premature linearisation of dynamic equations is developed. The adopted method is capable of reproducing nonlinear dynamic effects, such as beam stiffening due to centrifugal and Coriolis forces induced by rotation of the joints. Structural damping effects on the model's dynamic behaviour are also shown. A reliable model for a viscoelastic beam is subsequently presented. The governing equations of motion are derived using Hamilton's principle, and using the finite difference method, nonlinear partial differential equations are reduced to ordinary differential equations. For the purpose of flexible manipulator control, the standard form of state space equations for the viscoelastic link and the actuator is obtained. Simulation results indicate substantial improvements in dynamic behaviour, and a parameter sensitivity study is carried out to investigate the effect of structural damping on the vibration amplitude.

---

OPSOMMING

'n Voigt-Kelvin saamgestelde model word in 'n nie-reglynige dinamiese oplossingsruimte gebruik om 'n soepel hanteringsrobot te ontleed. Die modellering is daartoe instaat om nie-linêre dinamiese vertoning van styfheid wat ontstaan wees sentrifugale en Corioliskragte te ontleed. Die modelvergelykings is gebaseer op die Hamiltonbeginsel. Simulasieresultate toon betekenisvolle verbetering van sisteemgedrag. Verdere ontledings handel oor strukturele demping en vibrasie-amplitude.

---

 

 

1. INTRODUCTION

Many materials exhibit simultaneous elastic and viscous behaviour. Viscoelastic materials are widely used in vehicle technology, spacecraft attachments, and composite material development to reduce undesired vibrations. Developing realistic flexible multi-body system models requires the inclusion of the effect of damping that characterises most existing materials. Viscoelastic materials, when vibrating, dissipate energy; and as a consequence, the effect of initial disturbances or shock loadings dies out after a certain period of time depending on the amount of material damping. The inclusion of damping in computational models may also be necessary in order to damp out high frequency oscillations that do not have a significant effect on the solution. Dynamic analysis of viscoelastically damped systems is not easy to implement in engineering applications. Accurate modelling of a material's behaviour is essential to the study and accurate design of the correspondent system. Viscoelastic materials have been characterised by classical linear integer-order derivative models such as the Kelvin-Voigt and Maxwell models [1]. Many published research articles in the recent literature cover areas of dynamic modelling, vibration control, and stability analysis of flexible beams. In spite of the flexibility of a robotic manipulator, many modelling methodologies - analytical, assumed mode, finite element, and finite difference methods - have been proposed for the flexible manipulators. Theodore et al. [2], Kopacek et al. [3], and Ge et al. [4] compared the assumed mode method and the finite element method for modelling the flexible multi-link manipulators. To solve a large set of differential equations derived using the finite element method, a lot of boundary conditions should be considered; these are, in most situations, uncertain for the flexible manipulator. Dynamic characteristics of rotating beams have been investigated by Chung and Yoo [5]. They formulated a three-dimensional vibration of a rotating beam using a finite element method (FEM), and obtained the time responses and distribution of the deformations and stresses at a given rotating speed. The stability analysis for the flap-wise motion of a cantilever beam, subjected to a rotary oscillation, was studied by Chung et al. [6], who obtained the transition curves for a sinusoidal varying angular velocity. Similar analysis was carried out for a multi-layered composite beam by Yoo et al. [7]. The vibration and instability analyses of a rotating non-uniform beam with frequency-dependent structural damping, equivalent viscous damping, and root damping were studied by Lin and Leen [8]. The dynamic behaviour and instability of a rotating sandwiched beam with a constrained damping layer, subjected to an axial periodic load, using FEM, was studied by Lin and Chen [9].

A new rotating beam finite element model was developed by Gunda and Ganguli [10] in which the basic functions were obtained by an exact solution of the governing static homogeneous differential equation of a stiff string. The method proposed by them was found to be applicable to non-uniform rotating beams. Abolghasemi and Jalali [11] obtained the attractors of a rotating viscoelastic beam using a perturbation technique, and found that stable limit-cycles exist for the system in both the presence and the absence of structural damping. Dynamic modelling and control of vibration for a uniform rotating cantilever beam was carried out by Yang et al. [12]. The transverse vibration of the rotating beam was suppressed using a combinational technique that included positive position feedback control and the momentum exchange feedback method. Vibration behaviour and control of a clamped-free rotating flexible cantilever arm, with fully covered active constrained layer damping treatment using FEM, was investigated by Fung and Yau [13]. Dasa et al. [14] investigated the free out-of-plane vibration of a rotating beam with a non-linear spring mass system using a perturbation method. The non-linear coupled longitudinal and bending vibration of a variable speed rotating beam was modelled by Younesian and Esmailzadeh [15]. The solution of the nonlinear partial differential equation of motion was assumed to be a series of multiplications of time and position functions in accordance with the Galerkin modal summation technique. The multiple scales method was used to obtain the exact first-order solution. Roy et al. [16] studied the dynamic behaviour of a composite beam, formed by concentric layers of steel and aluminium. Two efficient modelling techniques used coupled (thermo-mechanical) ATF (Augmenting Thermodynamic Fields) and ADF (Anelastic Displacement Fields) displacements, to represent the constitutive relationship in time domain by using certain viscoelastic parameters. Dynamic investigations of an axially moving beam subject to a constant axial stress were carried out by Marynowski [17]. The Kelvin-Voigt and Maxwell models were considered as beam material models. The numerical solutions of full nonlinear and linearised equations were compared. The effects of axial travelling speed and internal damping on the dynamic stability of the axially moving beam were studied in detail. Temel [18] examined the transient analysis of viscoelastic helical rods subject to time-dependent loads in the Laplace domain. The governing equations for naturally twisted and curved spatial rods obtained using the Timoshenko beam theory were rewritten for cylindrical helical rods. The curvature of the rod axis, the effect of rotary inertia, and shear and axial deformations were considered in the formulation. The material of the rod was assumed to be homogeneous, isotropic, and linear viscoelastic.

This paper is organised follows. In Section 2, a mechanical description of the case study is presented; in Section 3, mathematical modelling of a viscoelastic beam is presented; in Section 4, the discretisation procedure for the continuous model of a viscoelastic beam is described; in Section 5 numerical results are presented; in Section 6 an actuator system dynamic model is presented; in Section 7 state space representation for the overall system, including the viscoelastic arm and actuator, is derived; and the last section considers the conclusions.

 

2. MECHANICAL DESCRIPTION

The coordinate systems used in developing the model are shown in Figures 1 and 2 where is the inertial reference frame, is a body coordinate system attached to the motor shaft, and is another body coordinate system attached to the hub such that its x₂ axis is directed along the initial situation of the arm. The rotation of frame relative to frame is described by the angle ^θ(t), and ^α denotes rotation of framerelative to frame.

 

 

 

 

Consider a viscoelastic beam as an Euler-Bernoulli beam with length L fixed on the hub in the horizontal plane. The joint is taken to have small torsional deformations. The hub is assumed to be rigid, and the viscoelastic link is radially attached to the hub. The motion of the viscoelastic rotating beam is described by the rigid body rotation (θ), the torsional deformation angle ^(α) that is measured relative to the motor coordinate system, and transverse displacement w(x,t) of the beam. The rotatory inertia is also considered. The beam is divided into N - 1 elements with the ⁱth element having two nodes i and i+1, and the transverse component of the displacement vector of an infinitesimal link element is described by ^{w(x t)}.

 

3. MATHEMATICAL MODELLING

In this section, we will consider dynamic modelling of a viscoelastic beam without joint flexibility with a clamped-free condition together with a Voigt-Kelvin constitutive model. The position of an infinitesimal element of the link located at a distance x from the frame origin in the X₂, Y₂ direction, relative to the {O₀,X₀Y₀Z₀^} reference frame, is given by Martins et al. [19]:

where e₁, e₂, e₃ are the unit vectors along the X₀, Y₀ axis and Z₀ axis respectively.

According to this expression, the kinetic energy involves kinetic energy of the hub and of the flexible link as follows:

wheredenotes time derivative of (_-) and (-).(-) denotes the dot product of vectors (-) .The potential energy is composed of internal energy due to the elastic deformation of the link from bending.

Potential energy due to gravity is not accounted for, since only the motion in the plane perpendicular to the gravitational field is considered. The non-conservative work is composed of the structural damping and torque applied at the hub. One could express the stress-strain relationship for a Kelvin-Voigt standard viscoelastic material as follows [20,21]:

in which _Eηs are respectively the module of elasticity and the damping factor and time derivative operator, and σ,ε represent the stress and strain values respectively. In order to arrive at the solution method, the stress-strain relationship is now converted into the following form

where for the Kelvin-Voigt standard viscoelastic material P = 1 and Q = E + η _s. The differential equation of the beam under viscoelastic material conditions without rotation of the hub is presented. The bending moment at each section of the beam is obtained as

Multiplying the above equation by ρ gives

where

in which ρ is the curvature radius due to bending. Using the above equation

where I and A respectively represent moment of inertia and the cross-sectional area of the beam. After two derivations of the last equation, this follows:

The non-conservative work done by the structural damping is given as

where _σd represents the stress of the viscous material of the Kelvin-Voigt model.

Using Eq. (12), the non-conservative work done by the Kelvin-Voigt structural damping is given by

and the non-conservative work done by the torque applied at the hub is

The Hamilton principle is used in conjunction with the known theorems of the calculus of variations to obtain the non-linear equation of motion for the Kelvin-Voigt model as follows:

where Eq. (15) is the system equation of motion, and Eq. (16) is the constraint equation expressing the total balance of angular momentum of the system. The first set of Eq. (17) expresses the two geometric boundary conditions at the clamped end of the link, and the second set presents natural boundary conditions at the free end of the link. Term is centrifugal force, is force due to tangential acceleration, is force due to Kelvin-Voigt structural damping model, is torque due to the Coriolis force, and is torque due to the force caused by tangential acceleration.

 

4. DISCRETISATION OF THE EQUATIONS OF MOTION

The continuous models consist of partial differential equations that are difficult to deal with both analytically and computationally. To obtain a set of ordinary differential equations, a discretisation of the equations developed previously must be performed. The finite difference method is used to approximate partial differential equations with a set of ordinary differential equations. In this method, the interval is divided into n-1 uniform segments with We substitute space derivatives with finite differences so that the second and fourth spatial derivatives are substituted by the second-order central difference approximation and the first spatial derivative by the first-order backward difference approximation.

And the third spatial derivative is approximated by

where abbreviated notations w_i(t) and W_o(t) are used to represent w(x.j) and w(rj) respectively. It must be noted that the fourth order derivative in every point depends on the information gathered from five points; therefor meshing must have at least four elements or five nodes. For discretisation of expressions that consist of integrals, the integral is replaced with summation and the beam is meshed from 1 to n 1 intervals such that the first node in the mesh is related to the hub and equals 0, and the final node relates to x = r + L and equals n.

In this section, equations of the Kelvin-Voigt model are considered and the set of ordinary differential equations are presented. By using the finite difference method, ordinary differential equations are expressed as follows:

Since the boundary value of the solution is obtained by substituting i = 0, this approximation is applied for i = 1,2,..., n, where

for i=1

Discretisation of the partial differential equation with respect to the spatial variable is performed by substituting a set of ordinary differential equations that can be written in matrix form as follows [19]:

where

where z₁ is the variables' vector with n +1 elements expressed as follows:

 

5. SIMULATION RESULTS

In order to illustrate the behaviour of the rotating viscoelastic beam with the Kelvin-Voigt standard model under excitation pulse, a beam with the geometrical and material parameters presented in Table 1 is considered. The pulse has a duration of two seconds and an amplitude of 0.55N_m respectively.

The nonlinear differential equation of motion (Eq.28) is numerically solved using the central difference method, and a computer program is developed using MATLAB software for numerical simulations. A parametric sensitivity study is also carried out to evaluate the effects of structural damping on time responses. Kinematic parameters of a rotating viscoelastic beam that demonstrates a pure elastic system are illustrated in Figure 4.

 

 

 

 

As shown in Figure 4, even for a pure elastic system (η = 0), the damping matrix is not null, and the centrifugal force produces oscillation responses for angular displacement, velocity, and tip deflection of the elastic beam. Kinematic responses of the viscoelastic beam are oscillatory, but they do not have a regular form. Responses after the end of the input pulse show sudden variation. It can easily be found that this event is relevant to the excitation pulse. It is also interesting to note that the centrifugal force behaves as a virtual damper, and so it reduces vibration amplitudes with an increase in time. Figure 5 illustrates angular displacement, angular velocity, and transverse displacement of the free end for the viscoelastic material system. It indicates that results obtained by different values of the damping factor are almost the same, but a lower damping factor produces larger vibration amplitudes compared with a higher damping factor. It is clear that the rate of decrease of vibration amplitude is increased with an increase in the damping value; also, the response of the viscoelastic system is more regular than the response of the pure elastic system. In this case, the centrifugal force behaves as a virtual damper as well. It can be seen that both centrifugal force and structural damping reduce vibration amplitudes considerably, but the damping factor has a higher reduction effect.

 

 

The effect of damping on the vibration amplitudes for a fixed time is illustrated in Figure 6. It can be seen that by increasing the viscosity factor up to N = 35 , displacement of the end effector of the beam almost decreases, but beyond this value, the amplitude of vibration increases. It was thus found that if the value of the damping is higher than 35, the tip deflection of a viscoelastic beam cannot be stable.

 

 

6. ACTUATOR SYSTEM DYNAMIC MODEL

The actuator dynamics are incorporated with the link state space equation into the overall system. A DC motor with a constant field - assumed to be the actuator for this study - is driven by the application of voltage to its armature terminals. Through a set of gears, the motor drives a load with moment of inertia J, subject to an external torque. Its dynamics are described by

Where T_e and T_m are the torque exerted on the motor shaft by the load, and transmitted respectively through the gears and produced torque by a DC motor. Thus

So the torque at the motor shaft is seen to be the torque generated by the motor, minus the torque required to accelerate the rotor. The torque exerted by the motor on the load shaft and transmitted through the gears, is NT_e. Newton's second law applied to the load is

Where ν is the gear ratio of the motor. Since this becomes

or

where J_e = J + ν²J_m is effective inertia seen at the load shaft. We need the current i by the Kirchhoff voltage law

where L' and R are the inductance and the resistance of the armature circuit. k_mω is the back emf and the equation in the matrix form is as follows:

By selection of θ and ω as the outputs, the output equation is

The simulation conditions are as follows: zero conditions for the state variables and

The time response of the DC motor is shown in Figure 7.

 

 

7. STATE SPACE REPRESENTATION

For the purpose of flexible manipulator control, we need a standard form of the state space equations for a viscoelastic link and actuator. The finite difference method is used to approximate partial differential equations with a set of ordinary differential equations as follows:

And from the previous section, the state space equation of the actuator is derived. In this section, state space equations of the total system are aimed at. To obtain the closed form of the state space equations of the system, state space representation of the DC motor (Eq. (37)) should be transformed into an equation of motion for the flexible beam of the same order.

Introducing variables with n+1 elements, i is the current of the motor and ; expansion of the state space representation of the driver motor is as follows:

With this separation,

Now the state space representation of equations of motion for the viscoelastic beam maintaining the total form of Eq. (39) must be extracted; where,

Where 0|M-¹K is the M^-1K matrix with one zero column added to it until its dimensions have become (n+1) _x (n+2). Using Eqs. (41) and (42), the total description of the state space equations of the overall system is:

And the closed form is

Selecting the angular position of the viscoelastic beam and the tip deflection of the beam as output variables, for the output equation,

Only Y (1,2) = Y (2, n + 2) = 1. Clearly, the state equation (Eq. (44)) has incorporated manipulator motion, including the hub rotation and link deflection and actuator dynamics, as well as sensor specification. The system matrix a is a function of link construction (mass and stiffness distributions) and actuator design (motor parameters), while output matrix C is a function of beam position and tip deflection.

 

8. CONCLUSION

In this paper, nonlinear dynamic modelling and state space approximation of a viscoelastic rotating beam was modelled and analysed; a complete dynamic model of a flexible beam without premature linearisation of dynamic equations was subsequently developed, while nonlinear dynamic effects, such as beam stiffening due to centrifugal and Coriolis forces induced by the rotation of the joints, were reproduced, and structural damping effects on the model's dynamic behaviour were investigated. The governing equations of motion were derived using the Hamilton principle, and the finite difference method was used to approximate partial differential equations of motion with a set of ordinary differential equations. Numerical simulations were carried out for both pure elastic and viscoelastic systems under Voigt-Kelvin material conditional systems. Simulating pure elastic systems showed that the centrifugal force behaves as a virtual damper, thus reducing vibration amplitudes with an increase in time. Viscoelastic simulation results indicated that the results obtained by different values of damping factor are almost the same; but for lower values of damping factor, the vibration amplitude was larger than when using higher values, and the response of the viscoelastic system was more regular than the response of the pure elastic system. The effect of structural damping on the vibration amplitude was investigated in a parametric sensitivity study. It was found that while increasing the damping factor, the amplitude of vibration almost decreased until the damping factor reached its critical value. Finally, for the purposes of control, the state space equations of a viscoelastic rotating beam and actuator were obtained.

 

REFERENCES

[1] Hughes, T.J.R. & Simo, J.C. -1998. -Computational inelasticity, Springer, New York.         [ Links ]

[2] Theodore, R.J. & Ghosal, A. -1995. -Comparison of the assumed modes method and finite element models for flexible multi link manipulators, Journal Robotics Research, 14 (2), pp 91111.         [ Links ]

[3] Kopacek, P., Desoyer, K. & Lugner, P. -1989. -Modeling of flexible robots- An introduction, IFAC Proc. Series (10), pp 21-27        [ Links ]

[4] Ge, S.S., Lee, T.H. & Zhu, G. -1995. -Comparison studies between assumed modes method and finite element method in modeling a single-link-flexible manipulator, Proceedings of the IEEE International Conference on Intelligent Control and Instrumentation, Singapore, pp 376-381.         [ Links ]

[5] Chung, J. & Yoo, H.H. -2002. -Dynamic analysis of a rotating cantilever beam by using the finite element method. Journal Sound Vibration, 249, pp 147-164.         [ Links ]

[6] Chung, J., Jung, D. & Yoo, H.H. -2004. -Stability analysis for the flap wise motion of a cantilever beam with rotary oscillation, Journal Sound Vibration, 273, pp 1047-1062.         [ Links ]

[7] Yoo, H.H., Lee, S.H. & Shin, S.H. -2005. -Flap wise bending vibration analysis of rotating multi- layered composite beams. Journal Sound Vibration, 286, pp 745-761.         [ Links ]

[8] Lin, S.M. & Leen, S.Y. -2004. -Prediction of vibration and instability of rotating damped beams with an elastically restrained root, Journal Mechanical Science, 46, pp 1173-1194.         [ Links ]

[9] Lin, C.Y. & Chen, L.W. -2003. -Dynamic stability of a rotating beam with a constrained damping layer,Journal Sound Vibration, 267, pp 209-225.         [ Links ]

[10] Gunda, J.B. & Ganguli, R. -2008. -Stiff-string basis functions for vibration analysis of high speed rotating beams, Journal Applied Mechanics, 75, pp 024502-1-024502-5.         [ Links ]

[11] Abolghasemi, M. & Jalali, M.A. -2003. -Attractors of a rotating viscoelastic beam,Journal Non Linear Mechanics, 38, pp 739-751.         [ Links ]

[12] Yang, J.B., Jiang, L.J. & Chen, D.C.H. -2004. -Dynamic modeling and control of a rotating Euler- Bernoulli beam, Journal Sound vibration, 274, pp 863-875.         [ Links ]

[13] Fung, E.H.K. & Yau, D.T.W. -2004. -Vibration characteristics of a rotating flexible arm with ACLD treatment, Journal Sound Vibration, 269, pp 165-182.

[14] Dasa, S.K., Rayb, P.C. & Pohit, G. -2007. -Free vibration analysis of a rotating beam with nonlinear spring and mass system, Journal Sound Vibration, 301, pp 165-188.         [ Links ]

[15] Younesian, D. & Esmailzadeh E. -2010. -Non-linear vibration of variable speed rotating viscoelastic beams, Nonlinear Dynamics, 60, pp 193-205.         [ Links ]

[16] Roy, H. Dutt, J.K. & Datta, P.K. -2009. -Dynamics of multilayered viscoelastic beams, Structural Engineering and Mechanics, 33(4), pp 391-406.         [ Links ]

[17] Marynowski, K. -2002. -Nonlinear dynamic analysis of an axially moving viscoelastic beam, Journal of Theoretical and Applied Mechanics, 40.         [ Links ]

[18] Temel, B. -2004. -Transient analysis of viscoelastic helical rods subject to time-dependent loads, Journal of Solids and Structures, 41, pp 1605-1624.         [ Links ]

[19] Martins, J., Ayala, M. & Costa, J.-2002, -Modeling of flexible beams for robotic manipulators, Multibody System Dynamics, 7, 79-100.         [ Links ]

[20] Esfandiar, H & Daneshmand, S. -2012. -Complete dynamic modeling and approximate state space equations of the flexible link manipulator, Journal of Mechanical Science and Technology, 26(9), pp 2845-2856.         [ Links ]

[21] Rodriguz, M. -2006. -Analysis of structural damping, Master's thesis, Lulea University of Technology.         [ Links ]

 

 

* Corresponding author
1 The author was enrolled for an M Eng degree in the Department of Mechanical Engineering, Firuzkooh Branch, Islamic Azad University, Firuzkooh, Iran
2 The author was enrolled for an M Eng degree in the Department of Mechanical Engineering, Majlesi Branch, Islamic Azad University, Isfahan, Iran.