- Email: [email protected]
Fuzzy parameter identification of direct drive robots
J. FranklinInst. Vol. 333(B), No. 2, pp. 289-310, 1996
~ ) Pergamon
S0016--0032(96)00026-9
Copyright © 1996 The Franklin Institute Published by Elsevier Science Ltd Printed in Great Britain 0016-0032/96 $15.00+0.00
Fuzzy ParameterIdentification of Direct Drive Robots by CHIA-JU
WU
and
HUANG
CHING-HUO
Department o f Electrical Engineering, National Yunlin Institute o f Technology, Touliu, Yunlin, 640 Taiwan (Received 22 August 1994," accepted 31 January 1996)
ABSTRACT: ,4fuzzy identification method is proposed to estimate parameters in dynamics of direct drive robots. In this method, linear membership functions are used and the output of each fuzzy rule is assumed to be the linear combination of input variables. The input-output data for identification are generated in a recursive manner. The number and the ranges of memberships functions are determined by a heuristic searching method and a statistical training method, respectively, in the sense of least squared output errors. Once the parameters in dynamics are estimated, existing control strategies, such as the well-known computed-torque method, can be applied to make the robot follow a desired trajectory. Simulation examples are included to illustrate the validity of the proposed method. Copyright © 1996 Published by Elsevier Science Ltd L Introduction
Research of robots has been prompted in recent years by a wide variety of goals. One of the goals is to develop control algorithms to make a robot move along a desired trajectory. Among the trajectory-tracking algorithms that have been proposed previously (1, 2), the computed-torque method provides a firm mathematical frame and yields a family of easy-to-understand control schemes. These schemes perform well provided that the parameters in dynamics of robots are known fairly accurately. However, in practice, there are uncertainties in the dynamics of robots. This may arise because the exact evaluation of dynamics is too costly or some parameters in dynamics are difficult to compute or measure (1, 2). Therefore, in order to implement the computed-torque control schemes, a parameter identification method for robots must be developed. This motivates the research in this paper. Since the successful applications of fuzzy model in identification of complex systems (3, 4), this paper will use fuzzy model for parameter identification of direct drive robots. However, considering that the dynamics of robots are usually highly nonlinear and coupled, the parameters in the dynamics will be estimated in a sequential and item-byitem manner, that is the parameters in the inertia term and the gravity term will be estimated first. Then with the estimated inertia term and gravity term, the parameters in the Coriolis/centripetal term are estimated. During the identification procedure, it will be assumed that the joint velocities, joint accelerations, and joint torques can be measured without error. 289
290
Chia-Ju Wu and Huang Ching-Huo
The estimation of the parameters in the friction term will not be discussed in this paper. This is because the direct drive robots have no gearing such that the joint friction is very small and can be neglected (5-8). If a robot is characterized by high gear ratios, like most of the commercial robots do, then the joint friction must be taken into account. In addition to eliminating joint friction, direct drive robots also allow the control of joint torques ( ~ 8 ) . This is required in many of the proposed controller though commercial robots usually do not have this feature. The proposed identification method is a two-stage approach. In the first stage, the input-output data used to estimate parameters in the inertia term, the Coriolis/centripetal term, and the gravity term are generated in a recursive manner. Then in the second stage, each parameter is estimated individually. The strategy is to represent the input-output relation of each parameter by fuzzy rules and the output of each rule as a linear combination of input variables (3, 4). The coefficients associated with each rule is then determined by the least squares method (9) to minimize the squared output errors, which means the differences between the output data of the original system and those of the fuzzy model. In order to determine the number and the ranges of membership functions that minimize the sum of squared errors, the breadth-first method (9) and the simulated annealing method (10, 11) are used, respectively. The former is a well-known searching method in the field of data structure and the latter is a statistical training method that has been used in neural networks for years. After estimating the parameters in dynamics of direct drive robots, the control schemes that are developed based on the computed-torque method will be applicable to the robots. A two-link planar robot will be used in the simulation examples for illustration.
IL Review of the Computed-Torque Method In general, the dynamics of an n-joint drive robot can be written as (1, 2)
g(o)O+ V(O, O) + G(O) = z
(1)
where 0 = [0~ 0z ... 0,]T~ ~" is the joint position vector, M(O)~ ~"×" is the inertia matrix, V(O,O)E ~" is the Coriolis/centripetal vector, G(O)~ ~" is the gravity vector, and z s R" is the joint torque/force vector. If every term in the above equation is known fairly accurately, then a family of control schemes called the computed-torque method (1, 2) can be applied to solve the trajectory-tracking problem of the above robot. For example, the PID computedtorque control law is given as (1, 2)
z = M(O) (Od+ Kv~ + Kpe + Kite" dt) + V(O, O) + G(O)
(2)
with the tracking error defined as e = 0d- 0
(3)
where K,,, Kp, and Ki are diagonal constant matrices, and 0d is the desired trajectory. Substituting Eq. (2) into Eq. (1) yields
(+K,~+Kpe+Ki~e'dt = 0.
(4)
Fuzzy Parameter Identification of Direct Drive Robots
291
From the above equation, one can find that the tracking error will approach zero as long as K,,, Kp and K~ are chosen such that K,, = k,,. L×,
(5)
Kp = kt,- I. xn
(6)
K, = k,./~ ×.
(7)
and
k, < k,
.k~.
(8)
The above analysis shows that the computed-torque controller provides a firm mathematical frame for the solution of trajectory-tracking problem of robots. However, as mentioned above, this family of control schemes is applicable only when the parameters in the dynamics are known fairly accurately. Therefore, these parameters should be estimated.
IlL Generation of Data for Identification Before presenting the procedure of parameter estimation, the input-output data for identification must be generated first. In this section, it will be shown how to generate identification data for parameters in the inertia term, the Coriolis/centripetal term, and the gravity term in a recursive manner. The first step is to write Eq. (1) in another form
(1, 2) M(O)O + Vco,(O)[O0] +
Vcen(O)[O2] -[- G(O)
= ~
(9)
where V,.or(O) is a matrix of Coriolis coefficients with dimension n x n(n-1)/2, [00] is an n ( n - 1)/2 x 1 vector of joint velocity products given by
0,03
[001=[0102
...
0,_,0,] T
(10)
V,.e,(0) is an n x n matrix of centripetal coefficients, and [02] is an n x 1 given by
[02]=[0, 2
07
...
0.2]L
(11)
To estimate parameters in M(O) and G(O), one can set 0 = Oinit
(12)
and
0=0
(13)
where Oi,i, is an arbitrarily chosen joint position vector. Then apply (n 2+ n + 2)/2 different joint torque vectors, r(l), ~(2) . . . . . ~((n 2+ n + 2)/2) to the system (the choice of the number (n2+n+2)/2 will become clear later). In this manner vectors [00] and [02] will be zero and the following set of simultaneous equations will be obtained
M(O,,,)O(i) + G(Oi,u) = r(/)
for i = 1,2 . . . . . (n 2 + n + 2)/2.
(14)
Chia-Ju Wu and Huan9 Ching-Huo
292
Since M(O~,,.) and G(0~..) are the same for the (n2+n+2)/2 equations in (14), by subtracting one equation from the next one, one has
M(Og,.t)(O(i+l)-O(i)) = z ( i + l ) - z ( i )
f o r / = 1,2 . . . . . (n2+n)/2.
(15)
Since the inertia matrix is symmetric (1, 2), there are (n2+n)/2 parameters to be determined. Therefore, if the joint accelerations and the joint torques can be measured without error, there will be (n2+ n)/2 unknown variables in the (n 2+n)/2 equations of (15). The value of M(Og..) can then be obtained by solving the set of simultaneous equations in (15). After obtaining the value of M(Oi..), by substituting it into any equation in (14), the value of G(O~,i,)is obtained. By changing the value of 0~,, in Eq. (12) recursively while keeping the value of 0 be zero, the inertia matrix M(O) and the gravity term G(O) for different value of 0 can be obtained. Therefore, a sufficient number of identification data for M(O) and G(O) can be generated without difficulty. To generate identification data for Vc~,(O),the first step is to fix joint 2 through joint n such that 0=[0,
0
0
...
O]T.
(16)
Substituting Eq. (16) into Eq. (9) yields Vcen(O)[O 2
0
O]r=z--M(O)O-G(O).
..,
(17)
If l)1,/), and z can be measured without error, and M(O) and G(O) have been estimated (the identification procedure will be shown in the next section), then the parameters in the first column of Vc~,(O)will be determined. Similarly, if 0=[0
...
0
0,
0
...
01T
(18)
then the parameters in the i-th column of Vce,(O)can be obtained. In this column-bycolumn manner, different identification data for Vce,(O)are generated. The procedure to generate identification data for V~o,(O)is very similar to that of Vce,(O). The only difference is that two elements in 0 are not equal to zero at one time. For example, if joint 3 through joint n are fixed, then one has 0=[0,
02
0
...
0] T.
(19)
Substituting Eq. (19) into Eq. (9) yields Vcor(O)[O 1
02
0
...
O]V=r-M(O)O-G(O)-V~,(O)[O~
O2 0 ...
0] r.
(20) If 01, 0z, 0, and z can be measured without error, and M(O), G(O), and V~e,(0) have been estimated, then the parameters in the first column of V~o~(0)will be obtained. In this column-by-column manner, different identification data for V~o~(O)are generated.
IV. Parameter Identification After generating identification data for M(O), G(O), Vce.(O), and Nor(O), a fuzzy model will be proposed to estimate parameters in these four terms. Since all these four
293
F u z z y P a r a m e t e r Identification o f D i r e c t Drive R o b o t s
terms are functions of 0, the fuzzy input variables will be 01, 02. . . . . 0,. The parameter to be estimated will be denoted by y, which could be any parameter in M(O), G(O), Vce,(O), and V~or(0). The fuzzy model to estimate the parameter is assumed to be in the following form (3, 4) Rule 1:If01 isA] . . . . . and 0, is A], theny = b ~ + b l "O1 + ' " + b ~ , Rulek:If
"0,.
01 isA] . . . . . and 0, isA,k , t h e n y = b 0 +k b l ' 0k n + . . . + b ~ - 0 n .
where A~, 1 ~< i ~< n, 1 ~
i ,
0,)
i=l
Y-
k
(21)
Y,(A'I(OI) ^ ' " ^
A'n(On))
i--1
where i (W~(O~) A " " A A .i( O . ) ) = m m" ( A ~ i( O ~ ) , A 2 i( 0 2 ) . . . . ,A.(O.)).
(22)
Let f~ be fl, = k
" .4;(0,) ^...^
Z (A',(OI)
i A.(O.)
(23)
i
A.(O.))
^...^
i=1
Then k
Y = ~ fi "(bid "~-btl ° Ol -~-'"--~ hi- On).
(24)
i=l
When m pairs of identification data Ol(j), 02(j), . . ., O,(j) --* y ( j ) ( j = 1,2 . . . . . m) are given, one can obtain (25)
Y = QB
where Y = [y(1) 8=[bo
[
f,
Q=
(1)
Lfll(m)
]~1 ( 1 ) 0 1
~ (1)
fll(m)Ol(m)
...
b .~
...
fl
...
b o~
...
y(2) ...
(1)0,(1)
fll(m)On(m)
b .~
y(rn)] T ..
bko
(26) ...
b~] T
(27)
...
pig(l)
...
fg(m) fg(m)O,(m) ...
fk(1)O,(1)
...
(28)
294
Chia-Ju Wu and Huan# Chin#-Huo fli(j) :
k
A'I (01 ( j ) ) ^ ' "
^ A~.(O.(j))
(29)
~, (Ai~(O~(j)) ^ . . . A A~(O.(j)) i=l
The coefficients vector B associated with each rule can then be obtained by the least squares method using Eq. (25), i.e. B = Q+Y
(30)
where Q+ is the pseudo-inverse of Q (13). In this manner, the sum of squared output errors (y_QB)r(y_QB)
= ( y _ Q Q + y ) r ( y _ Q Q + y)
(31)
will be minimized.
V. Number and Ranges of Membership Functions It should be noted that the procedure described in the above section can be performed only after the number of membership functions, the number of fuzzy rules, and the ranges of membership functions have been chosen. Since different choices will result in different estimated parameters, it will be shown how to determine these values such that the sum of squared output errors is globally minimal. In general, there seems no theoretical approach to determine the number of membership functions associated with each input variable. Therefore, a directed graph (9) will be created and the breadth-first method (9) will be used. Each node in the directed graph will denote a set of numbers, in which each one corresponds to the number of membership functions of an input variable. The direct graph of a fuzzy model with two input variables is shown in Fig. 1. The searching strategy is to make the number of membership functions as small as possible provided that the sum of squared output errors is within a predetermined value. If the number of membership functions associated with each input variable is determined, then the number of fuzzy rules will be the product of these numbers. For example, if there are two input variables and each input variable has two membership functions associated with it, then there will be four fuzzy rules. Since the number of fuzzy rules can be determined via the number of membership functions associated with each input variable, its determination will not be discussed in this paper. The membership functions in this paper are assumed to have linear shapes (3, 4). Therefore, for each membership function, there are two range values, which give the greatest grade 1 and the least grade 0, to be determined. By defining the performance index as the sum of squared output errors, a learning mechanism is presented to find the optimal range of values to minimize the performance index. This mechanism is developed based on a concept similar to the 'simulated annealing' method (10, 11) in the neural networks theory. The first step is to select a set of nominal range values. Then large random adjustments of range values are made. Only those adjustments that reduce the performance index or have the probability greater than a predetermined
Fuzzy Parameter Identification of Direct Drive Robots
295
(1,1)
/11-.... (I 2)
(2,1)
(3,1)
(2,2)
(I ,3)
f
(4,1)
(3,2)
(2,3)
(1,4)
FIG. 1. The direct graph of a fuzzy model with two input variables. In each node, the two numbers correspond to the numbers of membership functions of the first and the second input variables, respectively. value will be retained. The average size of the adjustment is then gradually reduced and optimal range values will eventually be reached. Since the learning time is always a finite value in practice, the learning mechanism will stop as the range values converge to fixed values or the learning time is longer than a predetermined value. In the latter case, it is possible that the performance index will be trapped in a local minimum. Therefore, in the following algorithm, the range values and their corresponding performance index after each random adjustment will be stored in a set. The range values which minimize the performance index is found by comparing all the elements in the set after the learning process. The details of the learning mechanism are described as follows.
Algorithm A (Learning mechanism for range values) Step Step
1. Let S be the set in which each element represents a set of range values and their corresponding performance index. Initially, this is an empty set. 2. Choose the time increment value At, the final learning time t~na/, and the probability comparison value 7.
296 Step Step Step Step Step Step Step Step
Step
Step Step
Chia-Ju W u and H u a n 9 C h i n g - H u o
Choose the nominal range values for each membership function. 4. According to procedure described in Section IV, find the performance index corresponding to the range values in Step 3. . Put the range values in Step 3 and their corresponding performance index into the set S as an element. 6. Set time t --- 0. 7. t = t-t-At. 8. Make random change in magnitude of each range value and add this magnitude to each corresponding range value. 9. Find the performance index corresponding to the range values in Step 8. 10. If the performance index is reduced, then: (1) retain the change in Step 8; (2) put the range values and their corresponding performance index into the set S as an element; and (3) go to Step 12; otherwise continue. 11. Calculate the probabilities for each magnitude change in Step 8. If the minimal value of these probabilities is greater than y, then: (1) retain the change in Step 8; and (2) put the range values and their corresponding performance index into the set S as an element; otherwise, return to their previous range values. 12. If t >>, t~,a~, then continue; otherwise go to Step 7. 13. Find the optimal range values which minimize the performance index by comparing all the elements in S. .
In Step 11 of the above algorithm, the probabilities are calculated by the Cauchy distribution as follows T
P(x) =
T 2+ x 2
(32)
where P ( x ) is the probability of a random change x, and T=
1 l+t"
(33)
The reason for choosing the Cauchy distribution in Eq. (32) is to lead the learning process to a linear convergence rate. A detailed discussion can be found in (14). Substituting Eq. (33) into Eq. (32) yields P(x) -
(1 +t) . 1 +x2(1 + t ) 2
(34)
For a random change of range value to be retained, if it does not reduce the performance index, the value of P ( x ) must be greater than the predetermined value y, i.e.
(1 +t) 1 +x~(1 + 8 ) 2
. . ~.
(35)
Solving for x in Eq. (35) yields
(1+8-~)
x2 < - -
~(l+t) 2
(36)
Fuzzy Parameter Identification of Direct Drive Robots
297
From this equation one can see that as time t increases, the values of x must be decreased to satisfy Eq. (36). Therefore, the probabilities of large random change are reduced gradually. In Step 8 of the above algorithm, there are many ways to determine the random change. For example, one can select a random number P(x) from a uniform distribution between zero and one. Then according to the gamma distribution
P(x) = e x p ( - x)
(37)
where P(x) is the probability of a random change with magnitude x, the magnitude x can be obtained as x = -ln(P(x)).
(38)
The corresponding random change is determined by adding randomly a positive or negative sign to the magnitude. Before presenting simulation examples to illustrate the proposed method, the process of parameter identification is summarized in the following algorithm.
Algorithm B (Parameter identification) Step Step Step
Step Step
Step Step Step Step Step
1. Generate identification data for M(O), G(O) according to the procedure described in Section III. 2. For every parameter in M(O) and G(O), perform the identification procedure in Steps 3-5. 3. Create a directed graph and start the identification procedure from the root node of the directed graph. Every node in the directed graph will denote a set of numbers, in which each one corresponds to the number of membership functions of an input variable. The number of fuzzy rules is also determined as the product of this set of numbers. 4. The optimal range values for each membership function are determined by the procedure described in Algorithm A. 5. If the sum of squared output errors obtained in Step 4 is less than a predetermined value or the number of fuzzy rules exceeds a predetermined number, then the identification procedure is completed. Otherwise, go to the next node in the directed graph according to the searching strategy of the breadth-first method. For the new node, go to Step 4 to find the optimal range values. 6. Generate identification data for V,.e,(O) according to the procedure described in Section III. 7. For every parameter in Vce,(O),repeat the identification procedure in Steps 3-5. 8. Generate identification data for V,or(O) according to the procedure described in Section III. 9. For every parameter in Vcor(O),repeat the identification procedure in Steps 3-5. 10. End.
After parameters in M(O), G(O), Vce,(O), and V,.o,(O) have all been estimated, the
Chia-Ju Wu and Huan9 Ching-Huo
298 Y-axis
nlz
,¢
/ /
¢ e2
/
/
/ /
O, X-axis FIG. 2. A two-link planar robot with rotary joints, where 0j and 02 denote the joint positions and
l~ and/2 denote the lengths of links. The dynamics equations are given in Eqs (39) and (40). computed-torque method can then be applied to make the robot move along a desired trajectory.
VI. Simulation Examples In the following simulation examples, the two-link planar robot shown in Fig. 2 will be used for illustration. The dynamics of this robot, which are assumed to be unknown, are given as follows (1, 2) Zl = [(ml
+m2)l~+m21~+2m21112cos 02]01 +(m21~+m:lll2 cos 02)02
--m21112(20,02+0~)sin 02 +
(m,
+m2)gll cos01 +m2g12cos(01 +02)
(39)
r2 = (m21z,+m21112cos 02)01 +m21202+m2111202 sin 02 +mzffl2cos(01 71-02). (4o) Comparing Eqs (39) and (40) with Eq. (9), one can find that the parameters to be estimated are
M,2]=[(m, +m2)l~+m21Z~+ 2m21112cos02 m21~+m21l12cos02] LM21 M22] m212+m21112COS02 m212
M=~M11
(41)
Fuzzy Parameter Identification of Direct Drive Robots
299 (42)
L Vco,2J
L
Vcen =[Vcenl 1 Vcenl2] =
0
0
--m2lll2sinOz I
LVcen21 Vcen22J [mzlll2 sin02 G = G2 =
[m+m2'cosO+m2'2osO+02,] m2ffl 2 cos(01 _t_02)
(43)
0
.
(44)
Through the simulation, the values of ml, m2, 1~, 12, and 9 are chosen to be 1 kg, 1 kg, 1 m, 1 m and 9.8 m/s 2, respectively. The ranges of 0t and 02 are assumed to be - r t ~< 0~ ~< rt and - r e ~< 02 ~< re.
Example 1 This example will show how to estimate the parameters in M(0), G(0), V,,n(O), and V,.or(O) according to the procedure described in Algorithm B. In applying Algorithm A to find the optimal range values, the values of (r, At, and ? are chosen to be 5 0 , 0 . 1 , and 0.5 s, respectively. The nominal range values for each membership function are chosen such that the membership functions are evenly distributed in the range [ - ~z, 7r] as shown in Fig. 3.
(a) Estimation of Mll and M12 During the procedure of generating identification data, one can find that: (1) M~ and M12 are functions of 02; (2) M21 is equal to Mj2 (since M is symmetric); and (3) Membership Functions
/
n u n
0 FIG. 3. Plots of the nominal membership functions when applying Algorithm A to find the optimal range values. If there are k membership functions associated with 01(02), then the membership functions have grade 0 at points 01(02) = - ~ + ( i - l ) ' 2 n / k , i = 1, 2, ..., k, and grade 1 at points 07(02) = - r~+ i" 2rc/k, i = 1, 2 . . . . . k.
Ch&-Ju Wu and Huan 9 Ching-Huo
300
TABLE I
Minimal values of Y,)°=°l(Error,) 2, maxi(Errori) 2 and min,(Errori)2 for k = 1,2 . . . . . 7, where k is the number of membership functions of 02 and Error, denotes [MH(i) --3;Ill(,)]for simplicity of notation k
E) °° (Error;) 2
maxi(Errorl) 2
1 2 3 4 5 6 7
129.2268 3.8635 1.1785 0.6543 0.2772 0.0956 0.0043
5.7896 1 0.2336 0.0284 0.0240 0.0072 1.7487" 10 -4
mini(Errori) 2 1.2693" 10 -6
2.3383" 10 -8 4 . 0 3 3 9 " 10 -6
1.7912" 10 -8 1.9075" 10 -7 4 . 3 5 2 9 " 10 -7 3.9178" 10 -9
M22 is a constant item. Therefore, only two parameters, namely M n and M~2, need to be estimated in the inertia matrix. Since MH is a function of 02, the fuzzy model to estimate M l l will be
Rule 1:
I f 02 is AI,
1 l then 3~ll = bo+bl "02
Rule k:
I f 02 is A k,
then AI11 = bok + b lk • 02
where _~t~l denotes the estimate of MI1 and k denotes both the number of fuzzy rules and the number of membership functions of 02. Suppose that 100 pairs of identification data for M~I and M~2 (02(/) ~ M . ( / ) , 02(/) -~ M12(/), i = 1,2 . . . . . 100, 02(0 = - n + [ ( i - 1)/99]" 2n) have been generated according to the procedure described in Section III. Then following the procedure in Algorithm B, the minimal values of Z~°° [ M ~ ( / ) - - M I I ( i ) ] 2 for k -- 1, 2 . . . . . 7, are found and shown in Table I. The plot of ( M l l - - M n ) for the case o f k = 7 is shown in Fig. 4. Following the same procedure, the minimal values of Z) °° [M~2(/)- ~ t l 2 ( i ) ] 2 for k = l, 2 . . . . . 7, are found and shown in Table II. The plot of (M~2- M~2) for the case of k = 7 is shown in Fig. 5.
(b) Estimation of Vcen21, V~12, and V~orl During the procedure of generating identification data, one can find that: (1) Vcen21,
Vcenl2, and V~orl are functions of 02, and (2) Vcenl 1 = Vcen2 2 = Vcor2 = O, Therefore, three parameters, namely V~e.2~, V~e.~2, and Vcor~, need to be estimated. Since the three parameters are functions of 02, their fuzzy models will be the same as that of Mll. Suppose that 100 pairs of identification data for Vce.21, Vcen~2, and V~or~ (02(/) Vcen21(/), 02(/) ~ Vcenl2(/), 02(/') ~ Veor~(/),i = 1,2 . . . . . 100,02(/) = --7r+[(i--1)/99] • 2r0 have been generated. Then following the procedure in Algorithm B, the minimal values of Z] °° [Vce.21(i)--l~c~.2~(/)] 2, E) °° [Vc~,,~z(/)--~,,~2(/)] 2, and E] °° [V~or~(i)--12~or~(/)] 2 for k = 1, 2 . . . . . 7 are found and shown in Tables III-V, respectively. The plots of (V~.2~ - 12~.2~), (V~e.,2- 12~.,2), and (V~or~- 12co.~) for the case of k = 7 are shown in Figs 6-8, respectively.
Fuzzy Parameter Identification of Direct Drive Robots
301
0.015] 0.01
0.005
o
i
-o.oo5
-0.01t -0.015 -4
-3
-2
2
02 FIG. 4. Plot of (Mjt--h4j~) with 7 membership functions associated with 02.
TABLE II
Minimal values of E~°°(Errori)2, maxi(Errorl) 2 and mini(Error32 for k = 1 , 2 , . . . , 7, where k is the number of membership functions of 02 and Errori denotes [Mi2(i)--~tl2(i)] for simplicity of notation k
E]°=° (Errori) 2
1 2 3 4 5 6 7
19.5107 0.6165 0.1831 0.1594 0.0839 0.0348 0.0013
maxi(Errori) 2 1.4516 0.0334 0.0233 0.0141 0.0091 0.0019 7.1864" 10 -5
mini(Errori)2 0 0 5.0728" 10 7 1.1613" ]0 -9
2.3873" 1O-11 1.1262" 10 - 7 3.1633" 10 1o
(c) Estimation o f GI and G2 D u r i n g the p r o c e d u r e o f g e n e r a t i n g identification data, one can find t h a t G1 a n d G2 are functions o f 01 a n d 02. Therefore, the fuzzy m o d e l o f G1 will be in the following form
302
Chia-Ju W u a n d H u a n 9 Chin#-Huo X 10 .3
8
6
4
2
I
=-
0
-2
-4
-6
-8
-10
-4
-3
-2
-1
0
1
2
3
O2 FIG. 5. Plot of (M~2-~1'I12) with 7 membership functions associated with 0a.
TABLEII1 Minimal values of El°=° (Errori) z, maxi(Errorl) 2 and mini(Errori) 2for k = 1,2 . . . . . 7, where k is the number of membership functions of 02 and Errori denotes [ Vc.e,Zl ( i) -- 12~e,21(i)] for simplicity of notation k
E) °° (Errorl) 2
1 2 3 4 5 6 7
15.8913 3.2934 0.4131 0.1861 0.0481 0.0137 0.0016
Rule l:
Rule (kl" k2):
If0~isA11
If01isA~ l
and
and
maxi(Errori)2 1.1337 0.2923 0.0495 0.0128 0.0030 0.0012 2.5404.10 -4
02 isA21,
02 isA~ z,
mini(Errori)2 0 0 6.3683" 10 7 4.3544.10 -8 1.5603.10 -7 1.8888.10 -9 0
then G1 = bo' + b l ' 0 ,
+b~'02
then Gl = bokl"k2+blkl'k2 . 01 .k_b~l.k2. 02
F u z z y P a r a m e t e r Identification o f Direct Drive Robots
303
TABLE IV
Minimal values of £ ] °=° (Err or i) 2, maxi( Error i) 2 and mini( Error i)2for k = l, 2 . . . . . 7, where k is the number of membership functions of 02 and Errori denotes [Vc,nl2(i) - I?~enl2(O]for simplicity of notation k
~,~°°(Errori) 2
1 2 3 4 5 6 7
14.1730 1.2100 0.4013 0.0863 0.0225 0.0193 0.0027
maxi(Errori) 2 0.9118 0.1847 0.0337 0.0089 0.0017 0.0023 3.2705" 10 4
mini(Errori) 2 0 0 7.4998" 1.3867" 1.2030" 5.2056" 2.3567"
l0 -7 10 7 10 s 10 7 10 ~
TABLE V
Minimal values of £ ]°° (Error i) 2, maxi( Error i)2 and mini( Error y for k = 1,2 . . . . . 7, where k is the number of membership functions of 02 and Errorl denotes [ V,,o,~(i) - V,.,,rl(i)1 for simplicity of notation k l 2 3 4 5 6 7
£]°_°(Errori) 2 41.0945 22.0589 4.2394 1.6625 0.3346 0.0880 0.0015
maxi(Errori) 2 2.5990 2.4096 0.3410 0.2254 0.0421 0.0069 1.2218" 10 -4
mini(Errori) 2 0 0 3.0397. l0 5 5.1873" 10 -7 3.6509" 10 -9
3.6089" 10 -7 3.7890" 10 -9
where 6~1 d e n o t e s the e s t i m a t e o f Gi, k l d e n o t e s the n u m b e r o f m e m b e r s h i p functions o f 01, a n d k2 d e n o t e s the n u m b e r o f m e m b e r s h i p functions o f 02. S u p p o s e t h a t 400 p a i r s o f identification d a t a for G1 a n d G2 ((0~(0, 02(j)) --* Gl(i x j ) , (Ol(i), 02(j)) ~ G2(i × j ) , i = 1, 2 . . . . . 2 0 , j = 1, 2 . . . . . 20, Ol(i) = - - n + [(i-- 1)/19]" 2n, 02(j) = - g - q - [ ( j - - 1 ) / 1 9 ] ' 2 n ) have been g e n e r a t e d p r i o r to a p p l y i n g A l g o r i t h m B. T h e n following the p r o c e d u r e in A l g o r i t h m B, the m i n i m a l values o f ~400 i= 1 [(G1 (i) -- (~1 (i))/y]2 for different c o m b i n a t i o n s o f k l a n d k2 are f o u n d a n d s h o w n in T a b l e VI. T h e p l o t o f ( G 1 - t~0/g for the case o f k l = 8 a n d k2 = 8 is s h o w n in Fig. 9. F o l l o w i n g the s a m e p r o c e d u r e , the m i n i m a l values o f w400 ,-,i=l [(G2 (i) - G2 (0)/912 for different c o m b i n a t i o n s o f k l a n d k2 are f o u n d a n d shown in T a b l e VII. T h e p l o t o f ( G 2 - G 2 ) / g for the case o f k l = 8 a n d k2 = 8 is shown in Fig. 10.
Example 2 This e x a m p l e shows h o w to a p p l y the P I D c o m p u t e d - t o r q u e c o n t r o l l e r b a s e d on the e s t i m a t e d d y n a m i c s to m a k e the r o b o t m o v e a l o n g the desired t r a j e c t o r y (1)
Chia-Ju Wu and Huang Ching-Huo
304 0.02
0.015
0.01
~d 0.005
I ~o
-0.005
-0.01 -4
~ -3
~ -2
~ -1
2
0
3
02
FIG. 6. Plot of (Vcen2I -- Veen21) with 7 membership functions associated with 02.
01a = 0.1 • sin(rtt)
(45)
02a = 0.1" cos(rtt).
(46)
As described in Section II, the PID computed-torque control law is given as
z -- M(O)(Oa+Kv~+Kpe+K~Se'dt)+ 12,.or(O)[O0]+Vcen(O)[Oz]--~a(o)
(47)
where/14(0), 12cor(0), and ~rcen(O ) a r e estimated with 7 membership functions associated with 02, and (~(0) is estimated with 8 membership functions associated with both 01 and 02. If the initial condition of [01 02] is chosen to be [0.1 0] and Kv, Kp, and K~ are chosen as
K,,=
Kp
and
=
50
oI
625
Fuzzy Parameter Identification of Direct Drive Robots
305
°°2 /
0.0151 0.01
/
0.005
I
-0.005 -0.01 /
-0.01 5 I
I
0 02 -4
t -3
L -2
, -1
~ 0
~ 1
r 2
3
02 FIG. 7. Plot of (Vce.12- I?,.e.12)with 7 membership functions associated with 02.
K; =
1200
then the plots of (O,d--0,) and (02d--02) are as shown in Figs 11 and 12, respectively.
VII. Conclusions
This work showed how to generate identification data in a recursive manner and how to apply a fuzzy model to estimate parameters in dynamics of a direct drive robot. From the simulation results shown in Tables I-VII, one can find that AIr(0), l?ce,(0), 12cot(0), and 6~(0) approximate the original parameters well. It is also expected that one can increase the number of training pairs to reduce the differences between the estimates and the original parameters. Based on the estimated dynamics, it was shown in Example 2 that the PID computed-torque controller gave a good result in trajectory-tracking problem of robots. Though the applications of the proposed method in the simulation examples are successful, it should be pointed out that heavy computation burden is inevitable in the identification procedure. Most of the computation time will spend in generating identification data and in searching for the optimal range values of each membership
Chia-Ju Wu and Huan9 Chiny-Huo
306 0.01
0.005
-0.005
-0.01
-0.015
' -3
-4
P -2
'
' 0
-1
' 1
' 2
3
O2
FIG. 8. Plot of
(V,.o,j-
I2,-o,0 with 7 membership functions associated with 02.
TABLE VI
Minimal values of 5~4°_°(Errori)2, maxi(Errori) 2 and mini(Errori) 2 for different combinations of kl and k2, where kl and k2 are the numbers of membership functions of O~ and 02, respectively, and Errori denotes [GI(i) - Gl(i)]/g for simplicity of notation kl
k2
Z]°- ° (Errori) 2
maxi(Errorl) 2
mini(Errori) 2
1 2 3 4 5 6 7 8
1 2 3 4 5 6 7 8
481.3954 165.8560 114.4186 40.0132 3.5234 2.3358 0.7215 0.0724
8.9184 4.4778 5.3172 1.6219 0.1778 0.0946 0.0263 0.0017
3.8526.10 5 1.2495 • l0 -5 4.1347.10 -8 1.0243" 10 -6
1.9926" 10 8 7.0697" l0 -9 2.6583" 10 -11 1.3458" 10 -28
function. Unless the c o m p u t a t i o n complexity o f the identification procedure is reduced, the p r o p o s e d m e t h o d can only be used in a off-line m a n n e r . Therefore, how to reduce the c o m p u t a t i o n complexity of the p r o p o s e d m e t h o d is suggested as a topic for further
Fuzzy Parameter Identification o f Direct Drive Robots
307
(Gt-al)/g 0.040.02o -0.02-0.04-
2
4
-4
-4
01
FIG. 9. Plot of (G~- (~1)/# with 8 membership functions associated with both 0~ and 02.
TABLE VII
Minimal values of Y?°_°(Errori) 2, maxi(Errori): and mini(Errori) 2 for difJerent combinations of kl and k2, where ks and k2 are the numbers of membership functions of 01 and 02, respectively, and Error, denotes [G2(i) - G2(i)]/g for simplicity of notation kl
k2
1 2 3 4 5 6 7 8
1 2 3 4 5 6 7 8
E]°°(Errori) z 162.6653 77.2357 19.9503 7.3790 1.4301 0.2837 0.0987 0.0329
maxi(Errori) z 1.6900 2.1638 1.6308 0.2935 0.0686 0.0217 0.0016 0.0007
min,(Errori) 2 1.5618 2.0618 1.5097 1.3840 3.2439 1.1431 7.0626" 3.8807"
10 -5 10 7 10 -7 10 -6 10 v 10 9 10 lo 10 9
research. I n a d d i t i o n , since it is a s s u m e d that the j o i n t velocities, j o i n t accelerations, a n d j o i n t torques c a n be m e a s u r e d w i t h o u t error, the m e a s u r e m e n t of j o i n t variables is suggested as a n o t h e r research topic.
308
Chia-Ju Wu and Huan9 Chiny-Huo
0.030.02-
0.01
-
A
0-0.01
-
-0.02-
-0.03> 4
o
/ -4
-4 FIG. 10. Plot of
-2 01
(G2-G2)/gwith 8 membership functions associated with both
01 and 02.
0.1
0.05
-0.05
j i 0 5
~ 1
15 1.
I 2
I 2.5
I 3
3.5
I 4
Time(sec.)
FIG. 11. Plot of (0~-- 0j), where 0~d is given in Eq. (45).
t
4.5
Fuzzy Parameter Identification o f Direct Drive Robots
309
0.02,
-0.02
I
-0.04
,~a~ "~1
-0.06
-0.08
r
"0"10
I
I
I
I
I
[
I
I
I
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Time(see.) FIG. 12. Plot of (02d--02), where 02d is given in Eq. (46).
A cknowledgements This work was supported by the National Science Council, Taiwan, Republic of China, under grant No. NSC83-0422-E-224-036.
References (1) F. L. Lewis, C. T. Abdallah and D. M. Dawson, "Control of Robot Manipulators", Macmillan, New York, NY, 1993. (2) J. J. Craig, "Introduction to Robotics: Mechanics and Control", 2nd Edn, AddisonWesley, Reading, MA, 1988. (3) M. Sugeno and K. Tanaka, "Successive identification of a fuzzy model and its applications to prediction of a complex system", Fuzzy Sets Syst., Vol. 42, pp. 315-334, 1991. (4) T. Takagi and M. Sugeno, "Fuzzy identification of systems and its applications to modeling and control", IEEE Trans. Syst. Man Cybernet., Vol. SMC-15, No. 1, pp. 116-132, 1985. (5) H. Asada, T. Kanade and I. Takeyama, "Control of a direct-drive arm", ASME J. Dynam. Syst., Meas., Contr., Vol. 105, pp. 136-142, 1983. (6) H. Asada and K. Youcef-Toumi, "Analysis and design of a direct-drive arm with a fivebar-link parallel drive mechanism", ASME J. Dynam. Syst. Meas. Contr., Vol. 106, pp. 225-230, 1984. (7) H. Asada and K. Youcef-Toumi, "Direct-Drive Robots: Theory and Practice", MIT Press, Cambridge, MA, 1987. (8) K. Youcef-Toumi and A. T. Y. Kuo, "High-speed trajectory control of a direct-drive manipulator", IEEE Trans. Robot. Automat., Vol. 9, pp. 102-108, 1993.
310
Chia-Ju Wu and Huan9 Ching-Huo
(9) R. L. Kruse, B. P. Leung and C. L. Tondo, "Data Structures and Program Design in C", Prentice-Hall, Englewood Cliffs, N J, 1991. (10) B. Kosko, "Neural Networks and Fuzzy Systems", Prentice-Hall, Englewood Cliffs, N J, 1992. (11) P. D. Wassermann, "Neural Computing: Theory and Practice", Vann Nostrand, Reinhold, NY, 1989. (12) C. C. Lee, "Fuzzy logic in control systems: fuzzy logic controller, part II", IEEE Trans. Syst. Man Cybernet., Vol. 20, pp. 419-435, 1990. (13) T. N. E. Greville, "The pseudoinverse of a rectangular or singular matrix and its applications to the solutions of systems of linear equations", SIAM Rev., Vol. 1, pp. 38-43, January, 1959. (14) H. Szu and R. Hartley, "Fast simulated annealing", Phys. Lett., Vol. 1222 (3,4), pp. 157162, 1987.
~ ) Pergamon
S0016--0032(96)00026-9
Copyright © 1996 The Franklin Institute Published by Elsevier Science Ltd Printed in Great Britain 0016-0032/96 $15.00+0.00
Fuzzy ParameterIdentification of Direct Drive Robots by CHIA-JU
WU
and
HUANG
CHING-HUO
Department o f Electrical Engineering, National Yunlin Institute o f Technology, Touliu, Yunlin, 640 Taiwan (Received 22 August 1994," accepted 31 January 1996)
ABSTRACT: ,4fuzzy identification method is proposed to estimate parameters in dynamics of direct drive robots. In this method, linear membership functions are used and the output of each fuzzy rule is assumed to be the linear combination of input variables. The input-output data for identification are generated in a recursive manner. The number and the ranges of memberships functions are determined by a heuristic searching method and a statistical training method, respectively, in the sense of least squared output errors. Once the parameters in dynamics are estimated, existing control strategies, such as the well-known computed-torque method, can be applied to make the robot follow a desired trajectory. Simulation examples are included to illustrate the validity of the proposed method. Copyright © 1996 Published by Elsevier Science Ltd L Introduction
Research of robots has been prompted in recent years by a wide variety of goals. One of the goals is to develop control algorithms to make a robot move along a desired trajectory. Among the trajectory-tracking algorithms that have been proposed previously (1, 2), the computed-torque method provides a firm mathematical frame and yields a family of easy-to-understand control schemes. These schemes perform well provided that the parameters in dynamics of robots are known fairly accurately. However, in practice, there are uncertainties in the dynamics of robots. This may arise because the exact evaluation of dynamics is too costly or some parameters in dynamics are difficult to compute or measure (1, 2). Therefore, in order to implement the computed-torque control schemes, a parameter identification method for robots must be developed. This motivates the research in this paper. Since the successful applications of fuzzy model in identification of complex systems (3, 4), this paper will use fuzzy model for parameter identification of direct drive robots. However, considering that the dynamics of robots are usually highly nonlinear and coupled, the parameters in the dynamics will be estimated in a sequential and item-byitem manner, that is the parameters in the inertia term and the gravity term will be estimated first. Then with the estimated inertia term and gravity term, the parameters in the Coriolis/centripetal term are estimated. During the identification procedure, it will be assumed that the joint velocities, joint accelerations, and joint torques can be measured without error. 289
290
Chia-Ju Wu and Huang Ching-Huo
The estimation of the parameters in the friction term will not be discussed in this paper. This is because the direct drive robots have no gearing such that the joint friction is very small and can be neglected (5-8). If a robot is characterized by high gear ratios, like most of the commercial robots do, then the joint friction must be taken into account. In addition to eliminating joint friction, direct drive robots also allow the control of joint torques ( ~ 8 ) . This is required in many of the proposed controller though commercial robots usually do not have this feature. The proposed identification method is a two-stage approach. In the first stage, the input-output data used to estimate parameters in the inertia term, the Coriolis/centripetal term, and the gravity term are generated in a recursive manner. Then in the second stage, each parameter is estimated individually. The strategy is to represent the input-output relation of each parameter by fuzzy rules and the output of each rule as a linear combination of input variables (3, 4). The coefficients associated with each rule is then determined by the least squares method (9) to minimize the squared output errors, which means the differences between the output data of the original system and those of the fuzzy model. In order to determine the number and the ranges of membership functions that minimize the sum of squared errors, the breadth-first method (9) and the simulated annealing method (10, 11) are used, respectively. The former is a well-known searching method in the field of data structure and the latter is a statistical training method that has been used in neural networks for years. After estimating the parameters in dynamics of direct drive robots, the control schemes that are developed based on the computed-torque method will be applicable to the robots. A two-link planar robot will be used in the simulation examples for illustration.
IL Review of the Computed-Torque Method In general, the dynamics of an n-joint drive robot can be written as (1, 2)
g(o)O+ V(O, O) + G(O) = z
(1)
where 0 = [0~ 0z ... 0,]T~ ~" is the joint position vector, M(O)~ ~"×" is the inertia matrix, V(O,O)E ~" is the Coriolis/centripetal vector, G(O)~ ~" is the gravity vector, and z s R" is the joint torque/force vector. If every term in the above equation is known fairly accurately, then a family of control schemes called the computed-torque method (1, 2) can be applied to solve the trajectory-tracking problem of the above robot. For example, the PID computedtorque control law is given as (1, 2)
z = M(O) (Od+ Kv~ + Kpe + Kite" dt) + V(O, O) + G(O)
(2)
with the tracking error defined as e = 0d- 0
(3)
where K,,, Kp, and Ki are diagonal constant matrices, and 0d is the desired trajectory. Substituting Eq. (2) into Eq. (1) yields
(+K,~+Kpe+Ki~e'dt = 0.
(4)
Fuzzy Parameter Identification of Direct Drive Robots
291
From the above equation, one can find that the tracking error will approach zero as long as K,,, Kp and K~ are chosen such that K,, = k,,. L×,
(5)
Kp = kt,- I. xn
(6)
K, = k,./~ ×.
(7)
and
k, < k,
.k~.
(8)
The above analysis shows that the computed-torque controller provides a firm mathematical frame for the solution of trajectory-tracking problem of robots. However, as mentioned above, this family of control schemes is applicable only when the parameters in the dynamics are known fairly accurately. Therefore, these parameters should be estimated.
IlL Generation of Data for Identification Before presenting the procedure of parameter estimation, the input-output data for identification must be generated first. In this section, it will be shown how to generate identification data for parameters in the inertia term, the Coriolis/centripetal term, and the gravity term in a recursive manner. The first step is to write Eq. (1) in another form
(1, 2) M(O)O + Vco,(O)[O0] +
Vcen(O)[O2] -[- G(O)
= ~
(9)
where V,.or(O) is a matrix of Coriolis coefficients with dimension n x n(n-1)/2, [00] is an n ( n - 1)/2 x 1 vector of joint velocity products given by
0,03
[001=[0102
...
0,_,0,] T
(10)
V,.e,(0) is an n x n matrix of centripetal coefficients, and [02] is an n x 1 given by
[02]=[0, 2
07
...
0.2]L
(11)
To estimate parameters in M(O) and G(O), one can set 0 = Oinit
(12)
and
0=0
(13)
where Oi,i, is an arbitrarily chosen joint position vector. Then apply (n 2+ n + 2)/2 different joint torque vectors, r(l), ~(2) . . . . . ~((n 2+ n + 2)/2) to the system (the choice of the number (n2+n+2)/2 will become clear later). In this manner vectors [00] and [02] will be zero and the following set of simultaneous equations will be obtained
M(O,,,)O(i) + G(Oi,u) = r(/)
for i = 1,2 . . . . . (n 2 + n + 2)/2.
(14)
Chia-Ju Wu and Huan9 Ching-Huo
292
Since M(O~,,.) and G(0~..) are the same for the (n2+n+2)/2 equations in (14), by subtracting one equation from the next one, one has
M(Og,.t)(O(i+l)-O(i)) = z ( i + l ) - z ( i )
f o r / = 1,2 . . . . . (n2+n)/2.
(15)
Since the inertia matrix is symmetric (1, 2), there are (n2+n)/2 parameters to be determined. Therefore, if the joint accelerations and the joint torques can be measured without error, there will be (n2+ n)/2 unknown variables in the (n 2+n)/2 equations of (15). The value of M(Og..) can then be obtained by solving the set of simultaneous equations in (15). After obtaining the value of M(Oi..), by substituting it into any equation in (14), the value of G(O~,i,)is obtained. By changing the value of 0~,, in Eq. (12) recursively while keeping the value of 0 be zero, the inertia matrix M(O) and the gravity term G(O) for different value of 0 can be obtained. Therefore, a sufficient number of identification data for M(O) and G(O) can be generated without difficulty. To generate identification data for Vc~,(O),the first step is to fix joint 2 through joint n such that 0=[0,
0
0
...
O]T.
(16)
Substituting Eq. (16) into Eq. (9) yields Vcen(O)[O 2
0
O]r=z--M(O)O-G(O).
..,
(17)
If l)1,/), and z can be measured without error, and M(O) and G(O) have been estimated (the identification procedure will be shown in the next section), then the parameters in the first column of Vc~,(O)will be determined. Similarly, if 0=[0
...
0
0,
0
...
01T
(18)
then the parameters in the i-th column of Vce,(O)can be obtained. In this column-bycolumn manner, different identification data for Vce,(O)are generated. The procedure to generate identification data for V~o,(O)is very similar to that of Vce,(O). The only difference is that two elements in 0 are not equal to zero at one time. For example, if joint 3 through joint n are fixed, then one has 0=[0,
02
0
...
0] T.
(19)
Substituting Eq. (19) into Eq. (9) yields Vcor(O)[O 1
02
0
...
O]V=r-M(O)O-G(O)-V~,(O)[O~
O2 0 ...
0] r.
(20) If 01, 0z, 0, and z can be measured without error, and M(O), G(O), and V~e,(0) have been estimated, then the parameters in the first column of V~o~(0)will be obtained. In this column-by-column manner, different identification data for V~o~(O)are generated.
IV. Parameter Identification After generating identification data for M(O), G(O), Vce.(O), and Nor(O), a fuzzy model will be proposed to estimate parameters in these four terms. Since all these four
293
F u z z y P a r a m e t e r Identification o f D i r e c t Drive R o b o t s
terms are functions of 0, the fuzzy input variables will be 01, 02. . . . . 0,. The parameter to be estimated will be denoted by y, which could be any parameter in M(O), G(O), Vce,(O), and V~or(0). The fuzzy model to estimate the parameter is assumed to be in the following form (3, 4) Rule 1:If01 isA] . . . . . and 0, is A], theny = b ~ + b l "O1 + ' " + b ~ , Rulek:If
"0,.
01 isA] . . . . . and 0, isA,k , t h e n y = b 0 +k b l ' 0k n + . . . + b ~ - 0 n .
where A~, 1 ~< i ~< n, 1 ~
i ,
0,)
i=l
Y-
k
(21)
Y,(A'I(OI) ^ ' " ^
A'n(On))
i--1
where i (W~(O~) A " " A A .i( O . ) ) = m m" ( A ~ i( O ~ ) , A 2 i( 0 2 ) . . . . ,A.(O.)).
(22)
Let f~ be fl, = k
" .4;(0,) ^...^
Z (A',(OI)
i A.(O.)
(23)
i
A.(O.))
^...^
i=1
Then k
Y = ~ fi "(bid "~-btl ° Ol -~-'"--~ hi- On).
(24)
i=l
When m pairs of identification data Ol(j), 02(j), . . ., O,(j) --* y ( j ) ( j = 1,2 . . . . . m) are given, one can obtain (25)
Y = QB
where Y = [y(1) 8=[bo
[
f,
Q=
(1)
Lfll(m)
]~1 ( 1 ) 0 1
~ (1)
fll(m)Ol(m)
...
b .~
...
fl
...
b o~
...
y(2) ...
(1)0,(1)
fll(m)On(m)
b .~
y(rn)] T ..
bko
(26) ...
b~] T
(27)
...
pig(l)
...
fg(m) fg(m)O,(m) ...
fk(1)O,(1)
...
(28)
294
Chia-Ju Wu and Huan# Chin#-Huo fli(j) :
k
A'I (01 ( j ) ) ^ ' "
^ A~.(O.(j))
(29)
~, (Ai~(O~(j)) ^ . . . A A~(O.(j)) i=l
The coefficients vector B associated with each rule can then be obtained by the least squares method using Eq. (25), i.e. B = Q+Y
(30)
where Q+ is the pseudo-inverse of Q (13). In this manner, the sum of squared output errors (y_QB)r(y_QB)
= ( y _ Q Q + y ) r ( y _ Q Q + y)
(31)
will be minimized.
V. Number and Ranges of Membership Functions It should be noted that the procedure described in the above section can be performed only after the number of membership functions, the number of fuzzy rules, and the ranges of membership functions have been chosen. Since different choices will result in different estimated parameters, it will be shown how to determine these values such that the sum of squared output errors is globally minimal. In general, there seems no theoretical approach to determine the number of membership functions associated with each input variable. Therefore, a directed graph (9) will be created and the breadth-first method (9) will be used. Each node in the directed graph will denote a set of numbers, in which each one corresponds to the number of membership functions of an input variable. The direct graph of a fuzzy model with two input variables is shown in Fig. 1. The searching strategy is to make the number of membership functions as small as possible provided that the sum of squared output errors is within a predetermined value. If the number of membership functions associated with each input variable is determined, then the number of fuzzy rules will be the product of these numbers. For example, if there are two input variables and each input variable has two membership functions associated with it, then there will be four fuzzy rules. Since the number of fuzzy rules can be determined via the number of membership functions associated with each input variable, its determination will not be discussed in this paper. The membership functions in this paper are assumed to have linear shapes (3, 4). Therefore, for each membership function, there are two range values, which give the greatest grade 1 and the least grade 0, to be determined. By defining the performance index as the sum of squared output errors, a learning mechanism is presented to find the optimal range of values to minimize the performance index. This mechanism is developed based on a concept similar to the 'simulated annealing' method (10, 11) in the neural networks theory. The first step is to select a set of nominal range values. Then large random adjustments of range values are made. Only those adjustments that reduce the performance index or have the probability greater than a predetermined
Fuzzy Parameter Identification of Direct Drive Robots
295
(1,1)
/11-.... (I 2)
(2,1)
(3,1)
(2,2)
(I ,3)
f
(4,1)
(3,2)
(2,3)
(1,4)
FIG. 1. The direct graph of a fuzzy model with two input variables. In each node, the two numbers correspond to the numbers of membership functions of the first and the second input variables, respectively. value will be retained. The average size of the adjustment is then gradually reduced and optimal range values will eventually be reached. Since the learning time is always a finite value in practice, the learning mechanism will stop as the range values converge to fixed values or the learning time is longer than a predetermined value. In the latter case, it is possible that the performance index will be trapped in a local minimum. Therefore, in the following algorithm, the range values and their corresponding performance index after each random adjustment will be stored in a set. The range values which minimize the performance index is found by comparing all the elements in the set after the learning process. The details of the learning mechanism are described as follows.
Algorithm A (Learning mechanism for range values) Step Step
1. Let S be the set in which each element represents a set of range values and their corresponding performance index. Initially, this is an empty set. 2. Choose the time increment value At, the final learning time t~na/, and the probability comparison value 7.
296 Step Step Step Step Step Step Step Step
Step
Step Step
Chia-Ju W u and H u a n 9 C h i n g - H u o
Choose the nominal range values for each membership function. 4. According to procedure described in Section IV, find the performance index corresponding to the range values in Step 3. . Put the range values in Step 3 and their corresponding performance index into the set S as an element. 6. Set time t --- 0. 7. t = t-t-At. 8. Make random change in magnitude of each range value and add this magnitude to each corresponding range value. 9. Find the performance index corresponding to the range values in Step 8. 10. If the performance index is reduced, then: (1) retain the change in Step 8; (2) put the range values and their corresponding performance index into the set S as an element; and (3) go to Step 12; otherwise continue. 11. Calculate the probabilities for each magnitude change in Step 8. If the minimal value of these probabilities is greater than y, then: (1) retain the change in Step 8; and (2) put the range values and their corresponding performance index into the set S as an element; otherwise, return to their previous range values. 12. If t >>, t~,a~, then continue; otherwise go to Step 7. 13. Find the optimal range values which minimize the performance index by comparing all the elements in S. .
In Step 11 of the above algorithm, the probabilities are calculated by the Cauchy distribution as follows T
P(x) =
T 2+ x 2
(32)
where P ( x ) is the probability of a random change x, and T=
1 l+t"
(33)
The reason for choosing the Cauchy distribution in Eq. (32) is to lead the learning process to a linear convergence rate. A detailed discussion can be found in (14). Substituting Eq. (33) into Eq. (32) yields P(x) -
(1 +t) . 1 +x2(1 + t ) 2
(34)
For a random change of range value to be retained, if it does not reduce the performance index, the value of P ( x ) must be greater than the predetermined value y, i.e.
(1 +t) 1 +x~(1 + 8 ) 2
. . ~.
(35)
Solving for x in Eq. (35) yields
(1+8-~)
x2 < - -
~(l+t) 2
(36)
Fuzzy Parameter Identification of Direct Drive Robots
297
From this equation one can see that as time t increases, the values of x must be decreased to satisfy Eq. (36). Therefore, the probabilities of large random change are reduced gradually. In Step 8 of the above algorithm, there are many ways to determine the random change. For example, one can select a random number P(x) from a uniform distribution between zero and one. Then according to the gamma distribution
P(x) = e x p ( - x)
(37)
where P(x) is the probability of a random change with magnitude x, the magnitude x can be obtained as x = -ln(P(x)).
(38)
The corresponding random change is determined by adding randomly a positive or negative sign to the magnitude. Before presenting simulation examples to illustrate the proposed method, the process of parameter identification is summarized in the following algorithm.
Algorithm B (Parameter identification) Step Step Step
Step Step
Step Step Step Step Step
1. Generate identification data for M(O), G(O) according to the procedure described in Section III. 2. For every parameter in M(O) and G(O), perform the identification procedure in Steps 3-5. 3. Create a directed graph and start the identification procedure from the root node of the directed graph. Every node in the directed graph will denote a set of numbers, in which each one corresponds to the number of membership functions of an input variable. The number of fuzzy rules is also determined as the product of this set of numbers. 4. The optimal range values for each membership function are determined by the procedure described in Algorithm A. 5. If the sum of squared output errors obtained in Step 4 is less than a predetermined value or the number of fuzzy rules exceeds a predetermined number, then the identification procedure is completed. Otherwise, go to the next node in the directed graph according to the searching strategy of the breadth-first method. For the new node, go to Step 4 to find the optimal range values. 6. Generate identification data for V,.e,(O) according to the procedure described in Section III. 7. For every parameter in Vce,(O),repeat the identification procedure in Steps 3-5. 8. Generate identification data for V,or(O) according to the procedure described in Section III. 9. For every parameter in Vcor(O),repeat the identification procedure in Steps 3-5. 10. End.
After parameters in M(O), G(O), Vce,(O), and V,.o,(O) have all been estimated, the
Chia-Ju Wu and Huan9 Ching-Huo
298 Y-axis
nlz
,¢
/ /
¢ e2
/
/
/ /
O, X-axis FIG. 2. A two-link planar robot with rotary joints, where 0j and 02 denote the joint positions and
l~ and/2 denote the lengths of links. The dynamics equations are given in Eqs (39) and (40). computed-torque method can then be applied to make the robot move along a desired trajectory.
VI. Simulation Examples In the following simulation examples, the two-link planar robot shown in Fig. 2 will be used for illustration. The dynamics of this robot, which are assumed to be unknown, are given as follows (1, 2) Zl = [(ml
+m2)l~+m21~+2m21112cos 02]01 +(m21~+m:lll2 cos 02)02
--m21112(20,02+0~)sin 02 +
(m,
+m2)gll cos01 +m2g12cos(01 +02)
(39)
r2 = (m21z,+m21112cos 02)01 +m21202+m2111202 sin 02 +mzffl2cos(01 71-02). (4o) Comparing Eqs (39) and (40) with Eq. (9), one can find that the parameters to be estimated are
M,2]=[(m, +m2)l~+m21Z~+ 2m21112cos02 m21~+m21l12cos02] LM21 M22] m212+m21112COS02 m212
M=~M11
(41)
Fuzzy Parameter Identification of Direct Drive Robots
299 (42)
L Vco,2J
L
Vcen =[Vcenl 1 Vcenl2] =
0
0
--m2lll2sinOz I
LVcen21 Vcen22J [mzlll2 sin02 G = G2 =
[m+m2'cosO+m2'2osO+02,] m2ffl 2 cos(01 _t_02)
(43)
0
.
(44)
Through the simulation, the values of ml, m2, 1~, 12, and 9 are chosen to be 1 kg, 1 kg, 1 m, 1 m and 9.8 m/s 2, respectively. The ranges of 0t and 02 are assumed to be - r t ~< 0~ ~< rt and - r e ~< 02 ~< re.
Example 1 This example will show how to estimate the parameters in M(0), G(0), V,,n(O), and V,.or(O) according to the procedure described in Algorithm B. In applying Algorithm A to find the optimal range values, the values of (r, At, and ? are chosen to be 5 0 , 0 . 1 , and 0.5 s, respectively. The nominal range values for each membership function are chosen such that the membership functions are evenly distributed in the range [ - ~z, 7r] as shown in Fig. 3.
(a) Estimation of Mll and M12 During the procedure of generating identification data, one can find that: (1) M~ and M12 are functions of 02; (2) M21 is equal to Mj2 (since M is symmetric); and (3) Membership Functions
/
n u n
0 FIG. 3. Plots of the nominal membership functions when applying Algorithm A to find the optimal range values. If there are k membership functions associated with 01(02), then the membership functions have grade 0 at points 01(02) = - ~ + ( i - l ) ' 2 n / k , i = 1, 2, ..., k, and grade 1 at points 07(02) = - r~+ i" 2rc/k, i = 1, 2 . . . . . k.
Ch&-Ju Wu and Huan 9 Ching-Huo
300
TABLE I
Minimal values of Y,)°=°l(Error,) 2, maxi(Errori) 2 and min,(Errori)2 for k = 1,2 . . . . . 7, where k is the number of membership functions of 02 and Error, denotes [MH(i) --3;Ill(,)]for simplicity of notation k
E) °° (Error;) 2
maxi(Errorl) 2
1 2 3 4 5 6 7
129.2268 3.8635 1.1785 0.6543 0.2772 0.0956 0.0043
5.7896 1 0.2336 0.0284 0.0240 0.0072 1.7487" 10 -4
mini(Errori) 2 1.2693" 10 -6
2.3383" 10 -8 4 . 0 3 3 9 " 10 -6
1.7912" 10 -8 1.9075" 10 -7 4 . 3 5 2 9 " 10 -7 3.9178" 10 -9
M22 is a constant item. Therefore, only two parameters, namely M n and M~2, need to be estimated in the inertia matrix. Since MH is a function of 02, the fuzzy model to estimate M l l will be
Rule 1:
I f 02 is AI,
1 l then 3~ll = bo+bl "02
Rule k:
I f 02 is A k,
then AI11 = bok + b lk • 02
where _~t~l denotes the estimate of MI1 and k denotes both the number of fuzzy rules and the number of membership functions of 02. Suppose that 100 pairs of identification data for M~I and M~2 (02(/) ~ M . ( / ) , 02(/) -~ M12(/), i = 1,2 . . . . . 100, 02(0 = - n + [ ( i - 1)/99]" 2n) have been generated according to the procedure described in Section III. Then following the procedure in Algorithm B, the minimal values of Z~°° [ M ~ ( / ) - - M I I ( i ) ] 2 for k -- 1, 2 . . . . . 7, are found and shown in Table I. The plot of ( M l l - - M n ) for the case o f k = 7 is shown in Fig. 4. Following the same procedure, the minimal values of Z) °° [M~2(/)- ~ t l 2 ( i ) ] 2 for k = l, 2 . . . . . 7, are found and shown in Table II. The plot of (M~2- M~2) for the case of k = 7 is shown in Fig. 5.
(b) Estimation of Vcen21, V~12, and V~orl During the procedure of generating identification data, one can find that: (1) Vcen21,
Vcenl2, and V~orl are functions of 02, and (2) Vcenl 1 = Vcen2 2 = Vcor2 = O, Therefore, three parameters, namely V~e.2~, V~e.~2, and Vcor~, need to be estimated. Since the three parameters are functions of 02, their fuzzy models will be the same as that of Mll. Suppose that 100 pairs of identification data for Vce.21, Vcen~2, and V~or~ (02(/) Vcen21(/), 02(/) ~ Vcenl2(/), 02(/') ~ Veor~(/),i = 1,2 . . . . . 100,02(/) = --7r+[(i--1)/99] • 2r0 have been generated. Then following the procedure in Algorithm B, the minimal values of Z] °° [Vce.21(i)--l~c~.2~(/)] 2, E) °° [Vc~,,~z(/)--~,,~2(/)] 2, and E] °° [V~or~(i)--12~or~(/)] 2 for k = 1, 2 . . . . . 7 are found and shown in Tables III-V, respectively. The plots of (V~.2~ - 12~.2~), (V~e.,2- 12~.,2), and (V~or~- 12co.~) for the case of k = 7 are shown in Figs 6-8, respectively.
Fuzzy Parameter Identification of Direct Drive Robots
301
0.015] 0.01
0.005
o
i
-o.oo5
-0.01t -0.015 -4
-3
-2
2
02 FIG. 4. Plot of (Mjt--h4j~) with 7 membership functions associated with 02.
TABLE II
Minimal values of E~°°(Errori)2, maxi(Errorl) 2 and mini(Error32 for k = 1 , 2 , . . . , 7, where k is the number of membership functions of 02 and Errori denotes [Mi2(i)--~tl2(i)] for simplicity of notation k
E]°=° (Errori) 2
1 2 3 4 5 6 7
19.5107 0.6165 0.1831 0.1594 0.0839 0.0348 0.0013
maxi(Errori) 2 1.4516 0.0334 0.0233 0.0141 0.0091 0.0019 7.1864" 10 -5
mini(Errori)2 0 0 5.0728" 10 7 1.1613" ]0 -9
2.3873" 1O-11 1.1262" 10 - 7 3.1633" 10 1o
(c) Estimation o f GI and G2 D u r i n g the p r o c e d u r e o f g e n e r a t i n g identification data, one can find t h a t G1 a n d G2 are functions o f 01 a n d 02. Therefore, the fuzzy m o d e l o f G1 will be in the following form
302
Chia-Ju W u a n d H u a n 9 Chin#-Huo X 10 .3
8
6
4
2
I
=-
0
-2
-4
-6
-8
-10
-4
-3
-2
-1
0
1
2
3
O2 FIG. 5. Plot of (M~2-~1'I12) with 7 membership functions associated with 0a.
TABLEII1 Minimal values of El°=° (Errori) z, maxi(Errorl) 2 and mini(Errori) 2for k = 1,2 . . . . . 7, where k is the number of membership functions of 02 and Errori denotes [ Vc.e,Zl ( i) -- 12~e,21(i)] for simplicity of notation k
E) °° (Errorl) 2
1 2 3 4 5 6 7
15.8913 3.2934 0.4131 0.1861 0.0481 0.0137 0.0016
Rule l:
Rule (kl" k2):
If0~isA11
If01isA~ l
and
and
maxi(Errori)2 1.1337 0.2923 0.0495 0.0128 0.0030 0.0012 2.5404.10 -4
02 isA21,
02 isA~ z,
mini(Errori)2 0 0 6.3683" 10 7 4.3544.10 -8 1.5603.10 -7 1.8888.10 -9 0
then G1 = bo' + b l ' 0 ,
+b~'02
then Gl = bokl"k2+blkl'k2 . 01 .k_b~l.k2. 02
F u z z y P a r a m e t e r Identification o f Direct Drive Robots
303
TABLE IV
Minimal values of £ ] °=° (Err or i) 2, maxi( Error i) 2 and mini( Error i)2for k = l, 2 . . . . . 7, where k is the number of membership functions of 02 and Errori denotes [Vc,nl2(i) - I?~enl2(O]for simplicity of notation k
~,~°°(Errori) 2
1 2 3 4 5 6 7
14.1730 1.2100 0.4013 0.0863 0.0225 0.0193 0.0027
maxi(Errori) 2 0.9118 0.1847 0.0337 0.0089 0.0017 0.0023 3.2705" 10 4
mini(Errori) 2 0 0 7.4998" 1.3867" 1.2030" 5.2056" 2.3567"
l0 -7 10 7 10 s 10 7 10 ~
TABLE V
Minimal values of £ ]°° (Error i) 2, maxi( Error i)2 and mini( Error y for k = 1,2 . . . . . 7, where k is the number of membership functions of 02 and Errorl denotes [ V,,o,~(i) - V,.,,rl(i)1 for simplicity of notation k l 2 3 4 5 6 7
£]°_°(Errori) 2 41.0945 22.0589 4.2394 1.6625 0.3346 0.0880 0.0015
maxi(Errori) 2 2.5990 2.4096 0.3410 0.2254 0.0421 0.0069 1.2218" 10 -4
mini(Errori) 2 0 0 3.0397. l0 5 5.1873" 10 -7 3.6509" 10 -9
3.6089" 10 -7 3.7890" 10 -9
where 6~1 d e n o t e s the e s t i m a t e o f Gi, k l d e n o t e s the n u m b e r o f m e m b e r s h i p functions o f 01, a n d k2 d e n o t e s the n u m b e r o f m e m b e r s h i p functions o f 02. S u p p o s e t h a t 400 p a i r s o f identification d a t a for G1 a n d G2 ((0~(0, 02(j)) --* Gl(i x j ) , (Ol(i), 02(j)) ~ G2(i × j ) , i = 1, 2 . . . . . 2 0 , j = 1, 2 . . . . . 20, Ol(i) = - - n + [(i-- 1)/19]" 2n, 02(j) = - g - q - [ ( j - - 1 ) / 1 9 ] ' 2 n ) have been g e n e r a t e d p r i o r to a p p l y i n g A l g o r i t h m B. T h e n following the p r o c e d u r e in A l g o r i t h m B, the m i n i m a l values o f ~400 i= 1 [(G1 (i) -- (~1 (i))/y]2 for different c o m b i n a t i o n s o f k l a n d k2 are f o u n d a n d s h o w n in T a b l e VI. T h e p l o t o f ( G 1 - t~0/g for the case o f k l = 8 a n d k2 = 8 is s h o w n in Fig. 9. F o l l o w i n g the s a m e p r o c e d u r e , the m i n i m a l values o f w400 ,-,i=l [(G2 (i) - G2 (0)/912 for different c o m b i n a t i o n s o f k l a n d k2 are f o u n d a n d shown in T a b l e VII. T h e p l o t o f ( G 2 - G 2 ) / g for the case o f k l = 8 a n d k2 = 8 is shown in Fig. 10.
Example 2 This e x a m p l e shows h o w to a p p l y the P I D c o m p u t e d - t o r q u e c o n t r o l l e r b a s e d on the e s t i m a t e d d y n a m i c s to m a k e the r o b o t m o v e a l o n g the desired t r a j e c t o r y (1)
Chia-Ju Wu and Huang Ching-Huo
304 0.02
0.015
0.01
~d 0.005
I ~o
-0.005
-0.01 -4
~ -3
~ -2
~ -1
2
0
3
02
FIG. 6. Plot of (Vcen2I -- Veen21) with 7 membership functions associated with 02.
01a = 0.1 • sin(rtt)
(45)
02a = 0.1" cos(rtt).
(46)
As described in Section II, the PID computed-torque control law is given as
z -- M(O)(Oa+Kv~+Kpe+K~Se'dt)+ 12,.or(O)[O0]+Vcen(O)[Oz]--~a(o)
(47)
where/14(0), 12cor(0), and ~rcen(O ) a r e estimated with 7 membership functions associated with 02, and (~(0) is estimated with 8 membership functions associated with both 01 and 02. If the initial condition of [01 02] is chosen to be [0.1 0] and Kv, Kp, and K~ are chosen as
K,,=
Kp
and
=
50
oI
625
Fuzzy Parameter Identification of Direct Drive Robots
305
°°2 /
0.0151 0.01
/
0.005
I
-0.005 -0.01 /
-0.01 5 I
I
0 02 -4
t -3
L -2
, -1
~ 0
~ 1
r 2
3
02 FIG. 7. Plot of (Vce.12- I?,.e.12)with 7 membership functions associated with 02.
K; =
1200
then the plots of (O,d--0,) and (02d--02) are as shown in Figs 11 and 12, respectively.
VII. Conclusions
This work showed how to generate identification data in a recursive manner and how to apply a fuzzy model to estimate parameters in dynamics of a direct drive robot. From the simulation results shown in Tables I-VII, one can find that AIr(0), l?ce,(0), 12cot(0), and 6~(0) approximate the original parameters well. It is also expected that one can increase the number of training pairs to reduce the differences between the estimates and the original parameters. Based on the estimated dynamics, it was shown in Example 2 that the PID computed-torque controller gave a good result in trajectory-tracking problem of robots. Though the applications of the proposed method in the simulation examples are successful, it should be pointed out that heavy computation burden is inevitable in the identification procedure. Most of the computation time will spend in generating identification data and in searching for the optimal range values of each membership
Chia-Ju Wu and Huan9 Chiny-Huo
306 0.01
0.005
-0.005
-0.01
-0.015
' -3
-4
P -2
'
' 0
-1
' 1
' 2
3
O2
FIG. 8. Plot of
(V,.o,j-
I2,-o,0 with 7 membership functions associated with 02.
TABLE VI
Minimal values of 5~4°_°(Errori)2, maxi(Errori) 2 and mini(Errori) 2 for different combinations of kl and k2, where kl and k2 are the numbers of membership functions of O~ and 02, respectively, and Errori denotes [GI(i) - Gl(i)]/g for simplicity of notation kl
k2
Z]°- ° (Errori) 2
maxi(Errorl) 2
mini(Errori) 2
1 2 3 4 5 6 7 8
1 2 3 4 5 6 7 8
481.3954 165.8560 114.4186 40.0132 3.5234 2.3358 0.7215 0.0724
8.9184 4.4778 5.3172 1.6219 0.1778 0.0946 0.0263 0.0017
3.8526.10 5 1.2495 • l0 -5 4.1347.10 -8 1.0243" 10 -6
1.9926" 10 8 7.0697" l0 -9 2.6583" 10 -11 1.3458" 10 -28
function. Unless the c o m p u t a t i o n complexity o f the identification procedure is reduced, the p r o p o s e d m e t h o d can only be used in a off-line m a n n e r . Therefore, how to reduce the c o m p u t a t i o n complexity of the p r o p o s e d m e t h o d is suggested as a topic for further
Fuzzy Parameter Identification o f Direct Drive Robots
307
(Gt-al)/g 0.040.02o -0.02-0.04-
2
4
-4
-4
01
FIG. 9. Plot of (G~- (~1)/# with 8 membership functions associated with both 0~ and 02.
TABLE VII
Minimal values of Y?°_°(Errori) 2, maxi(Errori): and mini(Errori) 2 for difJerent combinations of kl and k2, where ks and k2 are the numbers of membership functions of 01 and 02, respectively, and Error, denotes [G2(i) - G2(i)]/g for simplicity of notation kl
k2
1 2 3 4 5 6 7 8
1 2 3 4 5 6 7 8
E]°°(Errori) z 162.6653 77.2357 19.9503 7.3790 1.4301 0.2837 0.0987 0.0329
maxi(Errori) z 1.6900 2.1638 1.6308 0.2935 0.0686 0.0217 0.0016 0.0007
min,(Errori) 2 1.5618 2.0618 1.5097 1.3840 3.2439 1.1431 7.0626" 3.8807"
10 -5 10 7 10 -7 10 -6 10 v 10 9 10 lo 10 9
research. I n a d d i t i o n , since it is a s s u m e d that the j o i n t velocities, j o i n t accelerations, a n d j o i n t torques c a n be m e a s u r e d w i t h o u t error, the m e a s u r e m e n t of j o i n t variables is suggested as a n o t h e r research topic.
308
Chia-Ju Wu and Huan9 Chiny-Huo
0.030.02-
0.01
-
A
0-0.01
-
-0.02-
-0.03> 4
o
/ -4
-4 FIG. 10. Plot of
-2 01
(G2-G2)/gwith 8 membership functions associated with both
01 and 02.
0.1
0.05
-0.05
j i 0 5
~ 1
15 1.
I 2
I 2.5
I 3
3.5
I 4
Time(sec.)
FIG. 11. Plot of (0~-- 0j), where 0~d is given in Eq. (45).
t
4.5
Fuzzy Parameter Identification o f Direct Drive Robots
309
0.02,
-0.02
I
-0.04
,~a~ "~1
-0.06
-0.08
r
"0"10
I
I
I
I
I
[
I
I
I
0.5
1
1.5
2
2.5
3
3.5
4
4.5
5
Time(see.) FIG. 12. Plot of (02d--02), where 02d is given in Eq. (46).
A cknowledgements This work was supported by the National Science Council, Taiwan, Republic of China, under grant No. NSC83-0422-E-224-036.
References (1) F. L. Lewis, C. T. Abdallah and D. M. Dawson, "Control of Robot Manipulators", Macmillan, New York, NY, 1993. (2) J. J. Craig, "Introduction to Robotics: Mechanics and Control", 2nd Edn, AddisonWesley, Reading, MA, 1988. (3) M. Sugeno and K. Tanaka, "Successive identification of a fuzzy model and its applications to prediction of a complex system", Fuzzy Sets Syst., Vol. 42, pp. 315-334, 1991. (4) T. Takagi and M. Sugeno, "Fuzzy identification of systems and its applications to modeling and control", IEEE Trans. Syst. Man Cybernet., Vol. SMC-15, No. 1, pp. 116-132, 1985. (5) H. Asada, T. Kanade and I. Takeyama, "Control of a direct-drive arm", ASME J. Dynam. Syst., Meas., Contr., Vol. 105, pp. 136-142, 1983. (6) H. Asada and K. Youcef-Toumi, "Analysis and design of a direct-drive arm with a fivebar-link parallel drive mechanism", ASME J. Dynam. Syst. Meas. Contr., Vol. 106, pp. 225-230, 1984. (7) H. Asada and K. Youcef-Toumi, "Direct-Drive Robots: Theory and Practice", MIT Press, Cambridge, MA, 1987. (8) K. Youcef-Toumi and A. T. Y. Kuo, "High-speed trajectory control of a direct-drive manipulator", IEEE Trans. Robot. Automat., Vol. 9, pp. 102-108, 1993.
310
Chia-Ju Wu and Huan9 Ching-Huo
(9) R. L. Kruse, B. P. Leung and C. L. Tondo, "Data Structures and Program Design in C", Prentice-Hall, Englewood Cliffs, N J, 1991. (10) B. Kosko, "Neural Networks and Fuzzy Systems", Prentice-Hall, Englewood Cliffs, N J, 1992. (11) P. D. Wassermann, "Neural Computing: Theory and Practice", Vann Nostrand, Reinhold, NY, 1989. (12) C. C. Lee, "Fuzzy logic in control systems: fuzzy logic controller, part II", IEEE Trans. Syst. Man Cybernet., Vol. 20, pp. 419-435, 1990. (13) T. N. E. Greville, "The pseudoinverse of a rectangular or singular matrix and its applications to the solutions of systems of linear equations", SIAM Rev., Vol. 1, pp. 38-43, January, 1959. (14) H. Szu and R. Hartley, "Fast simulated annealing", Phys. Lett., Vol. 1222 (3,4), pp. 157162, 1987.