Comparison of Descriptor Models and Reduced Dynamic Models for Constrained Robots

Copyright © IFAC Robot Control. Vienna. Austria. 1991

MODELLING AND IDENTIFICATION

COMPARISON OF DESCRIPTOR MODELS AND REDUCED DYNAMIC MODELS FOR CONSTRAINED ROBOTS K.A. Tahboub*, Th. Schmidt, R. Schiipphaus and P.C. Miiller Safety Control Engineering. University ofWuppertal. p.a. Box 100127. D-5600 Wuppertal1. Germany

·On leave from the Hcbron Technical Engineering College, West Bank

Abstract. There are numerous applications where the end effector of a robot manipulator comes in contact with a rigid external object. These applications include the use of manipulators for carrying out assemb ling and machining tasks. In this paper we study two methods of modeling constrained robots: t he descriptor model where the constrained robot is co nsidered as a system of differential eq uation s with a lgebraic constraints and the reduced dynamic model where the minimal number of variab les. velocity and contact force, are used. The models were derived for a PUII·IA 560 robot in context with two benchmark tasks. Different spaces, the joint space and the task space, were considered for defi ning the model in the most suitable form. A comparitive study based on the dynamic simulation of the two methods in terms of the model complex ity, the required computational cost and the fl exibility of varying the task is presented.

Keywords . Robots; force control; modeling; descriptor systems; simulation ; assembling; machining.

1

2

Introduction

There are numerous applications where the end effector of a robot manipulator comes in contact with a rigid external object. These applications include the use of manipulators for carrying out assembling and machining tasks. A convenient coordinate frame for the desc ription of such tasks is the task frame which can be stationary as in the case of assembli ng tasks or non stationary as in the case of machining tas ks. Two methods of modeling constrained robots are pre· sented in this paper. The first is the descriptor model where the constrained robot is considered as a system of differential equations with algebraic constraints and the second is the redu ced dynamic model which ben efi ts from the orthogonality between the constraint force and the velocity directions in the task space. These two methods are compared with regard to two benchmark tasks corresponding to assembling and machining where a model of the six degree of freedom PUIIIA 560 robot is used for the simulation. The paper is organized as follows: section 2 introduces the different spaces which are normally considered in planning and controlling constrained tasks. The two benchmark tasks under study, namely, assembling and machining, are stated and illustrated in section 3. Section 4 introduces t he descriptor model and a method for finding the constraint equations while section 5 presents the redu ced dynamic model method. The possible numerical integration met hods for simulating the different models are discussed in sect ion 6; the results of the simulations are presented and interpreted in section 7 and th e final conclusions are given in section 8.

Task Frame, Reference Cartesian Frame, and Joint Space

When the probl em of constrained robots is considered , one talks about the motion of the end effector and the forces of co ntact b etween it and the constraining surface. The natural and logi cal way to specify these motions and forces is to as· sign values to the velocity (translational and rotational) and force (force and moment) vectors in the task frame where either a velocity or a force component is assigned along and around each axis of this frame. This means that the constraints on the motions and the directions of the forces are set in this frame ; this will be shown in section 3 through two examples . The robot equations of motion are normally written in the joint space in the form

M(q)q + C(q, q) =

r

(1)

where M( q) is the square mass matrix, C( q, q) represents the Coriolis, centrifugal, and gravity vector, q is the vector: of joint coordinates , and r is the vector of the generalized input control forces at the joints. One is normally interested in the resulting Cartesian motion beside that in the joint space. A suitable Cartesian frame is the frame of reference which is attached to the origin of link zero of the robot (the base of the robot). The relation between the joint position and the position and orientation of the end effector (written in the frame of reference) is gi"en by the forward and the inverse kinematic equat ions, whereas the relation between the joint velocities and the end effector translational and rotational velocities

9

is given by the basic lacobian matrix J(q)

Table 1: The Natural and Artificial Constraints for the Peg into the Hole Task

[ ~ ] = J(q)q

(2)

where V and n are the translational and rotational velocity vectors correspondingly. If constraint forces appear, then equation (1) becomes

M(q)q + C(q,q) =

r + JT(q) [ ~:

]

Natural constraints

Artificial constraints

It

=

Vt

Vb

=

In = 0.1 Ib = 0.2

Tt

=0

W,

= - 0.003(1 2

=

Tn

= 0.01

0 =0

Vn

(3)

\.4.'n

where F. and r. are the constraint force and torque vectors between the end effector and the constraining environment mapped to the joint space through the transpose of the la· cobian matrix. The joint acceleration q in (~) can be replaced by the end effector acceleration vector X where

0

0

= -0.0012(1 2

-

-

101)

101)

A General Machining Task

3.2

r..lachining tasks in general, require moving the cutting tool which is hold by the end effector on a surface while applying the necessary normal forces. To guarantee covering general desired motions, the differential geometry techniques were used in [6] to model the motion along a space curve. It was shown there how to find the natural and artificial constraints given the desired equation of motion. Here, the selected task consists of machining the internal surface of a big cylinder placed around the robot as shown in Fig. 2. The desired motion of the tool draws a helix

(4) This results in the equivalent form

(5) where

(6)

r(l) = OAcos(l)i

+ OAsin(l)j + 0.031k.

(8)

and

Cx(q,q) = C(q,q) - Mx(q)j(q)q.

(7)

The corresponding natural and artificial constraints are given in Table 2. The task space is found from the desired motion

The dynamic model in form (5) provides an effective means to describe the dynamic behaviour of the end effector and to integrate the position and force control as will be seen in section 5.

3

Table 2: The Natural and Artificial Constraints for a Machining Task

Benchmark Tasks

In this section, two benchmark tasks will be stated and ex· plained, they will serve as the core of the simulation in section 7. The two tasks are selected in such a way to represent a wide variety of robot tasks with constraints.

Natural constraints

Artificial constraints Vt

It

= 0

Vn

=

Vb

= 0

In = 0.1 Ib = 0.2

Tt

=0

w,

Wn

3.1

A Simple Assembling Task

Tb

Inserting a peg into a hole is considered as an example of a simple assembling task. Figure 1 shows the robot, the peg, and the hole as well as the frame of reference ~, and the task frame R t which is chosen such that its T axis has the same orientation as the hole. The task consists of inserting the peg into the hole while rotating it about the hole axis. To fulfil this task, 6 natural constraints have to be met which are imposed by the geometry of the mating parts, see Table I. Whereby I, T, v, and W stand for the force, the moment, the translational and rotational velocities correspondingly; the subscripts I, n, and b stand for the axis of the task frame; so It for example means the force along the I axis. The artificial constraints are interpreted as equations to be satisfied by the end effector, the control system has to enforce these equations as closely as possible. For this task, the robot is stationary at time 10 , It is desired to exert constant forces with the above given values (these values are selected only as an illustration, in fact, for such a task these forces are to be maintained zero to guarantee the right desired motion and orientation) and insert the peg 0.1 m with a 0.5 rad rotation within 10.0 seconds and ending with zero velocity. A parabola velocity function is selected to fulfil these requirements, the equations of these parabolas are listed in Table 1.

0

= 0

=0

Tn Wb

= 0.1011

= 0.0748 = 0.01 = 0.9972

equation where T points always to the tangent of the curve, N points to the normal, and B points to the binormal. The values of v" Wt , and Wb are calculated to fulfil equation (8) and maintain the desired orientation of the task frame. The values of In, Ib, and Tn are set as desired, here their values are selected only as an example.

4

Descriptor Model Space

III

the Joint

In the case of constrained motion due to assembling or machining tasks, for a complete description of the system dynamics the robot equations of motion in the joint space (1) and additional algebraic equations are needed. The latter ones describe the constraints as functions of the joint angles in the form

(q) = O.

(9)

The constraints also cause an additional term in the equations of motion containing the product of the Jacobian J ~ = 3(q)/3q and the contact force vector (Lagrangian multipliers). By introducing an augmented state vector, the so-called descriptor vector

IQ

Reduced Dynamic Model In the Task Space

5 (10)

Instead of writing explicit constraint equations, as it is necessary in the descriptor model method, one can consider the natural and artificial constraints, given for example in Table 1 and Table 2, to obtain a nonsingular dynamic model. In this case, an implicit or a verbal description of the constraints is sufficient for the derivation of the model. The derivation was presented in [7] and so only the results will be given here.

equations (1) and (9) form a system of differential-algebraic equations:

4.1

Constraint Equations for the Assembling Task

In this section a simple method to find the constraints equations of the first task is presented. The constraints c)( q) can be formulated in terms of the task frame, which describes the position and orientation of the end effector, and the initial task frame, which is stationary with its Tto-axis collinear to the centre line of the hole. These constraints simply read

Examining the constraints given in Table 1, it is obvious that the force is exerted along the N axis and the B axis and the moment about the N and the B axis; while the end effector moves and rotates along and around the Taxis. Hence, the force and the acceleration vectors can be written in one vector and equation (5) can then be written as

(12) (13) (14) (15)

N t · Tto = 0, Bt · T to =0, (Pt - Pto) · N to = 0 ,

(Pt - Pto)· Bto = 0,

Reduced Dynamic Model for the Assembling Task

5.1

A(q)Z

+ Cx(q,q)

=

r

( 19)

where A is the combined mass and Jacobian matrix given as A

where Pt denotes the position of the origin of the task frame R t in the frame of reference. The dot products of equations (12), (13) force the N t and the B t axis of the task-frame

[ MX123 T

_J T 4S6N

_JT 123N

_J T 45GB

_J T 123B

1

MXm T

(20)

to be perpendicular to the centre line of the hole. Thus, the rotational degrees of freedom around the B t - and the Nt-axis are constrained. Equations (14), (15) describe the translational constraints. Here the projections of the line segment from the origin of the initial task frame to the origin of the task frame on the N to - and on the Bto-axis are set to equal zero. These four equations, written in the task space, correspond to the four positional constraints given in Table 1. In order to obtain c)( q) , which is required by the model in (11), the forward kinematic matrix is employed to end up with explicit nonlinear constraints equations in the joint space.

where the subscript 123 means the first three columns of the subscripted matrix. The T, N, and B vectors in the above matrix are responsible for the transformation between the task and the reference frames. The vector Z contains the acceleration and force elements in the task space

Constraint Equations for the Machining Task

The assembling task of the previous section can be seen as an easy special case, where the translation and rotation take place along and around the stationary T axis. Here, we handle a more general case where the task frame is nonstationary and the path of the end effector is a general space curve with curvature and torsion components. The same reasoning and the same steps of obtaining the model in (19) are followed to obtain the reduced model for the machining task. Since the constraints are given in the task space while equation (5) is given in the frame of reference, a transformation is necessary. To switch from the velocity level of the constraints in Table 2 to the acceleration level of equation (5), the velocity vectors have to be differentiated once after the above mentioned transformation. This results in an additional normal acceleration term due to the curvature of the path. The reduced model is then given by

4.2

Here, it is important to note that matrix A is nonsingular.

5.2

Three equations describing the three given constraints on the motion are to be explicitly found . Beginning with equation (8), the first constraint is derived from the radius which is to be maintained constant, this means

(16) where P zt and P vt are the x and y components of the end effector position. The second equation is obtained by noting that the normal vector lies in the x - y plane, so

(17)

N, =0.

The third equation is necessary only if something like a groove in the surface exists, which constrains the motion, thus from (8)

P,t = 0.03 '" (arctan(PvdPx,)

+ n * 2TI),

(18)

Reduced Dynamic Model for the Machining Task

(22)

where

where

n = 0,1, ...

If this condition does not exist, then the end effector can move freely on the cylinderical surface.

11

A

T - J I23 N

[MX123T

-JTI23B

_J T 4S6N Mx,,,B] and

z = [v,

fn

fb

the BDF which are all supported by ANDECS [3]. Th e integration algorithm consists of solving the reduced model equat ion for the vector Z and the n integrating the acceleration. Two variations were experimented: in the first the state vector contains four/six (depending on the task) variables in the task space only while in the second one the state vector contains fourteen/fifteen ,·ariables in the reference and in the task spaces: namely, the position , the orientation, and the translational and rotational velocities.

Mx",T (23)

Wt

Tn

Wb

r

(24)

and k is the curvature of the path. Note that the angular acceleration has two components; the first corresponds to the curvature and the other to the torsion.

Numerical Integration

6 6.1

7

Numerical Integration of Descriptor Systems

A six degree of freedom PUMA 560 robot was considered for the simulation. The dynamic model of Armstrong [I] in the joint space was used without the centrifugal and Coriolis forces . In order to fulfil the simulation process of the different models, it was necessary to find symbolically the forward kinematic matrix, the inverse kinematic, the Jacobian matrix and its derivative; this was accomplished by using REDUCE [5] software package. To find the constraint equations and their derivatives, it was recommended to use a symbolic working software like MATHEMATICA [S]. To compare the two methods, the necessary input force vector r to generate the desired position and contact forces was calculated for the two tasks; to check the validity of the different models, r was calculated every time using all possible models (in the joint space and in the task space) and the same results were always obtained for the same task. Many simulation experiments were run varying the following

A parameter of significance for systems of differential-alge braic equations is the index. There are several definitions for this parameter [4]. For example, the differential index denotes the number of required differentiations of the system equations for obtaining differential equations for all elements of the descriptor vector (10). Thus, mechanical systems with holonomic constraints are of index 3. Another definition of the index appears in the field of numerics. In the analysis of roundoff errors the difference between the solutions of t.he exact system F(Y, Y) = 0 and a perturbed system

F(Y, Y) = 8(t) at the end of the integration interval can be estimated by an upper bound 11

Y(t) - Y(t)II S; C(II Y (O)

Y(O)II

+

+ oseSt max II.5(m +

max II.5(m-')(~),,),

ose9

(25)

• The state vector: the joint positions, velocities and the contact forces in the descriptor model; the position, orientation and the velocities of the end effector or the task frame variables in the reduced model.

where C depends only on F and the integration interval t. Here m denotes the perturbation index, which (for mechanical systems) equals the differential index. This index indicates, how many derivatives of the perturbation influence the solution for Y when initial-value-problem (IVP) solvers for ordinary differential equations are applied to this class of systems. Thus, special IVP solvers for index-3-systems which stabilize the solution, for example RADAU5 [4] or ODASSL [2], must be applied. We have chosen the latter one because of its implementation in the software package ANDECS ' [3]. In ODASSL, the stabilization of the solution is achieved by solving the overdetermined differentialalgebraic system of equations (ODAE)

=

[M-

1

(q)[r+;Tx-C(q,4)] ~(q)

• Integration code: Runge-Kutta, Adams - Bashforth - Moulton and BDF method for the reduced model; ODASSL for the descriptor model. • The sampling time: 0.1, 0.01 and 0.001 second. • The error tolerances: low and high. Because of the shortage in the space, only the simula· tion results of task 1 will be shown. It is important to note that no feedback was provided to handle any integration errors or deviations. Figure 3 shows the desired rotational and translational velocities of the end effector (w, and v,) around and along the tangent and the desired normal contact force; these variables were selected to be shown since they represent the translational and rotational motions as well as the contact forces. The simulation parameters for the two methods are given in Table 3.

(26)

~(q)

by the Gauss-Newton-iteration. This ODAE contains not only the index-3-constraints on position level, but also the index-2-constraints on the velocity level, which are obtained by differentiating the index-3-constraints with respect to time.

6.2

Simulation and Results

Table 3: The Simulation Parameters

Numerical Integration of Reduced Dynamic Models

Since the derived models in (19) and (22) are non singular, it i, possible to use any available ordinary differential equation solver like Runge-Kutta, Adams-Bashforth-Moulton or

Descriptor Model

Reduced !>.!odel

Integration method

ODASSL

LSODE2

Sampling time [sec]

0.01

0.01

Relative tolerance

1.0E - 08 q: 1.0E - 10

1.0E - 06

Absolute tolerance

ANDECS (ANalysis and DEsign of Controlled Systems) is a registered trademark of DLR. I

12

q: 1.0E - 07

1.0E - 08

CPU time [sec]

380.0

11.1

No. of function calls

2732

152

No. of Jacohian ca\ls

1276

15

The differences between the obtained variables after the simulation and the desired ones are shown in Figures 4, 5, and 6 where the curves with the stars correspond to the reduced dynamic model and the solid curves correspond to the descriptor model. The simulation parameters of the reduced model were chosen to become the best possible results. For the descriptor model, due to the extensive numerical calculations required for the Gauss-Newton-iteration in ODASSL the best reachable accuracy is restricted; the corrector step did not con\'erge with absolute error tolerances less than 10- 7 for the \'elocities. Note also that, due to the above reason, the required CPU time is clearly longer than that required for the reduced model. The constrained equations for this task signifies that the T axis of the end effector remains constant. The simulation showed that this was the case in the descriptor model where the error was exactly zero and in the reduced model where the error was of order 10- 9 . It is senseless to show any figures in this regard.

8

[3] Gaus, N., and M. Otter (1990). Simulation of dynamic systems with ANDECS-module DSSIM. User's Guide and Reference Manual. ,B-Release of Version 2.0. TR R27-90, DLR, Institut fiir Dynamik der Flugsysteme, D-8031 Oberpfaffenhofen. [4] Hairer, E., Ch . Lubich, and M. Roche (1989). The numerical solution of differential-algebraic systems by Runge- Kutta methods. Lecture Notes in Mathematics 1409, Springer- Veriag, Berlin. [5] Hearn, A. (1987). Reduce user's manual, Version 3.3. Rand Publication CP78. [6] Tahboub, K.A., and P.C. Miiller (1990). A unified modeling approach for assembly and machining robot control tasks. Submitted for publication in IEEE Trans. Robotics and Automation. [7] Tahboub, K.A., and P.C. Miiller (1991). A reduced dynamic model for constrained robots in the task frame. Robotersysteme, Vol. 7, pp. 45 -52.

Conclusions

The two methods were presented and illustrated through two examples, all the necessary steps to realize them were carried out. The differences between the two methods are basic. The descriptor method employs the robot model unchanged and requires in addition algebraic equations to describe the geometric constraints, these equations have to be found for every task explicitly and symbolically; the contact forces are calculated as a result of satisfying the constraint equations. In the reduced dynamic method, the constraints are observed in the structure of the model, they do not need to be explicitly written and no equations need to be solved to obtain the reduced model. It suffices to know in which direction can the end effector move, as a result the contact forces have a more active meaning as those in the descriptor model. One is limited to use the few available DAE codes and has to prepare the necessary equations in the descriptor model. The simulation results show that the accuracy in maintainig the constraints was high. Even in the presence of modeling errors the constraints are observed, so the descriptor model is more suitable for assembling tasks where exact contact forces are not needed. For the reduced dynamic model , one has the advantage of using well established and developed integration codes to end up with more accurate results. Moreover, it is possible to cover a large class of tasks without changing the structure of the model. The aims of the work were reached; the two benchmark tasks were developed and realized, it was also shown how it is possible to use the descriptor method to handle constrained robots where the constraints were not only on the position but also on the orientation. The strong and weak points of the two considered methods were also discussed.

[8] Wolfram, S. (1988). Mathematica, a system for doing mathematics by computer. Addison- Wesley, Redwood City.

Robot

Ra

Figure 1: A peg into a hole task

Cylinder

ilelix

References

--- ---

[1] Armstrong, B., O. Khatib, and J. Burdick (1986). The

Robot

explicit dynamic model and inertial parameters of the Pl'~!A 560 arm. IEEE Int. Conf. Robotics and Automation , pp. 510-518. [2] Fuhrer, C. (1988). Differential-algebraische Gleichungssytcmc in mechanischcn Mehrkorpcrsystcmen. Theorie, numerische Ansiitze und Anwendungen. Dissertation TU Munchen.

Figure 2: A cylinder machining task

13

Assembling task

Assembling task

l'

E-1 1.0

E-5 O.SO 0.25

0.8 06

'" 0.00

1:

~-0 .25

0.4 0.2

~ -O .SO

0.0

;;; -0 .75

.,>

.~

.. -1.00

-0.2

u; c:

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 time ! s

~

E1

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 time!s

-7

- - desired translational velocity desired rotational veloc~y rm I f r

- - error of the descriptor model - - error of the reduced model

Figure 3: The desired translational and rotational velocities of the end effector along and around the tangent , and the desired normal contact force

Figure 5: The error in the translational velocity

Assembling task

t

Assembling task

t

E-4 0.0

z

-0.8

.,

E-3 0.2 0.0

~

.Q

-0.2

!! -2.4 3.. -3.2

;;; -0.4 §

-4.0

-0.8

-4.8

-1.0

iij

c:

.2 ftj

!!

~

Ql

-1.6

'"" ~>

E1 -7

0

c:

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9 time!s

E1

-0.6

0.0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.9

-7

time!s

- - error of the descriptor model - - error of the reduced model

E1 -7

- - error of the descriptor model - - error of the reduced model

Figure 4: The error in the rotational velocity

Figure 6: The error in the contact normal force

14