- Email: [email protected]
Identification and adaptive control of robots
Mechatronics.Vol. 3, No. 2, pp. 149-166, 1993 Printed in Great Britain
0957-4158/93 $6.00+0.00 Pergamon Press Ltd
IDENTIFICATION AND ADAPTIVE CONTROL OF ROBOTS B. LAy'ros Technical Universityof Budapest. Departmentof ProcessControl H-1111 Budapest,Mfiegyeternrkp. 9., Hungary (Received 30 June 1992)
Abstract - The paper deals with the model identification of industrial robots as the basis for control algorithm design. The first part discusses the classical decentralized joint control using internal current and velocity loops and external position control, where the plant to be identified is the black-box containing the intero.,tl loops and an integrator. The paper presents the identification of a SCARA robot using LS, GLS and IV methods and real closed loop data under normal working conditions. The second part discusses the convergence test of self'tuning adaptive control for SCARA and PUMA robots using the advanced robot control algorithm of Slotine and Li to find the unknown parameters of the robot dynamic model. The convergence has been tested by using a transputer based simulation system and the symbolic form of the robot dynamic model.
1. INTRODUCTION Robot control systems have usually a hierarchical structure and can be divided into specification level (robot programming language, graphic simulation, path design), algorithmic level and torque control (power electronic) level. From the theoretical point of view it seems to be acceptable to divide the control algorithms into decentralized joint control and mole advanced robot control algorithms. Decentralized joint control attempts to control every joint axis separately without considering the coupling effects between the different joint axes. To this class belong the three loop current, velocity and position contro!, the PID control of the position with motor voltage controller output, etc. To the class of advanced control algorithms belong such methods like computed torque techniques, resolved motion acceleration control, hybrid position and force control, model reference adaptive control, selftuning adaptive control, fuzzy control, etc. [1, 2, 3]. Characteristic for this second class of control algorithms is that they are based on the nonlinear dynamic model of the robot. The paper discusses parameter estimation methods to find the average linear model for decentralized joint control, and selftuning adaptive control to find the nonlinear dynamic model during control.
2. IDENTIFICATION OF THE LINEARIZED MODEL FOR DECENTRALIZED POSITION CONTROL
Industrial robot control systems often apply internal current (torque) control, central velocity (revolution) control and outer position (joint variable) control for the robot axes. It is assumed that the actuators are DC motors with gear reductions.The controllers for the current and velocity loops are analog PI con149
150
B. LANTOS
trollers, while the position controllers are DDC controllers containing PID or PI controllers and a feedforward path for the desired velocity. If the mechanical and electrical parameters of the motors and the gear reductions are known then the following method can be suggested for the design of the decentralized joint controllers. For the current control loop we can consider the backward electromotor force as a disturbance, and the design problem reduces to the design of a PI controller for a first order system with the electrical time constant. The integrator assures that the gravity load is compensated in steady state. After having designed the current controller, the backward electromotor force can be taken into consideration using root locus method, and the dominant low frequency pole can be determined for the velocity control design, whose place is very critical for the quality of the whole control system. This dominant pole can be cancelled by the zero of the PI controller of the velocity loop, and the controller gain can be chosen to result nearly critically damped closed loop. The resulting system after velocity control can be approximated by a time invariant third order linear system and an extra integrator. This is only an average model, because a nonlinear system linearized along a time varying path yields always a time variant linear system. An excellent treatment of identification methods can be found in [4]. The identification problem in discrete time is now the following. Let z -~ denote the backward shift operator. The linear, time invariant approximation of the black-box containing the current and velocity control loops of a robot axis can be assumed in the form of a discrete time SISO (single input x, single output y) system
B(z-~)
y(t) = A - - ~
z-a x(t),
alz-l+...Wam Z-m,
:
1
=
b~z-l+...+b,,z -m,
+
y(,) = -aly( t - 1 ) . . . . . amY( t
(1) - -
m).4_
+ blx(t - d - 1)+'"+bmx(t - d - m), where r is the sampling time, d . r is the time delay (if any) and t is an integer multiple of the sampling time. Let p and hr (t) denote the system parameter and the observation respectively, where p = (a,, o~ ..... a,,, b~, b~ ..... bin) ~, hT(t) = ( - y ( t - 1) ..... - y ( t - m), x ( t - d - 1) ..... x ( t - d - m)), then
y(t) = hr (t)p.
(2)
In the case of the least-squares method (LS) the system parameter p is the analytical solution of the optimum problem N
2
V ( p ) : ,=,,~d+l { Y ( t ) - h r ( t ) P } - - 4
min,
(3)
where V(p) is the sum of the squares of the equation errors, and N is the number of observations. Using the notations
Identification and adaptive control o f robots Y=(y(m+d+l)
15 i
. . . . . y ( N ) ) T,
,,__ Ibm(m+ =
- y ( m + d)
....
!
:
-y(N)
....
l
y ( d + 1)
x(,,,)
...
x(1)
x(U-d-1)
...
x(U-d-m)
]
i
y ( N - m)
!
!
1'
the LS optimum can be written in the form /3 = (H/-H) -1HTY.
(4)
If the system output contains an additive noise v(t), then y(t) = h ~ ( t ) p + v(t). Let us introduce the notations |
N
q)~-N-d-m
y~ h ( t ) h r ( t ) = 1 HrH, t=m+d+l N - d - m 1
qghv - U _ d _
N
m
y~ h ( t ) v ( t ) t=m+d+l
N
-
1 HTv, d - m
then ¢p~ and ¢Phv are convergent to the mean values ¢Phh ---) E { h ( ' ) hT (')} and ¢Phv ~ E { h ( t ) v ( t ) } respectively, and p satisfies b = P + ~0~,]¢h~•
(5)
Hence/3 --~ p asymptotical if h(t) and v(t) are uncorrelated. This is always satisfied if v(t) is white noise. Otherwise v(t) may be a coloured noise, and the estimate/3 is not necessarily convergent to the system parameter p. The error is dependent on the signal to noise ratio. Since the identification of the robot is based on closed loop data, we can expect coloured noise and biased estimation. Hence the application of instrumental variables (IV) or generalized least squares (GLS) can be suggested. The method of instrumental variables substitutes the LS estimation by the relation /3 = ( Z z H ) - I Z r Y .
(6)
The instrumental variable z(t) should satisfy that ¢Pzh is invertable and ~Pzvhas zero mean value, where 1
¢Pzh - N - d - m _
qg~v
1
N- d - m
N
y~ z ( t ) h r ( t ) = 1 ZrH, t=m+d+l N - d - m N
• z(t)v(t) = 1 ZTV. t=m+d+l N - d - m
152
B. LANTOS
The instrumental variables can be generated by the auxiliary model
f~= ~z Z--
x(t),
z~(m + d + i
1)] ]=
(7)
Z~"(N) + d) =
....
:
:
-~(N-1)
...
+ 1)
41)
]
:
....
fJ(N-m) x(N-d-1) :
: x(n-d-m) :
...
J"
The auxiliary model can be initialled by using LS estimation. However the model is not exact, thus and v(t ) and similarly y(t) and v(t) are correlated. Hence a sequence of estimations should be applied:
z(t ) = h(t, ~, x)
1. LS estimation---) A,/}.
2. Generating the instrumental variables:/]33(t) =
Bx(t ).
3. IV estimation---> ,~,/~.
4. Jump to 2. if the parameters are not convergent, otherwise stop. The method of the generalized least squares assumes coloured noise v(t) produced from the white noise e(t) by the filter C(z -1 ). The system model can be assumed in the form
a(z-1)y(t) = B(z-l)z-dx(t)+C(z-~)e(t) a(z-~)~(t) = B(z-1)z-aYc(t)+e(t), ~(t) = c-l(z-1)y(t), Yc(t)=C-l(z-1)x(l).
or (8)
Since C(z -1) is usually unknown, thus C(z-1), a(z -1) and B(z -1) have to be identified simultaneously. The algorithm has been implemented assuming C(z -1 ) = 1 / F ( z -l): 1. Let/~(z-1) = 1. 2. Determine ~(t) and ~(t) using
F(z-' ): ~:(t) =/~(z-' )y(t)
and Yc(t) = l:(z-'
)x(t).
3. Determine the parameters o f / t ( z - ' ) and/~(z-') from :~(,) and 5'(t) by using LS estimation.
A(z-1)y(t)-B(z-~)x(t-d). 5. Apply LS estimation for the model F(z -1)v(t) = e(t), ( a : = F, B:= 0, y(t): = v(t), and no x(t) is present),
4. Determine v(t) =
where e(t) is a white noise to find the new value of/6.
6. Jump to 2. if the parameters are not convergent otherwise stop.
Identification and adaptive control of robots
153
During the structure estimation the purpose is to find the system order (m) and the time delay (d). In the case of LS estimation the loss function can be used for this aim. For increasing system order LS estimations are performed for the same observations {xn, Yn}, and the equation errors normalized to one step (the residuals) are determined for every system order. That system order has to be found, for which the residual can not be essentially improved by increasing the system order. In the case of IV or GLS estimation the whiteness of the residuals is necessary to the unbiased estimation, that can be tested using the autocorrelation function. It may be another practical method to compare the sum of the squares of the differences between measured system output and computed model output for different model orders. The efficiency of the identification of the linear approximation will be demonstrated for a SCARA robot. The robot repeats a composite motion cycle, so that it is useful to design new controller parameters for the position controllers of the different axes. The base for the controller design may be the approximated linear model of the resulting subsystem containing the internal velocity and current control loops of the robot axis. The subsystem begins at the D/A converter (DAin) after the position controllers and is finished after the encoder and counter converting the joint position into large resolution integer values (00 or qhO). The robot has been sent through the path formulated by the robot programming language and using the actual (not necessarily suitable) controller parameters. During this normal work 4000 closed loop data x(t) := Dain(t) and y(t):= 0(t) have been collected with 1 ms sampling time and stored for later identification. The identification has been performed in batch mode in a host computer using LS, IV and GLS methods. The measured output of the subsystem (qhO)and the answer of the identified models (qhidO)have been computed and compared for different system orders and methods, and the best result has been accepted. The best approximation has been achieved for GLS method and order 4 without time delay [5]. For the second robot axis Fig. 1 shows the subsystem input, the subsystem output and the answer of the identified model for the same input. It can be seen that the approximation is excellent, although the energy content of the input signal is not high. The discrete time model of the subsystem and the equivalent continuous time model for step inputs (for zero order hold) are respectively
Gtz_l~ = G(s) =
--0.007195z -l + 0. 005618z -2 + 0.004719z -3 + 0.008259z -4 1 - 0. 828761z -I - 0. 281923z -2 - 0. 481865z -3 + 0. 592652z -4 ' 24.3666 - 0.049844s + 1. 24168.10 4 s 2 - 1.23455.10 -8 s 3
(9) (10)
s(I+Tls)(l+2VTzs+Tz2s 2) T1 = 0.00528626 sec,
D = 0.0756408,
Tz = 0.000052803 sec.
The continuous time subsystem is of non-minimum phase (there is a zero at s = 469), and has an approximately cancelled conjugate complex pole (-167 + j2202) and zero (-18 4 + j2043) pair in the high frequency domain. To avoid numerical instability caused by the large amount of data, the Householder method has been applied for matrix inversion. Moreover, it was necessary to test the poles of the discrete time model in different steps of the algorithm, and correct the pole in the neighborhood of z = 1 if it became numerically unstable. A built in program has supported the controller design [5]. The position controller design can be elaborated immediately in discrete time using Tustin (w) transformation, root locus method or state feedback. In this case, it is less critical if the discrete time model contains a negative real pole in z with an odd number of multiplicity. However the preferred control algorithms for position control were PIF or PIDF control, where F denotes the feedforward path from the desired path velocity 0a to the set point input of the velocity control. The controller design has been performed in continuous time. Since in this case the open loop contains
154
B. LANTOS DAin [2] qhO [2] qhldO [2]
0 "A C]
2~
20
.;J C 0 tO O.
,,~
t5
io
(.. t,-
o W e.-
-5
e-I
0
5
~0
15 t
20 E5 [ E -:L see I
30
35
40
Fig. 1. ldentifiLation of the hm~arized model
two integrators, hence one zero of the controller should be placed in the low frequency domain and a large phase margin has to be designed. The continuous time PID controller has been implemented in discrete time using its step equivalent form. Although the identification and controller design has been performed in batch mode on a host computer, the whole process can be implemented in real time on a multiprocessor control architecture.
3. SELFTUNING ADAPTIVE CONTROL The base for advanced robot control algorithms is the robot dynamic model, which gives the relation between driving torque (z), inertia matrix (H) and the centripetal, Coriolis and friction effects (CO), and the gravitational part (D):
H(q)it + h(q, (t) = H(q)~ + C(q, dl){l + D(q)
= z.
(11)
Using the mass mi, the center of mass Pci and inertia moment Kei belonging to a coordinate system having origin in the center of mass but axes parallel to the Denavit-Hartenberg frame Ki, the quantities in the dynamic model can be expressed using recurrence or direct formulas [2, 6]. Similarly, H and h can
Identification and adaptive control of robots
155
be determined also by using the inertia moment Ki expressed in the Denavit-Hartenberg frame Ki. The later case is preferable because the motion equation becomes linear in mi, m~Pci and Ki. The main problem of the realization of advanced robot control algorithms is caused by the fact that the driving torque should be computed by the controller within 1-10 ms. Hence the controller should be realized by using a multiprocessor architecture containing fast arithmetical processors which allow parallel computation. Preferable to parallel computation are nonrecursive algorithms for the computation of the dynamic model H(q)~ + C(q, Cl)Z + D(q) = "t"and Y(q, q, ~, ~)a = r so that the time consuming computations can be distributed among the processors [2]. A special list processing program has been developed in Pascal for IBM-PC/AT which determines H, h and o~ in symbolic form [8]. If the robot dynamic model is known in symbolic form where the inertia parameters (mass, center of mass, inertia moments) are also contained in symbolic form, then independent parameters ot1, ot2 .... can be chosen so that the dynamic model is linear in the parameters:
H(q)~ + C(q, (I)~ + O(q) = Y(q, q, ~, ~)a = "r.
(12)
The identification problem is the determination of the actual parameter value a (Slotine and Li [9]). The method can be generalized to Cartesian and hybrid control (Slotine and Li [10]). Let q,(t) be the reference signal defined later and s = 0r - q. The control law is chosen as "t:=/~/(q)qr + C(q, (l)¢r
+
D(q) + Kt,(q r - q) + KoS,
(13)
where/q, C, b is computed at the estimated & and measured q, q. The Lyapunov-like function
v=lsTH(q)s
1 + ~(qr - q)TKe(q, - q) + I ( ~ - a)T/"(t~ -- a )
(14)
Can be USed, where K e, K D, e are positive definite symmetric matrices. If a does not change during the adaptation process then the adaptation law
d_.._~= F_lyT (q, (l, glr, qr) s dt
(15)
can be chosen and
(/= -sT Kos <- 0 is satisfied. If the reference signal is t
qr := qa + A S(qa - q ) dr ~ s = ((la -dl)+ A(qa - q )
(16)
0
and the roots of det (sl + A) = 0 lay in the left half of the complex plane then the trajectory is convergent to the sliding surface s = 0 and the error qa - q --~ 0 exponentially on them. Hence the control (13) with the parameter estimation (15) is globally asymptotically stable. Since the mass, center of mass and inertia matrix of the robot links are usually unknown, the convergence test of the selftuning adaptive control method has been performed on simulated robot models, allowing the comparison of known model parameters and their estimations. Multitransputer based simulation system has been developed to examine a broad class of control algorithms containing usual industrial robot control (position, velocity and current control), computed torque technics and decentralized PID controllers, resolved motion acceleration control, hybrid position and force control, selftuning adaptive control and fuzzy control [3]. It is assumed that the control architecture contains low level torque (current) control and the actuators are DC motors with gear reductions.
156
B. LANTOS
The controllers for the current loops are analog PI controllers. The mass and inertia parameters of the link and the actuator (or load) mounted on the link are reduced to resulting parameters, and the effect of the kinetical energy of the rotor is approximated by ~ t g , i ( v~ is the gear reduction and Ori is the rotor inertia). The friction of the rotor and the gear reduction is approximated by v2f~i + fg,ar,,. Let L,i, Rri, qi and c2i denote the rotor inductivity, rotor resistance, electromotor force constant and torque constant respectively. The coupled mechanical and electrical equations have the form EDij(q)dlj + EEDijk(q)(lj~lk + D i ( q ) + ~iit)riqi +(V2fri + fgear,i)(]i = ViC2iiri, j j k
diri
Lri ~
~ .
+ Krilri = Ui -- Cli Vi~ti '
(17) (18)
where i~ is the rotor current and u~ is the motor voltage. The PI controllers are realized by operational amplifiers, therefore the set-points [a are voltages. The simulation is based on the robot dynamic model in symbolic form. Using the substitutions Hii = Dii := Dii + v26)ri, hi "=/1/ +(V2fri + fsear.i){ti, and the state variables q, q, i~, 1 the system can be written in the form of dynamic and static equations:
dq_
dt-q' d___~q= H_1{coll (vic2,i, i _/~)}, dt di__L~ = co. { - c. v,o, - R.,i.,)}. dt
(19)
d ! = 7~ - T = i~ - coil ( r, ir, ), dt u:coll
apii ]
1
The torque output r of the high level controller is computed in the sampling instants and converted to ~'e.After D/A conversion their components are sent to the reference signal inputs (set point functions) of the current loops. Further state equations have to be added to the system if the sensors contain time constants. Saturations in the D/A converters, the operational amplifiers, the motors and the sensors can be taken into consideration during simulation.
3.1. Convergence Test for SCARA Robot The dynamic model of the 4dof RRTR SCARA robot has been generated using symbolic computation. The functions Dij, Dok and D i of the robot model are sums of terms. Every term is the product of a constant, an expression of the parameters and an elementary function (a product of powers of Si = sin qi, 6",. = cos qi or qi, depending on the joint type). The parameters /'1 . . . . . Pn of the links should be corrected by the rotor inertia, gear reduction and friction:
P1 = m2 a2 + P2xzm2 + p2xlml + p2y2m2 + rn3a2 + 2pcxlrnlal + + 2pcxzm2a2 + mla~ + P~xsm3 + m4a 2 + p2y3m3 + Kczl + 2p~x3m3a2 + + gcz 2 + p2x4m 4 + p2y4m 4 + p2ylrrl1 + gcz 3 + rn:2a~ + gcz 4 + + rn3a2 + m4 a2 + Orl V2
Identification and adaptive control o f robots
157
P2 = Pcy4m4a2 P3 = ,Dcx4m4a2
P4 = P5 P6 P7 P8
(P~2rr~+ Pcy3ma)a,
= (Pcy2ra2 + m2a2 + Pcx3m3 +m3a2 + m4a2)al = Pcy4m4al = Pcx4m4al = p2x2m2 + p23m3 + 2Pcx2rn2a2 + 2Pcx3rn3a2 + p2y2m2 + P2yam3 +
+ m2a~ + p2x4m4 + m3a~ + p2y4m4 + Kcz2 + Kcz3 + Kcz4 + rl'14a~ P9 = p2x4m4 nt"p2cy4rl'14 -I- gcz 4
t'1o
(20)
= m3 + m4 + e,3 ~
Pll = (m4 + m3)gz
P12 = o,s P13 =
(~r4 V~4
P14 = e15 =
f2
P16 = f3
/'17
=
A
Appendix 1 contains the dynamic model of the SCARA robot. The simulation system makes it possible to distinguish the dynamic parameter values of the robot (for computations of the sensor signals using Runge-Kutta method; parameters P~) from the dynamic parameters of the selftuning adaptive controllers (the estimated parameters Bi). Commands allow to modify Ke, Ko, A and F for the adaptive controller during the identification. The convergence test has been performed for initial parameters Bi (0) scaled up by 500%. Because the robot is of open chain type, the identification has been elaborted in a backward order of the links. The influence of the parameters onto the driving torques of the links has been formulated in a parameter graph supporting the tuning strategy. Since the command set of the simulation system allows to switch off the tuning for some parameters in different (acceleration, constant velocity, steady state) path sections, it was possible to find the best tuning strategy by using the relations in the parameter graph. It is a real problem of error based adaptation, that the control error becomes small if the parameter adaptation converges, but the energy of the error signal is not enough to continue the tuning. Therefore, in every stage of the adaptation, the parameters have been divided into two groups: parameters for which the tuning has already been finished, and parameters to be tuned simultaneously. If the tuning for some parameters of the second group was finished, then the user could decide whether the actual, the initial or a previously stored value has to be used for a remaining parameter in the second group to continue the identification. The path, control error, set point for current control and parameter tuning are shown in Figs 2-5 for 5 ms sampling rate. The real (P~), the estimated (B,) parameters and the relative accuracy (Pi / Bi in %) are shown in Table 1. It can be seen that the accuracy of all the dominant parameters is satisfactory, the control error is small, and there is no saturation in the control system.
158
B. L A N T O S .do t~] qd0 {2] Qd0 tS] ado [4]
~ O ~ , A ~ O ~ X
7
8 B 4
3
|
o t
[ see ]
Fig. 2. Path signals/or the S C ~ A robot
O
dq [i] d. [a] dtl [a] a . t4]
~
x
[]
L
2.5 0.0
~ -2.5 ~-5.0 -7.5
/A
-~,-iO. 0
/
c~-t2.5 i 0
[ ne ] Fig. 3. Error signals of the SCARA robot
Identification and adaptive control o f robots Idh 11] lab [2] Idh[3] im~t4]
159
~ 0 ~ A ~ 0 ~ X
30 25 20
o
~5 -1o .o
-15
/o
-~0
o
1 t
[ see ]
Fig. 4. Adaptive controller outputs of the SCARA robot
st11] A [14] s [15] B [17]
\x -5
-1o
-20 -25 -30
0 t
[ see |
Fig. 5a. Parameter tuning of the SCARA robot
160
B. LANTOS n
El( t1 B[ 4)
~
tat 51
O (
"
]
s t el S tt2] e (~3] 22.5 20.0 :17.5 J.5.0 t2.5 t0.0 7.5 5.0 ~.5
/
0.0 t
[ aec )
Fig 5b Parameter tuning of the SCARA robot
e ['~6]
tB
t6
to
\ o
~
~
~
.~
5"
"~ t
~
~
~
[ aec ]
Fig. 5c. Parameter tuning of the SCARA robot
16
ti
~
Identification and adaptive control o f robots
161
TABLE 1 Parameter Convergence for S C A R A Robot real
estimated
P1 = 4.230 P2 = 0 P3 = 0 P4 = - 0.162 P5 = 0.947 P6 = 0 P7 = 0 Ps = 0.727 P9 = 0.0040 P10 = 130.252 P11 = - 6 7 . 1 9 9 P12 = 0.419 PI3 = 0.405 P14 = 14.500 P15 = 13.800 P16 = 3948.866 P17 = 13.400
B 1 = 4.494 B2=0 B3 = 0 B 4 = - 0.515 B 5 = 0.866 B6 = 0 B7=0 Bs = 0.765 B 9 = 0.0044 Blo = 131.577 Bn = - 6 7 . 1 6 5 B12 = 0.746 B13 = 0.522 Bla = 14.488 B15 = 13.788 B16 = 4082.003 BI7 = 13.402
accuracy 94.1%
31.5% 109.4%
95.0% 91.0% 99.0% 100.1% 56.1% 77.6% 100.1% 100.1% 96.7% 100.0%
3.2. Convergence Test for 3dof PUMA Arm The dynamic model of the first 3 degree of freedom of the PUMA 560 robot has been generated using symbolic computation [8]. The parameters P1. . . . . P21 of the links have been corrected by the rotor inertia, gear reduction and friction:
P1 = 2pcy3m3d3 + Kyl + rrt3d2 + Orl v? P2 = m3a~ + Ky2 + 2pcx2m2a~2 + m2a~ P3 = 2pcx3m3a3 + m3a~ + Kz3 P4 = Pcz3m3a3 P5 = Kx2 P6 =
K~3
P7 = Pcy2m2a2 Ps = (Pox3 + as)m3a2
P9 = Pcz3m3a2 r,o = (pcy3a3 + Pcx3d3 + a3a3 )m3 PH = (Pcy3a3 + m3d3 +Pcz2a2)a2 P12 =
Pcz3m3d3
P13 = m3 a2 + 2pcx3m3a3 + ,ocx2m2a2 + K~2 + Ky3 + m2a 2 + m3a 2 + 6)r2~ P14 = 2pcx3msa3 + Ky3 + rosa 2
P15 = 2#cx3m3a3 + m3a2 + Kz3 - Kx3
(21)
162
B. L A N T O S /°16 -- m34 + rr~a2 - Kx2 + 2pcx2m2a2 + Ky2 Ply = - ( Pcz3m3a3 - ,Ocyzm2a2) /°18 = (a3 +,Ocx3)m3gz P19 = (m3a2 + Pcx2m2 + m2a2)g~ P20 = ,Ocy2m2gz P21 = Pczam3gz
P22 = or3 v~ P23 = fl
P24 = f2 P25 = f3 The parameters P1. . . . . P25 have been computed from the geometrical, inertial and motor parameters of a real P U M A 560 robot (Tam et al. [11 ]). The convergence test has been performed for initial parameters B i (0) scaled up by 500%, and in backward order of the links. A suitable tuning strategy has been found by using the parameter graph of the 3dof P U M A ann. The convergence of the parameter tuning is shown in Table 2 for 1 ms sampling rate. The accuracy of the estimated dominant parameters is satisfactory [12]. It has been found, using simulation support to find suitable controller parameters, that the control error is also small without any saturation in the system. TABLE 2 Parameter Convergence for P U M A Robot real P1 = 1.151 P2 = 3.434 P3 = 0.014 P4 = - 0.020 P5 = 1.333 P6 = 0.546 P7 = 0.048 P8 = - 0.017 P9 = 0.458 P10 = - 0.007 P n = 2.449 P~2 = 0.160 Pi3 = 3.333 P14 = 0.547 P15 = - 0.532 P16=2.101 PI7 = 0.069 P18 = 0.383 P19 = - 53.115 P20 = - t.119 P21 = - 10.608 P2z = 0.583 P23 = 14.500 P24 = 14.300 P2s = 19.300
estimated B 1 = 1.092 B 2 = 3.278 B 3 = 0.016 B4 = - 0.029 B 5 = 1.646 B 6 = 0.523 B 7 = 0.027 B8 = - 0.033 B9 = 0.418 B10 = - 0.016 Bll = 2.393 B12 = 0.196 B13 = 3.160 BI4 = 0.497 B15 = - 0.693 n16 = 2.618 Bl7 = 0.043 n18 = 0.387 n19 = - 5 3 . 1 1 3 B2o = - 1 . t 1 8 B21 = - 10.607 822 = 0.569 B23 = 14.178 B24 = 14.295 B25 = 18.794
accuracy 105.4% t04.7% 83.9% 69.1% 81.0% 104.3% 175.7% 50.2% 109.5% 44.0% 102.4% 81.4% 105.5% 109.9% 76.8% 80.3% 159.7% 99.2% 100.0% 100.1% 100.0% 102.4% 102.2% 100.0% 102.7%
Identification and adaptive control of robots
163
For the 6dof case a multistep identification process can be suggested. The first step may be the identification of the parameters of the first 3dof of the robot for fixed q4, q5, q6. The parameters will not be unique because of the dependence on the fixed q4, qs, q6 configuration. The second step is the identification of the parameters of the 6dof robot using the previously identified dominant parameters for initial values. The first and second step can be performed in a laboratory. The third step should be repeated during the robot application if the influence of an unknown load (mass, center of mass and inertia moments) has to be identified. In this third step the identification of the robot parameters is switched off and only the identifiable load parameters are tuned. Appendix 2 summarizes some methods for model decomposition. The minimum dynamics parameters can be found by using the method of Kawasaki and Kanzaki [ 13], if the robot is described in frames based on the modified Denavit-Hartenberg notation. Determining the independent inertial parameters combined with model reduction can be performed by using the method of Sheu and Walker [14]. 4. CONCLUSIONS Model identification and selftuning adaptive control of industrial robots have been discussed in the paper. The identification of the average linearized model of the subsystem containing the internal current and velocity control loops has been presented using LS, IV and GLS parameter estimation methods for a SCARA robot based on real closed loop data under normal working conditions.The model can be applied to the design of the position controller in decentralized joint control systems. The selftuning adaptive control algorithm of Slotine and Li [9] has been applied with some practical modifications to the SCARA robot and the first 3 degree of freedom of the PUMA 560 arm. Simulation results have experimentally proved the convergence of parameter identification. The control errors were satisfactorily small during the adaptation process. Multitransputer based simulation system has supported the design of the controller parameters avoiding saturation. The simulation system turned out to be an easy tool to find the favourable tuning strategy for real-time implementation in the future. Acknowledgement - Support for the research of the model identification to decentralized joint control was provided by MIKROMATIKA Ltd. Budapest under grant No. BME--421.022/89. Support for the research of selftuning adaptive control is provided by the Hungarian National Research Program under grant No. OTKA-756. Appendix 1. Dynamic model of the SCARA robot Using the elementary functions El = 1, F 2 = $ 4 , F 3 = C 4, F 4 = S 2, F 5 = C 2, F 6 = $24, F 7 = C24, the symmetry properties of Di/ and Dijt, and the robot parameters '°1. . . . . P17, the nonzero elements in the robot dynamic model can be written in the following form: D n = F~P1 - 2F2P2 + 2F3P3 - 2F4P4 + 2FsP5 - 2F6P6 + 2FTP7 i)12 = F1Ps - F4P4 - 2F2P2 + 2F3P3 + FsPs - F6P6 + F7P7 014 = -&1"2 + Fa1"9 + F31"3 - F6e6 + FTP7 Dz2 = F1P8 - 2 F 2 Pz + 2F3P3 /)24 = - F 2 P 2 + E1P9 + F3P3 033 = ElelO 044 = F~e9 0112 = -0211 = 0122 = -F41"5 - F61"7 - F71"6 - &1"4 Oll 4 = -O411 ~-- O124 -~. O144 = - F 2 e 3 - F3P2 - f 6 e 7 - f T e 6 O214 = -O412 = 0224 = - 0 4 2 2 = & 4 4 = -&1"3 - 5 1 ' 2 03 = - E l e n
164
B. LANTOS El,1 = qrl El,2 = - F 2 ( 2 ( q r l
q-qr2)+Or4)-f3(t~lOr4
"t-q40rl "1"q20r4 + q40r2 "b q40r4)
Y1,3 = F3 (2(qrl + qr2 ) + qr4 ) -- F2 (qlOr4 + O,qrl 4- q20r4 + 04qr2 + 04qr, )
YL4 = - & (20,, + 0,2) - Fs(q,G2 + 020,~ + 0242) r~.5 = 6 ( 2 0 , , +O,2)-G(q~Or2 +020r~ + 020~2) Y~.6 = -F6(20,, + q~z + q~4)- F7(0~0,2 + qEG~ + 0~0,4 + 040,1 +
+ 02qr2 + 020r4 + 040,2 + ql4g4 ) + q20,2 + 020r4 + q40,2 + q40~) ~,9 =
0,4
rL14 = 0,1 Y2.2 =--F2(2(qr, +q,2)+qr4)-F3(qlqr4 + q40rl q" 02q, a -t-040,2 +040,4)
I2.3 = Fa(2(qrl + 0,2)+ 0,4)- F2(010,4 +040,1 "t-q20r4 + 1~4qr2 +040,4) r2.4 = -F40,1 + Fsqlq,1 Y2,5 =
Fsq,lF4ihOrl
r2.~ = FTq, l + F6010,1 ~,8 =
0,1 + 4,2
1"2.9 4,4 =
r2.12 q,2 r2,1s = 4,2 =
Y3,1O = qr3 Y3,11 =
-I
Y3,16 = qr3 r,,2 = -F2 (0n + 0,a ) + F3 (010,1 + 01qr2 + 020,1 + 02q,2 )
V4,3 = F3 (0,1 + qr2) + F2 (010r, + q10r2 + Y4,6 = -F6qrl + F7010~1 r4,7 = F70~1+ F6qlOrl
q20rl +
q20r2 )
I'4.9 = 0,1 + 0,2 + q,4 r4.13 = q~4 Yn,17 = 0,4
Appendix 2. Theoretical background for the decomposition Based on the Lagrange method [2,7], the functions D., Diik and Di in the dynamic model can be expressed using the partial angular velocities (itj_ l ) and par'tial velocities (idj_l)in the following form:
Identification and adaptive control of robots
' {ij} { Ks ( , )tj_1 ~ti-1+ m,'dy-1 "d DU=,: m~x ~ " ' i-1 --I'nsPcs "(
stj - I
X
"ai-1
+ "ti-1 X
165
-
"dj_,)},
Oij = Off, Dijk =
{ITrace(Ks)Sti-l'(Stj-lxst~-l)
~
-
s =max{i.k}
K" -~d
X Dt
.'t
+msp~,['tj_l(stk_l • Sdi_l)"d I't
"t
-',,_,(',,_,.'..i)]+
.(',.i
J -<
Oi# = Oiky, Di# =-Okj, if i, k > j, Diji = 0 if i > j , n D i = ~,.ms(-aTo,sg)(Sti_l
$=|
Xpc s
+'d
i-I ) •
To find the minimum dynamic parameters or the influence of the load, it is sometimes better to use the Newton-Euler method. The dynamic model is linear in the robot parameters, and the vectors ~ and def
•..)T.
Especially, the angular acceleration and the acceleration have the form: i
i j-i
iei = Fi4 + ¢I)i44 = ~ , i ty-lqj" + Y~ •
j=l
j=2 k=l
i ai = "Qiq + ~i{]q = i i
i t k-1 X i tj - l q y" q "k ,
• = ~ i aj_lqj- + Z 'tj_~ × ia:_,q).2 + 2 Y j=l
j=l
i
i
"" ~, 'ty_1 X idk-lqyqt.
j = l k=j+l
Let iEi-l.i and iNi-l.i denote the force and torque respectively, applied by link i - 1 to link i and expressed in the frame of link i. Then the Newton-Euler equations have the form: iF.i-,.i = iF.i.i+, + mi(a i _ ig)+ ([eixJ+[ogix][toix])mipc,, i N ,-l.i = i N i,i+l + [ iPi-l,i X] iF. i-l,i +
+ Kiei + Ogi x (Kitoi) + rniPci X(ai - ig), i-iF. iF; i-l,i = mi-l,i i-l.i ,
i-1N
i-l.i = mi-l.i
iN
i-l.i.
If the non modified Denavit-Hartenberg form is used, then the driving torque (force) has the form ~i = < i - I N i-Li, ti-1 > = < i N i-l.i, A Ti-l.i ti-1 > for revolute joints,
~i = < i-IF: i-Li, ti-1 > = < if. i-l.i,
A~_- 1.i ti - 1 > for sliding joints,
166
B. L A N T O S
Especially, if the f r a m e s o f the last link a n d the load are the s a m e ( K n = 9(c), t h e n n Fn,t, = mL(an '*N
rig) + ([~n X] + [O)n X][O)n X])mLPcL, X (Kt. con ) + ml.Pct
-
and t h e i n f l u e n c e o f the load on link i is
]
i-l j | =
|i-I N | L i-l,i _1
Ai_l.,
0
L[pi_,,nx]Ai
In - '
Ai_,. n
"F. nx
n,L
"
REFERENCES [1] S. G. Tzafestas: Adaptive• r•bus•• and fuzzy ru•e-based c•ntr•• •f r•b••ic manipu•a••rs. •n Inte••igem R•b••ic S•stems (Edjted by Tzafestas S. G.), pp. 313~,19. Marcel Dekker Inc., New-York (1991). [2] B. Lantos: Robot Control. Hungarian Academic Press. Budapest (1991), in Hungatian. [3] B. Lantos: Simulation environment for robot control design. IFAC/IFIP/IMACS Symposium on Robot Control SYROCO '91, Vienna, pp. 417-422 (1991) [41 L. Ljang: System Identification: Theory~or the User. Prentice-Hall Inc., Englewood Cliffs (1987). [5] B. Lan••s• •. Ha••sz and Zs. Juh•isz: •denti•cati•n and ••mr•••er design f•r SCARA r•b•ts. Res•arch Repor•' BME-42 •.•22/89• 173 Budapest (1989), in Hungarian. [61 M. Vukobratovic and V. Potkonjak: Dynamics #Manipulation Robots. Springer Verlag (1985). [7] R. P. ••u•: Rob•t Manipula••rs: Ma•hematics• Pr•graming• and C•mr••. The M•T press. Cambridge• M. A. ( •982). [8] P. Csala (1991): Robot model generating in symbolic form and its application in adaptive controls. Doctoral (techn. univ.) Thesis, TU Budapest (1991), in progress. [9] J.E.S~tineandW.Li:Adap~ivemanipu~at~rc~n~r~:acasestudy.~EEE~ru.C~Rob~ticsandAut~mati~n~Ra~eigh~pp. 13921400 (1987). [10] 3. E. Slotine. andW. Li: Adaptive strategies in constrained manipulation. IEEE Int. Con/. Robotics and Automation, Raleigh, pp. 595--601 (1987). [11] T. J. Tam, A. K. Bejczy, X. Yun. and X. Ding: Dynamic equations for six-link Puma 560 robot arm. Robotics Laborawo Report. SSM-RL-86--05. Washington University, St. Luis (1986). [12] B. Lant~s and ~. Csa~: C~nvergence test ~f se~ftuning r~bot c~ntr~. Au1~mation ~92~Budapest~ pp. 4~2~ ~3 ( ~992). [13l H. Kawasaki and K. Kanzaki: Minimum dynamics parameters of robot models. IFAC/IFIP/IMACS Symposium on Robot Control SYROCO '91, Vienna, pp. 111-116 (1991). [14] S. Y. Sheu and M. W. Walker: identifying the independent inertial parameter space o¢ robot manipulators. Int. J. Robotics Research 10, pp, 668-683 (1991).
0957-4158/93 $6.00+0.00 Pergamon Press Ltd
IDENTIFICATION AND ADAPTIVE CONTROL OF ROBOTS B. LAy'ros Technical Universityof Budapest. Departmentof ProcessControl H-1111 Budapest,Mfiegyeternrkp. 9., Hungary (Received 30 June 1992)
Abstract - The paper deals with the model identification of industrial robots as the basis for control algorithm design. The first part discusses the classical decentralized joint control using internal current and velocity loops and external position control, where the plant to be identified is the black-box containing the intero.,tl loops and an integrator. The paper presents the identification of a SCARA robot using LS, GLS and IV methods and real closed loop data under normal working conditions. The second part discusses the convergence test of self'tuning adaptive control for SCARA and PUMA robots using the advanced robot control algorithm of Slotine and Li to find the unknown parameters of the robot dynamic model. The convergence has been tested by using a transputer based simulation system and the symbolic form of the robot dynamic model.
1. INTRODUCTION Robot control systems have usually a hierarchical structure and can be divided into specification level (robot programming language, graphic simulation, path design), algorithmic level and torque control (power electronic) level. From the theoretical point of view it seems to be acceptable to divide the control algorithms into decentralized joint control and mole advanced robot control algorithms. Decentralized joint control attempts to control every joint axis separately without considering the coupling effects between the different joint axes. To this class belong the three loop current, velocity and position contro!, the PID control of the position with motor voltage controller output, etc. To the class of advanced control algorithms belong such methods like computed torque techniques, resolved motion acceleration control, hybrid position and force control, model reference adaptive control, selftuning adaptive control, fuzzy control, etc. [1, 2, 3]. Characteristic for this second class of control algorithms is that they are based on the nonlinear dynamic model of the robot. The paper discusses parameter estimation methods to find the average linear model for decentralized joint control, and selftuning adaptive control to find the nonlinear dynamic model during control.
2. IDENTIFICATION OF THE LINEARIZED MODEL FOR DECENTRALIZED POSITION CONTROL
Industrial robot control systems often apply internal current (torque) control, central velocity (revolution) control and outer position (joint variable) control for the robot axes. It is assumed that the actuators are DC motors with gear reductions.The controllers for the current and velocity loops are analog PI con149
150
B. LANTOS
trollers, while the position controllers are DDC controllers containing PID or PI controllers and a feedforward path for the desired velocity. If the mechanical and electrical parameters of the motors and the gear reductions are known then the following method can be suggested for the design of the decentralized joint controllers. For the current control loop we can consider the backward electromotor force as a disturbance, and the design problem reduces to the design of a PI controller for a first order system with the electrical time constant. The integrator assures that the gravity load is compensated in steady state. After having designed the current controller, the backward electromotor force can be taken into consideration using root locus method, and the dominant low frequency pole can be determined for the velocity control design, whose place is very critical for the quality of the whole control system. This dominant pole can be cancelled by the zero of the PI controller of the velocity loop, and the controller gain can be chosen to result nearly critically damped closed loop. The resulting system after velocity control can be approximated by a time invariant third order linear system and an extra integrator. This is only an average model, because a nonlinear system linearized along a time varying path yields always a time variant linear system. An excellent treatment of identification methods can be found in [4]. The identification problem in discrete time is now the following. Let z -~ denote the backward shift operator. The linear, time invariant approximation of the black-box containing the current and velocity control loops of a robot axis can be assumed in the form of a discrete time SISO (single input x, single output y) system
B(z-~)
y(t) = A - - ~
z-a x(t),
alz-l+...Wam Z-m,
:
1
=
b~z-l+...+b,,z -m,
+
y(,) = -aly( t - 1 ) . . . . . amY( t
(1) - -
m).4_
+ blx(t - d - 1)+'"+bmx(t - d - m), where r is the sampling time, d . r is the time delay (if any) and t is an integer multiple of the sampling time. Let p and hr (t) denote the system parameter and the observation respectively, where p = (a,, o~ ..... a,,, b~, b~ ..... bin) ~, hT(t) = ( - y ( t - 1) ..... - y ( t - m), x ( t - d - 1) ..... x ( t - d - m)), then
y(t) = hr (t)p.
(2)
In the case of the least-squares method (LS) the system parameter p is the analytical solution of the optimum problem N
2
V ( p ) : ,=,,~d+l { Y ( t ) - h r ( t ) P } - - 4
min,
(3)
where V(p) is the sum of the squares of the equation errors, and N is the number of observations. Using the notations
Identification and adaptive control o f robots Y=(y(m+d+l)
15 i
. . . . . y ( N ) ) T,
,,__ Ibm(m+ =
- y ( m + d)
....
!
:
-y(N)
....
l
y ( d + 1)
x(,,,)
...
x(1)
x(U-d-1)
...
x(U-d-m)
]
i
y ( N - m)
!
!
1'
the LS optimum can be written in the form /3 = (H/-H) -1HTY.
(4)
If the system output contains an additive noise v(t), then y(t) = h ~ ( t ) p + v(t). Let us introduce the notations |
N
q)~-N-d-m
y~ h ( t ) h r ( t ) = 1 HrH, t=m+d+l N - d - m 1
qghv - U _ d _
N
m
y~ h ( t ) v ( t ) t=m+d+l
N
-
1 HTv, d - m
then ¢p~ and ¢Phv are convergent to the mean values ¢Phh ---) E { h ( ' ) hT (')} and ¢Phv ~ E { h ( t ) v ( t ) } respectively, and p satisfies b = P + ~0~,]¢h~•
(5)
Hence/3 --~ p asymptotical if h(t) and v(t) are uncorrelated. This is always satisfied if v(t) is white noise. Otherwise v(t) may be a coloured noise, and the estimate/3 is not necessarily convergent to the system parameter p. The error is dependent on the signal to noise ratio. Since the identification of the robot is based on closed loop data, we can expect coloured noise and biased estimation. Hence the application of instrumental variables (IV) or generalized least squares (GLS) can be suggested. The method of instrumental variables substitutes the LS estimation by the relation /3 = ( Z z H ) - I Z r Y .
(6)
The instrumental variable z(t) should satisfy that ¢Pzh is invertable and ~Pzvhas zero mean value, where 1
¢Pzh - N - d - m _
qg~v
1
N- d - m
N
y~ z ( t ) h r ( t ) = 1 ZrH, t=m+d+l N - d - m N
• z(t)v(t) = 1 ZTV. t=m+d+l N - d - m
152
B. LANTOS
The instrumental variables can be generated by the auxiliary model
f~= ~z Z--
x(t),
z~(m + d + i
1)] ]=
(7)
Z~"(N) + d) =
....
:
:
-~(N-1)
...
+ 1)
41)
]
:
....
fJ(N-m) x(N-d-1) :
: x(n-d-m) :
...
J"
The auxiliary model can be initialled by using LS estimation. However the model is not exact, thus and v(t ) and similarly y(t) and v(t) are correlated. Hence a sequence of estimations should be applied:
z(t ) = h(t, ~, x)
1. LS estimation---) A,/}.
2. Generating the instrumental variables:/]33(t) =
Bx(t ).
3. IV estimation---> ,~,/~.
4. Jump to 2. if the parameters are not convergent, otherwise stop. The method of the generalized least squares assumes coloured noise v(t) produced from the white noise e(t) by the filter C(z -1 ). The system model can be assumed in the form
a(z-1)y(t) = B(z-l)z-dx(t)+C(z-~)e(t) a(z-~)~(t) = B(z-1)z-aYc(t)+e(t), ~(t) = c-l(z-1)y(t), Yc(t)=C-l(z-1)x(l).
or (8)
Since C(z -1) is usually unknown, thus C(z-1), a(z -1) and B(z -1) have to be identified simultaneously. The algorithm has been implemented assuming C(z -1 ) = 1 / F ( z -l): 1. Let/~(z-1) = 1. 2. Determine ~(t) and ~(t) using
F(z-' ): ~:(t) =/~(z-' )y(t)
and Yc(t) = l:(z-'
)x(t).
3. Determine the parameters o f / t ( z - ' ) and/~(z-') from :~(,) and 5'(t) by using LS estimation.
A(z-1)y(t)-B(z-~)x(t-d). 5. Apply LS estimation for the model F(z -1)v(t) = e(t), ( a : = F, B:= 0, y(t): = v(t), and no x(t) is present),
4. Determine v(t) =
where e(t) is a white noise to find the new value of/6.
6. Jump to 2. if the parameters are not convergent otherwise stop.
Identification and adaptive control of robots
153
During the structure estimation the purpose is to find the system order (m) and the time delay (d). In the case of LS estimation the loss function can be used for this aim. For increasing system order LS estimations are performed for the same observations {xn, Yn}, and the equation errors normalized to one step (the residuals) are determined for every system order. That system order has to be found, for which the residual can not be essentially improved by increasing the system order. In the case of IV or GLS estimation the whiteness of the residuals is necessary to the unbiased estimation, that can be tested using the autocorrelation function. It may be another practical method to compare the sum of the squares of the differences between measured system output and computed model output for different model orders. The efficiency of the identification of the linear approximation will be demonstrated for a SCARA robot. The robot repeats a composite motion cycle, so that it is useful to design new controller parameters for the position controllers of the different axes. The base for the controller design may be the approximated linear model of the resulting subsystem containing the internal velocity and current control loops of the robot axis. The subsystem begins at the D/A converter (DAin) after the position controllers and is finished after the encoder and counter converting the joint position into large resolution integer values (00 or qhO). The robot has been sent through the path formulated by the robot programming language and using the actual (not necessarily suitable) controller parameters. During this normal work 4000 closed loop data x(t) := Dain(t) and y(t):= 0(t) have been collected with 1 ms sampling time and stored for later identification. The identification has been performed in batch mode in a host computer using LS, IV and GLS methods. The measured output of the subsystem (qhO)and the answer of the identified models (qhidO)have been computed and compared for different system orders and methods, and the best result has been accepted. The best approximation has been achieved for GLS method and order 4 without time delay [5]. For the second robot axis Fig. 1 shows the subsystem input, the subsystem output and the answer of the identified model for the same input. It can be seen that the approximation is excellent, although the energy content of the input signal is not high. The discrete time model of the subsystem and the equivalent continuous time model for step inputs (for zero order hold) are respectively
Gtz_l~ = G(s) =
--0.007195z -l + 0. 005618z -2 + 0.004719z -3 + 0.008259z -4 1 - 0. 828761z -I - 0. 281923z -2 - 0. 481865z -3 + 0. 592652z -4 ' 24.3666 - 0.049844s + 1. 24168.10 4 s 2 - 1.23455.10 -8 s 3
(9) (10)
s(I+Tls)(l+2VTzs+Tz2s 2) T1 = 0.00528626 sec,
D = 0.0756408,
Tz = 0.000052803 sec.
The continuous time subsystem is of non-minimum phase (there is a zero at s = 469), and has an approximately cancelled conjugate complex pole (-167 + j2202) and zero (-18 4 + j2043) pair in the high frequency domain. To avoid numerical instability caused by the large amount of data, the Householder method has been applied for matrix inversion. Moreover, it was necessary to test the poles of the discrete time model in different steps of the algorithm, and correct the pole in the neighborhood of z = 1 if it became numerically unstable. A built in program has supported the controller design [5]. The position controller design can be elaborated immediately in discrete time using Tustin (w) transformation, root locus method or state feedback. In this case, it is less critical if the discrete time model contains a negative real pole in z with an odd number of multiplicity. However the preferred control algorithms for position control were PIF or PIDF control, where F denotes the feedforward path from the desired path velocity 0a to the set point input of the velocity control. The controller design has been performed in continuous time. Since in this case the open loop contains
154
B. LANTOS DAin [2] qhO [2] qhldO [2]
0 "A C]
2~
20
.;J C 0 tO O.
,,~
t5
io
(.. t,-
o W e.-
-5
e-I
0
5
~0
15 t
20 E5 [ E -:L see I
30
35
40
Fig. 1. ldentifiLation of the hm~arized model
two integrators, hence one zero of the controller should be placed in the low frequency domain and a large phase margin has to be designed. The continuous time PID controller has been implemented in discrete time using its step equivalent form. Although the identification and controller design has been performed in batch mode on a host computer, the whole process can be implemented in real time on a multiprocessor control architecture.
3. SELFTUNING ADAPTIVE CONTROL The base for advanced robot control algorithms is the robot dynamic model, which gives the relation between driving torque (z), inertia matrix (H) and the centripetal, Coriolis and friction effects (CO), and the gravitational part (D):
H(q)it + h(q, (t) = H(q)~ + C(q, dl){l + D(q)
= z.
(11)
Using the mass mi, the center of mass Pci and inertia moment Kei belonging to a coordinate system having origin in the center of mass but axes parallel to the Denavit-Hartenberg frame Ki, the quantities in the dynamic model can be expressed using recurrence or direct formulas [2, 6]. Similarly, H and h can
Identification and adaptive control of robots
155
be determined also by using the inertia moment Ki expressed in the Denavit-Hartenberg frame Ki. The later case is preferable because the motion equation becomes linear in mi, m~Pci and Ki. The main problem of the realization of advanced robot control algorithms is caused by the fact that the driving torque should be computed by the controller within 1-10 ms. Hence the controller should be realized by using a multiprocessor architecture containing fast arithmetical processors which allow parallel computation. Preferable to parallel computation are nonrecursive algorithms for the computation of the dynamic model H(q)~ + C(q, Cl)Z + D(q) = "t"and Y(q, q, ~, ~)a = r so that the time consuming computations can be distributed among the processors [2]. A special list processing program has been developed in Pascal for IBM-PC/AT which determines H, h and o~ in symbolic form [8]. If the robot dynamic model is known in symbolic form where the inertia parameters (mass, center of mass, inertia moments) are also contained in symbolic form, then independent parameters ot1, ot2 .... can be chosen so that the dynamic model is linear in the parameters:
H(q)~ + C(q, (I)~ + O(q) = Y(q, q, ~, ~)a = "r.
(12)
The identification problem is the determination of the actual parameter value a (Slotine and Li [9]). The method can be generalized to Cartesian and hybrid control (Slotine and Li [10]). Let q,(t) be the reference signal defined later and s = 0r - q. The control law is chosen as "t:=/~/(q)qr + C(q, (l)¢r
+
D(q) + Kt,(q r - q) + KoS,
(13)
where/q, C, b is computed at the estimated & and measured q, q. The Lyapunov-like function
v=lsTH(q)s
1 + ~(qr - q)TKe(q, - q) + I ( ~ - a)T/"(t~ -- a )
(14)
Can be USed, where K e, K D, e are positive definite symmetric matrices. If a does not change during the adaptation process then the adaptation law
d_.._~= F_lyT (q, (l, glr, qr) s dt
(15)
can be chosen and
(/= -sT Kos <- 0 is satisfied. If the reference signal is t
qr := qa + A S(qa - q ) dr ~ s = ((la -dl)+ A(qa - q )
(16)
0
and the roots of det (sl + A) = 0 lay in the left half of the complex plane then the trajectory is convergent to the sliding surface s = 0 and the error qa - q --~ 0 exponentially on them. Hence the control (13) with the parameter estimation (15) is globally asymptotically stable. Since the mass, center of mass and inertia matrix of the robot links are usually unknown, the convergence test of the selftuning adaptive control method has been performed on simulated robot models, allowing the comparison of known model parameters and their estimations. Multitransputer based simulation system has been developed to examine a broad class of control algorithms containing usual industrial robot control (position, velocity and current control), computed torque technics and decentralized PID controllers, resolved motion acceleration control, hybrid position and force control, selftuning adaptive control and fuzzy control [3]. It is assumed that the control architecture contains low level torque (current) control and the actuators are DC motors with gear reductions.
156
B. LANTOS
The controllers for the current loops are analog PI controllers. The mass and inertia parameters of the link and the actuator (or load) mounted on the link are reduced to resulting parameters, and the effect of the kinetical energy of the rotor is approximated by ~ t g , i ( v~ is the gear reduction and Ori is the rotor inertia). The friction of the rotor and the gear reduction is approximated by v2f~i + fg,ar,,. Let L,i, Rri, qi and c2i denote the rotor inductivity, rotor resistance, electromotor force constant and torque constant respectively. The coupled mechanical and electrical equations have the form EDij(q)dlj + EEDijk(q)(lj~lk + D i ( q ) + ~iit)riqi +(V2fri + fgear,i)(]i = ViC2iiri, j j k
diri
Lri ~
~ .
+ Krilri = Ui -- Cli Vi~ti '
(17) (18)
where i~ is the rotor current and u~ is the motor voltage. The PI controllers are realized by operational amplifiers, therefore the set-points [a are voltages. The simulation is based on the robot dynamic model in symbolic form. Using the substitutions Hii = Dii := Dii + v26)ri, hi "=/1/ +(V2fri + fsear.i){ti, and the state variables q, q, i~, 1 the system can be written in the form of dynamic and static equations:
dq_
dt-q' d___~q= H_1{coll (vic2,i, i _/~)}, dt di__L~ = co. { - c. v,o, - R.,i.,)}. dt
(19)
d ! = 7~ - T = i~ - coil ( r, ir, ), dt u:coll
apii ]
1
The torque output r of the high level controller is computed in the sampling instants and converted to ~'e.After D/A conversion their components are sent to the reference signal inputs (set point functions) of the current loops. Further state equations have to be added to the system if the sensors contain time constants. Saturations in the D/A converters, the operational amplifiers, the motors and the sensors can be taken into consideration during simulation.
3.1. Convergence Test for SCARA Robot The dynamic model of the 4dof RRTR SCARA robot has been generated using symbolic computation. The functions Dij, Dok and D i of the robot model are sums of terms. Every term is the product of a constant, an expression of the parameters and an elementary function (a product of powers of Si = sin qi, 6",. = cos qi or qi, depending on the joint type). The parameters /'1 . . . . . Pn of the links should be corrected by the rotor inertia, gear reduction and friction:
P1 = m2 a2 + P2xzm2 + p2xlml + p2y2m2 + rn3a2 + 2pcxlrnlal + + 2pcxzm2a2 + mla~ + P~xsm3 + m4a 2 + p2y3m3 + Kczl + 2p~x3m3a2 + + gcz 2 + p2x4m 4 + p2y4m 4 + p2ylrrl1 + gcz 3 + rn:2a~ + gcz 4 + + rn3a2 + m4 a2 + Orl V2
Identification and adaptive control o f robots
157
P2 = Pcy4m4a2 P3 = ,Dcx4m4a2
P4 = P5 P6 P7 P8
(P~2rr~+ Pcy3ma)a,
= (Pcy2ra2 + m2a2 + Pcx3m3 +m3a2 + m4a2)al = Pcy4m4al = Pcx4m4al = p2x2m2 + p23m3 + 2Pcx2rn2a2 + 2Pcx3rn3a2 + p2y2m2 + P2yam3 +
+ m2a~ + p2x4m4 + m3a~ + p2y4m4 + Kcz2 + Kcz3 + Kcz4 + rl'14a~ P9 = p2x4m4 nt"p2cy4rl'14 -I- gcz 4
t'1o
(20)
= m3 + m4 + e,3 ~
Pll = (m4 + m3)gz
P12 = o,s P13 =
(~r4 V~4
P14 = e15 =
f2
P16 = f3
/'17
=
A
Appendix 1 contains the dynamic model of the SCARA robot. The simulation system makes it possible to distinguish the dynamic parameter values of the robot (for computations of the sensor signals using Runge-Kutta method; parameters P~) from the dynamic parameters of the selftuning adaptive controllers (the estimated parameters Bi). Commands allow to modify Ke, Ko, A and F for the adaptive controller during the identification. The convergence test has been performed for initial parameters Bi (0) scaled up by 500%. Because the robot is of open chain type, the identification has been elaborted in a backward order of the links. The influence of the parameters onto the driving torques of the links has been formulated in a parameter graph supporting the tuning strategy. Since the command set of the simulation system allows to switch off the tuning for some parameters in different (acceleration, constant velocity, steady state) path sections, it was possible to find the best tuning strategy by using the relations in the parameter graph. It is a real problem of error based adaptation, that the control error becomes small if the parameter adaptation converges, but the energy of the error signal is not enough to continue the tuning. Therefore, in every stage of the adaptation, the parameters have been divided into two groups: parameters for which the tuning has already been finished, and parameters to be tuned simultaneously. If the tuning for some parameters of the second group was finished, then the user could decide whether the actual, the initial or a previously stored value has to be used for a remaining parameter in the second group to continue the identification. The path, control error, set point for current control and parameter tuning are shown in Figs 2-5 for 5 ms sampling rate. The real (P~), the estimated (B,) parameters and the relative accuracy (Pi / Bi in %) are shown in Table 1. It can be seen that the accuracy of all the dominant parameters is satisfactory, the control error is small, and there is no saturation in the control system.
158
B. L A N T O S .do t~] qd0 {2] Qd0 tS] ado [4]
~ O ~ , A ~ O ~ X
7
8 B 4
3
|
o t
[ see ]
Fig. 2. Path signals/or the S C ~ A robot
O
dq [i] d. [a] dtl [a] a . t4]
~
x
[]
L
2.5 0.0
~ -2.5 ~-5.0 -7.5
/A
-~,-iO. 0
/
c~-t2.5 i 0
[ ne ] Fig. 3. Error signals of the SCARA robot
Identification and adaptive control o f robots Idh 11] lab [2] Idh[3] im~t4]
159
~ 0 ~ A ~ 0 ~ X
30 25 20
o
~5 -1o .o
-15
/o
-~0
o
1 t
[ see ]
Fig. 4. Adaptive controller outputs of the SCARA robot
st11] A [14] s [15] B [17]
\x -5
-1o
-20 -25 -30
0 t
[ see |
Fig. 5a. Parameter tuning of the SCARA robot
160
B. LANTOS n
El( t1 B[ 4)
~
tat 51
O (
"
]
s t el S tt2] e (~3] 22.5 20.0 :17.5 J.5.0 t2.5 t0.0 7.5 5.0 ~.5
/
0.0 t
[ aec )
Fig 5b Parameter tuning of the SCARA robot
e ['~6]
tB
t6
to
\ o
~
~
~
.~
5"
"~ t
~
~
~
[ aec ]
Fig. 5c. Parameter tuning of the SCARA robot
16
ti
~
Identification and adaptive control o f robots
161
TABLE 1 Parameter Convergence for S C A R A Robot real
estimated
P1 = 4.230 P2 = 0 P3 = 0 P4 = - 0.162 P5 = 0.947 P6 = 0 P7 = 0 Ps = 0.727 P9 = 0.0040 P10 = 130.252 P11 = - 6 7 . 1 9 9 P12 = 0.419 PI3 = 0.405 P14 = 14.500 P15 = 13.800 P16 = 3948.866 P17 = 13.400
B 1 = 4.494 B2=0 B3 = 0 B 4 = - 0.515 B 5 = 0.866 B6 = 0 B7=0 Bs = 0.765 B 9 = 0.0044 Blo = 131.577 Bn = - 6 7 . 1 6 5 B12 = 0.746 B13 = 0.522 Bla = 14.488 B15 = 13.788 B16 = 4082.003 BI7 = 13.402
accuracy 94.1%
31.5% 109.4%
95.0% 91.0% 99.0% 100.1% 56.1% 77.6% 100.1% 100.1% 96.7% 100.0%
3.2. Convergence Test for 3dof PUMA Arm The dynamic model of the first 3 degree of freedom of the PUMA 560 robot has been generated using symbolic computation [8]. The parameters P1. . . . . P21 of the links have been corrected by the rotor inertia, gear reduction and friction:
P1 = 2pcy3m3d3 + Kyl + rrt3d2 + Orl v? P2 = m3a~ + Ky2 + 2pcx2m2a~2 + m2a~ P3 = 2pcx3m3a3 + m3a~ + Kz3 P4 = Pcz3m3a3 P5 = Kx2 P6 =
K~3
P7 = Pcy2m2a2 Ps = (Pox3 + as)m3a2
P9 = Pcz3m3a2 r,o = (pcy3a3 + Pcx3d3 + a3a3 )m3 PH = (Pcy3a3 + m3d3 +Pcz2a2)a2 P12 =
Pcz3m3d3
P13 = m3 a2 + 2pcx3m3a3 + ,ocx2m2a2 + K~2 + Ky3 + m2a 2 + m3a 2 + 6)r2~ P14 = 2pcx3msa3 + Ky3 + rosa 2
P15 = 2#cx3m3a3 + m3a2 + Kz3 - Kx3
(21)
162
B. L A N T O S /°16 -- m34 + rr~a2 - Kx2 + 2pcx2m2a2 + Ky2 Ply = - ( Pcz3m3a3 - ,Ocyzm2a2) /°18 = (a3 +,Ocx3)m3gz P19 = (m3a2 + Pcx2m2 + m2a2)g~ P20 = ,Ocy2m2gz P21 = Pczam3gz
P22 = or3 v~ P23 = fl
P24 = f2 P25 = f3 The parameters P1. . . . . P25 have been computed from the geometrical, inertial and motor parameters of a real P U M A 560 robot (Tam et al. [11 ]). The convergence test has been performed for initial parameters B i (0) scaled up by 500%, and in backward order of the links. A suitable tuning strategy has been found by using the parameter graph of the 3dof P U M A ann. The convergence of the parameter tuning is shown in Table 2 for 1 ms sampling rate. The accuracy of the estimated dominant parameters is satisfactory [12]. It has been found, using simulation support to find suitable controller parameters, that the control error is also small without any saturation in the system. TABLE 2 Parameter Convergence for P U M A Robot real P1 = 1.151 P2 = 3.434 P3 = 0.014 P4 = - 0.020 P5 = 1.333 P6 = 0.546 P7 = 0.048 P8 = - 0.017 P9 = 0.458 P10 = - 0.007 P n = 2.449 P~2 = 0.160 Pi3 = 3.333 P14 = 0.547 P15 = - 0.532 P16=2.101 PI7 = 0.069 P18 = 0.383 P19 = - 53.115 P20 = - t.119 P21 = - 10.608 P2z = 0.583 P23 = 14.500 P24 = 14.300 P2s = 19.300
estimated B 1 = 1.092 B 2 = 3.278 B 3 = 0.016 B4 = - 0.029 B 5 = 1.646 B 6 = 0.523 B 7 = 0.027 B8 = - 0.033 B9 = 0.418 B10 = - 0.016 Bll = 2.393 B12 = 0.196 B13 = 3.160 BI4 = 0.497 B15 = - 0.693 n16 = 2.618 Bl7 = 0.043 n18 = 0.387 n19 = - 5 3 . 1 1 3 B2o = - 1 . t 1 8 B21 = - 10.607 822 = 0.569 B23 = 14.178 B24 = 14.295 B25 = 18.794
accuracy 105.4% t04.7% 83.9% 69.1% 81.0% 104.3% 175.7% 50.2% 109.5% 44.0% 102.4% 81.4% 105.5% 109.9% 76.8% 80.3% 159.7% 99.2% 100.0% 100.1% 100.0% 102.4% 102.2% 100.0% 102.7%
Identification and adaptive control of robots
163
For the 6dof case a multistep identification process can be suggested. The first step may be the identification of the parameters of the first 3dof of the robot for fixed q4, q5, q6. The parameters will not be unique because of the dependence on the fixed q4, qs, q6 configuration. The second step is the identification of the parameters of the 6dof robot using the previously identified dominant parameters for initial values. The first and second step can be performed in a laboratory. The third step should be repeated during the robot application if the influence of an unknown load (mass, center of mass and inertia moments) has to be identified. In this third step the identification of the robot parameters is switched off and only the identifiable load parameters are tuned. Appendix 2 summarizes some methods for model decomposition. The minimum dynamics parameters can be found by using the method of Kawasaki and Kanzaki [ 13], if the robot is described in frames based on the modified Denavit-Hartenberg notation. Determining the independent inertial parameters combined with model reduction can be performed by using the method of Sheu and Walker [14]. 4. CONCLUSIONS Model identification and selftuning adaptive control of industrial robots have been discussed in the paper. The identification of the average linearized model of the subsystem containing the internal current and velocity control loops has been presented using LS, IV and GLS parameter estimation methods for a SCARA robot based on real closed loop data under normal working conditions.The model can be applied to the design of the position controller in decentralized joint control systems. The selftuning adaptive control algorithm of Slotine and Li [9] has been applied with some practical modifications to the SCARA robot and the first 3 degree of freedom of the PUMA 560 arm. Simulation results have experimentally proved the convergence of parameter identification. The control errors were satisfactorily small during the adaptation process. Multitransputer based simulation system has supported the design of the controller parameters avoiding saturation. The simulation system turned out to be an easy tool to find the favourable tuning strategy for real-time implementation in the future. Acknowledgement - Support for the research of the model identification to decentralized joint control was provided by MIKROMATIKA Ltd. Budapest under grant No. BME--421.022/89. Support for the research of selftuning adaptive control is provided by the Hungarian National Research Program under grant No. OTKA-756. Appendix 1. Dynamic model of the SCARA robot Using the elementary functions El = 1, F 2 = $ 4 , F 3 = C 4, F 4 = S 2, F 5 = C 2, F 6 = $24, F 7 = C24, the symmetry properties of Di/ and Dijt, and the robot parameters '°1. . . . . P17, the nonzero elements in the robot dynamic model can be written in the following form: D n = F~P1 - 2F2P2 + 2F3P3 - 2F4P4 + 2FsP5 - 2F6P6 + 2FTP7 i)12 = F1Ps - F4P4 - 2F2P2 + 2F3P3 + FsPs - F6P6 + F7P7 014 = -&1"2 + Fa1"9 + F31"3 - F6e6 + FTP7 Dz2 = F1P8 - 2 F 2 Pz + 2F3P3 /)24 = - F 2 P 2 + E1P9 + F3P3 033 = ElelO 044 = F~e9 0112 = -0211 = 0122 = -F41"5 - F61"7 - F71"6 - &1"4 Oll 4 = -O411 ~-- O124 -~. O144 = - F 2 e 3 - F3P2 - f 6 e 7 - f T e 6 O214 = -O412 = 0224 = - 0 4 2 2 = & 4 4 = -&1"3 - 5 1 ' 2 03 = - E l e n
164
B. LANTOS El,1 = qrl El,2 = - F 2 ( 2 ( q r l
q-qr2)+Or4)-f3(t~lOr4
"t-q40rl "1"q20r4 + q40r2 "b q40r4)
Y1,3 = F3 (2(qrl + qr2 ) + qr4 ) -- F2 (qlOr4 + O,qrl 4- q20r4 + 04qr2 + 04qr, )
YL4 = - & (20,, + 0,2) - Fs(q,G2 + 020,~ + 0242) r~.5 = 6 ( 2 0 , , +O,2)-G(q~Or2 +020r~ + 020~2) Y~.6 = -F6(20,, + q~z + q~4)- F7(0~0,2 + qEG~ + 0~0,4 + 040,1 +
+ 02qr2 + 020r4 + 040,2 + ql4g4 ) + q20,2 + 020r4 + q40,2 + q40~) ~,9 =
0,4
rL14 = 0,1 Y2.2 =--F2(2(qr, +q,2)+qr4)-F3(qlqr4 + q40rl q" 02q, a -t-040,2 +040,4)
I2.3 = Fa(2(qrl + 0,2)+ 0,4)- F2(010,4 +040,1 "t-q20r4 + 1~4qr2 +040,4) r2.4 = -F40,1 + Fsqlq,1 Y2,5 =
Fsq,lF4ihOrl
r2.~ = FTq, l + F6010,1 ~,8 =
0,1 + 4,2
1"2.9 4,4 =
r2.12 q,2 r2,1s = 4,2 =
Y3,1O = qr3 Y3,11 =
-I
Y3,16 = qr3 r,,2 = -F2 (0n + 0,a ) + F3 (010,1 + 01qr2 + 020,1 + 02q,2 )
V4,3 = F3 (0,1 + qr2) + F2 (010r, + q10r2 + Y4,6 = -F6qrl + F7010~1 r4,7 = F70~1+ F6qlOrl
q20rl +
q20r2 )
I'4.9 = 0,1 + 0,2 + q,4 r4.13 = q~4 Yn,17 = 0,4
Appendix 2. Theoretical background for the decomposition Based on the Lagrange method [2,7], the functions D., Diik and Di in the dynamic model can be expressed using the partial angular velocities (itj_ l ) and par'tial velocities (idj_l)in the following form:
Identification and adaptive control of robots
' {ij} { Ks ( , )tj_1 ~ti-1+ m,'dy-1 "d DU=,: m~x ~ " ' i-1 --I'nsPcs "(
stj - I
X
"ai-1
+ "ti-1 X
165
-
"dj_,)},
Oij = Off, Dijk =
{ITrace(Ks)Sti-l'(Stj-lxst~-l)
~
-
s =max{i.k}
K" -~d
X Dt
.'t
+msp~,['tj_l(stk_l • Sdi_l)"d I't
"t
-',,_,(',,_,.'..i)]+
.(',.i
J -<
Oi# = Oiky, Di# =-Okj, if i, k > j, Diji = 0 if i > j , n D i = ~,.ms(-aTo,sg)(Sti_l
$=|
Xpc s
+'d
i-I ) •
To find the minimum dynamic parameters or the influence of the load, it is sometimes better to use the Newton-Euler method. The dynamic model is linear in the robot parameters, and the vectors ~ and def
•..)T.
Especially, the angular acceleration and the acceleration have the form: i
i j-i
iei = Fi4 + ¢I)i44 = ~ , i ty-lqj" + Y~ •
j=l
j=2 k=l
i ai = "Qiq + ~i{]q = i i
i t k-1 X i tj - l q y" q "k ,
• = ~ i aj_lqj- + Z 'tj_~ × ia:_,q).2 + 2 Y j=l
j=l
i
i
"" ~, 'ty_1 X idk-lqyqt.
j = l k=j+l
Let iEi-l.i and iNi-l.i denote the force and torque respectively, applied by link i - 1 to link i and expressed in the frame of link i. Then the Newton-Euler equations have the form: iF.i-,.i = iF.i.i+, + mi(a i _ ig)+ ([eixJ+[ogix][toix])mipc,, i N ,-l.i = i N i,i+l + [ iPi-l,i X] iF. i-l,i +
+ Kiei + Ogi x (Kitoi) + rniPci X(ai - ig), i-iF. iF; i-l,i = mi-l,i i-l.i ,
i-1N
i-l.i = mi-l.i
iN
i-l.i.
If the non modified Denavit-Hartenberg form is used, then the driving torque (force) has the form ~i = < i - I N i-Li, ti-1 > = < i N i-l.i, A Ti-l.i ti-1 > for revolute joints,
~i = < i-IF: i-Li, ti-1 > = < if. i-l.i,
A~_- 1.i ti - 1 > for sliding joints,
166
B. L A N T O S
Especially, if the f r a m e s o f the last link a n d the load are the s a m e ( K n = 9(c), t h e n n Fn,t, = mL(an '*N
rig) + ([~n X] + [O)n X][O)n X])mLPcL, X (Kt. con ) + ml.Pct
-
and t h e i n f l u e n c e o f the load on link i is
]
i-l j | =
|i-I N | L i-l,i _1
Ai_l.,
0
L[pi_,,nx]Ai
In - '
Ai_,. n
"F. nx
n,L
"
REFERENCES [1] S. G. Tzafestas: Adaptive• r•bus•• and fuzzy ru•e-based c•ntr•• •f r•b••ic manipu•a••rs. •n Inte••igem R•b••ic S•stems (Edjted by Tzafestas S. G.), pp. 313~,19. Marcel Dekker Inc., New-York (1991). [2] B. Lantos: Robot Control. Hungarian Academic Press. Budapest (1991), in Hungatian. [3] B. Lantos: Simulation environment for robot control design. IFAC/IFIP/IMACS Symposium on Robot Control SYROCO '91, Vienna, pp. 417-422 (1991) [41 L. Ljang: System Identification: Theory~or the User. Prentice-Hall Inc., Englewood Cliffs (1987). [5] B. Lan••s• •. Ha••sz and Zs. Juh•isz: •denti•cati•n and ••mr•••er design f•r SCARA r•b•ts. Res•arch Repor•' BME-42 •.•22/89• 173 Budapest (1989), in Hungarian. [61 M. Vukobratovic and V. Potkonjak: Dynamics #Manipulation Robots. Springer Verlag (1985). [7] R. P. ••u•: Rob•t Manipula••rs: Ma•hematics• Pr•graming• and C•mr••. The M•T press. Cambridge• M. A. ( •982). [8] P. Csala (1991): Robot model generating in symbolic form and its application in adaptive controls. Doctoral (techn. univ.) Thesis, TU Budapest (1991), in progress. [9] J.E.S~tineandW.Li:Adap~ivemanipu~at~rc~n~r~:acasestudy.~EEE~ru.C~Rob~ticsandAut~mati~n~Ra~eigh~pp. 13921400 (1987). [10] 3. E. Slotine. andW. Li: Adaptive strategies in constrained manipulation. IEEE Int. Con/. Robotics and Automation, Raleigh, pp. 595--601 (1987). [11] T. J. Tam, A. K. Bejczy, X. Yun. and X. Ding: Dynamic equations for six-link Puma 560 robot arm. Robotics Laborawo Report. SSM-RL-86--05. Washington University, St. Luis (1986). [12] B. Lant~s and ~. Csa~: C~nvergence test ~f se~ftuning r~bot c~ntr~. Au1~mation ~92~Budapest~ pp. 4~2~ ~3 ( ~992). [13l H. Kawasaki and K. Kanzaki: Minimum dynamics parameters of robot models. IFAC/IFIP/IMACS Symposium on Robot Control SYROCO '91, Vienna, pp. 111-116 (1991). [14] S. Y. Sheu and M. W. Walker: identifying the independent inertial parameter space o¢ robot manipulators. Int. J. Robotics Research 10, pp, 668-683 (1991).