Dynamic identication of a Mitsubishi

PA-10 robotic manipulator

R. van der Aalst B.Sc.
DCT 2008.050
Traineeship report
Coach(es): Dr. Marcelo H. Ang Jr
(National University of Singapore)
Lim C.W. M.Eng.
(Singapore Institute of Manufacturing Technology)
Supervisor: Prof.dr. H. Nijmeijer
(Eindhoven University of Technology)
Technische Universiteit Eindhoven
Department Mechanical Engineering
Dynamics and Control Group
Eindhoven, April, 2008
i
Abstract
In this report the system parameters identication for a Mitsubishi PA-10 robotic manipulator is
considered. The Singapore Institute of Manufacturing Technology has, in collaboration with the
National University of Singapore, developed a specic method to do this system identication.
This report focusses on using that identication method.
The main advantages of the method are that it consists of a structured procedure to follow so
it is easily applied to dierent kind of robotic manipulators and furthermore no need of measured
acceleration information is required.
Before the identication procedure can be used the full dynamic model of the robotic manip-
ulator is derived. Reduction of system parameters is structurally applied because of the necessity
to reduce system complexity if to be used for identication purposes.
The identication procedure is based on achieving a fully undamped system by using mod-
eled friction compensation. Originally the friction compensation is done based on the assump-
tion of viscous friction being dominantly present in the robotic joints. This does not result in
undamped system behaviour. It is shown by means of experimental friction identication that
Coulomb friction is signicantly present in the robotic joints. Relating to this result the friction
compensation model is expanded with Coulomb friction using experimentally identied param-
eter values obtained from a least-square t through gathered data points. Despite using this
expended friction model as compensation, the system is still eected by friction inuences that
remain because of the real friction deviating from the modeled friction. It is shown that this
deviations between modeled, and thus compensated friction, and real friction result in specic
system behaviour that is clearly not nicely undamped.
It is concluded that a friction model, gained by experimental means, can result in control
improvements while used to compensate the friction but that it is not possible to succeed in
obtaining a fully undamped system. There will always be deviations between the modeled friction
used to compensate and the real friction present. Therefore the identication method should
be regarded impractical and it can not be used to identify system parameters as supported by
erroneous identication results.
ii
iii
Contents
1 Introduction 3
2 Robot modeling 5
2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
2.2 The Mitsubishi PA-10 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5
2.3 Kinematics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6
2.4 Dynamic rigid-body modeling . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7
2.5 Base parameter set . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 8
2.5.1 Denition of the inertial parameters . . . . . . . . . . . . . . . . . . . . . 8
2.5.2 Elimination of inertial parameters . . . . . . . . . . . . . . . . . . . . . . 9
2.5.3 Regrouping . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 9
2.6 Friction modeling. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 10
3 Literature review 13
3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 13
3.2 Some available identication methods . . . . . . . . . . . . . . . . . . . . . . . . 13
4 Experimental identication of the inertial parameters 15
4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15
4.2 Model reduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 15
4.3 Methodology . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 16
4.3.1 Experiment in practice . . . . . . . . . . . . . . . . . . . . . . . . . . . . 20
5 Experimental results 21
5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21
5.2 Practical problems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 21
5.3 Results original procedure . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 23
5.4 Results using identied friction . . . . . . . . . . . . . . . . . . . . . . . . . . . . 28
5.5 Friction inuence in practice . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 30
6 Conclusions and Recommendations 33
6.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
6.2 Contribution . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 33
6.3 Recommendations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
A Equations Base Parameter Set 37
1
B Experimental setup 39
C Friction identications 40
2
Chapter 1
Introduction
Nowadays robotic manipulators take more and more a signicant place in industry. They improve
quality and safety while reducing costs and processing time in fully automated production lines.
In these automated production lines robotic manipulators can perform a wide variety of tasks
ranging from welding, cutting, bonding and painting to material handling. As research pushes
performance limits more possibilities for application of robotic manipulators arise. A promising
technique that is currently topic of many researches to achieve robotic control improvement is
the usage of model-based control.
At the Singapore Institute for Manufacturing Technology (SIMtech) the robotic team is work-
ing on the implementation of force control on robotic manipulators. To achieve this, one of the
major challenges is to implement such a model-based controller in a well-structured manner.
This subgoal on its own, the implementation of model-based control, possesses challenges as
describing the complete dynamic behaviour of the robot using a mathematical model. To be able
to use such a mathematical model describing the robotic manipulator dynamics, it is necessary
to know the values of all system parameters. The robotic manipulator used is a Mitsubishi PA-10
as visible in gure 1.1.
The specic objective of this research is dened as follows; describe the complete dynamics of
the Mitsubishi PA-10 robotic manipulator and identify all unknown system parameters, used in
the obtained dynamic model, by experimental means. To do this identication an experimental
method as proposed by SIMtech
1
[14] is applied.
This report is organized as follows. First, in Chapter 2, the robotic dynamics will be addressed
using a mathematical model. Furthermore this model is prepared for usage in the experimental
identication of unknown parameters. Before the experimental identication method is fully
addressed in Chapter 4, some alternatives are briey mentioned in Chapter 3. The theory of
Chapter 4, regarding the identication experiment, is applied on the Mitsubishi PA-10. Chapter
5 presents and discusses the results of the experiment. Finally general discussion takes place in
Chapter 6 followed by recommendations for future research progress.
1
Developed in collaboration with the National University of Singapore (NUS).
3
Figure 1.1: The Mitshubishi PA-10 robotic manipulator.
4
Chapter 2
Robot modeling
2.1 Introduction
The general objective of SIMtech is to implement model-based control. To make this possible a
dynamic model has to be derived. After the general introduction of a Mitsubishi PA-10 robotic
manipulator in section 2.2, some coordinate frame conventions are dened in section 2.3. In
section 2.4 the equations of motion are stated. Finally the set of minimal inertial parameters
in which the equations of motion can be expressed are determined in section 2.5, besides some
other assumptions that simplify the model.
2.2 The Mitsubishi PA-10
Figure 2.1: Graphical representation of the Mitshubishi PA-10 robotic manipulator [5].
The Mitsubishi PA-10 is a robotic manipulator as depicted in gure 2.1. It possesses seven de-
grees of freedom, thus a so-called shoulder, elbow and wrist are present. Indicated are the limits
of the joint rotations and the maximum joint velocities.
5
2.3 Kinematics
To be able to describe the robotic kinematics a coordinate frame has to be assigned to each joint.
The coordinate frame assignment as used throughout this report is visible in gure (2.2) and
is conform the modied Denavit-Hartenberg (MDH) convention [1]. The corresponding MDH-
parameters are listed in table 2.1 and according to the MDH convention dened as:
a
i
The distance from z
i
to z
i+1
measured along x
i
,

i
The angle between z
i
and z
i+1
measured about x
i
,
d
i
The distance from x
i1
to x
i
measured along z
i
,

i
The angle between x
i1
and x
i
measured about z
i
,
where i = 1, . . . , 7. Furthermore the frames are placed in such a way that a positive sensor
feedback corresponds to the positive direction conform the MDH convention.
-
6
rh x
0
y
0
z
0
-
6
rh x
1
y
1
z
1
-
?
rh
x
2
z
2
y
2
-
6
rh x
3
y
3
z
3
-
?
rh
x
4
z
4
y
4
-
6
rh x
5
y
5
z
5
-
?
rh
x
6
z
6
y
6
-
6
rh
x
7
y
7
z
7
Figure 2.2: Links and frame assignment of the
PA-10.
link i a
i

i
d
i

i
1 0 0 0

1
2 0 90 0

2
3 0 90 d
3

3
4 0 90 0

4
5 0 90 d
5

5
6 0 90 0

6
7 0 90 d
7

6
Table 2.1: Modied Denavit-Hartenberg pa-
rameters
6
Using this frame assignment the kinematic relation between frame i and frame i + 1, and corre-
sponding joints, is described with the following homogeneous transformation matrix:
T
i
i+1
=
_

_
cos(
i
) sin(
i
) 0 a
i
sin(
i
)cos(
i
) cos(
i
)cos(
i
) sin(
i
) sin(
i
)d
i
sin(
i
)sin(
i
) cos(
i
)sin(
i
) cos(
i
) cos(
i
)d
i
0 0 0 1
_

_
(2.1)
Furthermore the PA-10 is an open kinematic chain and therefore the forward kinematics can be
computed as a product of homogenous transformations between circumjacent coordinate frames:
T
0
i
= T
0
1
T
1
2
. . . T
i1
i
(2.2)
2.4 Dynamic rigid-body modeling
A rigid-body model is used to describe the robot dynamics. The rigid-body model can be de-
rived using the Euler-Lagrange formulation [1,2,3]. The general form of the model with friction
included is as follows:
(t) = D(q(t)) q(t) +C(q(t), q(t)) +h(q(t)) +
f
(t) (2.3)
Where q is a vector representing the joint angles, q = [
1
...
7
]
T
, and q and q represent joint veloc-
ities and joint accelerations respectively. Furthermore D(q(t)) represents the inertia matrix and
C(q(t)), h(q(t)) and
f
(t) represent vectors that contain the Coriolis/centripetal, gravitational
and friction forces respectively. The elements of inertia matrix D(q(t)) are dened as:
d
j,k
=
7

i=max(j,k)
Tr
_
_

q
k
T
0
i
_
J
i
_

q
j
T
0
i
_
T
_
(2.4)
and the elements of C(q(t), q(t))) are dened using Christoel symbols:
c
j
=
7

k=1
7

l=1
c
j,k,l
q
k
q
l
, (2.5)
c
j,k,l
=
7

i=max(j,k,l)
Tr
_
_

2
q
k
q
l
T
0
i
_
J
i
_

q
k
T
0
i
_
T
_
j, k, l = 1, . . . , 7 (2.6)
Furthermore the gravitational eects are added in h(q(t)) as:
h
j
=
7

i=j
_
m
i
g
_

q
j
T
0
i
_
r
i
_
j = 1, . . . , 7 (2.7)
Here m
i
is the mass of the ith link, r
i
contains the homogeneous coordinates of the center of
mass of link i expressed in the ith coordinate frame:
r
i
= [x
i
y
i
z
i
1]
T
(2.8)
7
and g represents a gravity vector g = [0 0 g 0] where g = 9.8062m/s
2
. J
i
is the pseudo inertia
tensor described as:
J
i
=
_

_
I
xx
i
+I
yy
i
+I
zz
i
2
I
xy
i
I
xz
i
m
i
x
i
I
xy
i
I
xx
i
I
yy
i
+I
zz
i
2
I
yz
i
m
i
y
i
I
xz
i
I
yz
i
I
xx
i
+I
yy
i
I
zz
i
2
m
i
z
i
m
i
x
i
m
i
y
i
m
i
z
i
m
i
_

_
(2.9)
Here I
xx
i
, I
yy
i
and I
zz
i
denote the principal moments of inertia and I
xy
i
, I
yz
i
and I
xz
i
the prod-
ucts of inertia for link i. Furthermore m
i
x
i
, m
i
y
i
and m
i
z
i
represent the rst moment of link
i. All expressed with respect to the ith coordinate frame as dened in the modied Denavit-
Hartenberg convention.
2.5 Base parameter set
The derived dynamic model, equation (2.3), is because of its size very impractical to use for iden-
tication or real-time control purposes. To make is suitable the model size has to be reduced. It is
well known that such a dynamic model, excluding friction, can be linearly expressed in its inertial
parameters. Each link is described using 10 inertial parameters, so for the Mitsubishi PA-10 70
inertial parameters in total. However, some of these parameters have no eect on the dynamic
model and others may appear in linear combinations. By eliminating the inertial parameters that
have no eect on the dynamic model and by regrouping the linear combinations of the inertial
parameters one can get a set of minimum inertial parameters, the base parameter set (BPS), that
still fully describes the manipulators dynamics. To do this elimination and regrouping by hand
is a very tedious and time-consuming operation. Furthermore it is not sure that the optimal so-
lution of the BPS is obtained. Here a structured method is used as proposed by [15,16,17] to nd
the BPS.
2.5.1 Denition of the inertial parameters
Again m
i
denotes the mass and m
i
x
i
, m
i
y
i
and m
i
z
i
the rst moments of link i. The inertia
tensor of link i about the origin of coordinate frame i is denoted as:
I
i
i
=
_
_
I
xx
i
I
xy
i
I
xz
i
I
xy
i
I
yy
i
I
yz
i
I
xz
i
I
yz
i
I
zz
i
_
_
(2.10)
The inertial parameters of the ith link are collected in a vector X
i
dened as:
X
i
= [I
xx
i
I
xy
i
I
xz
i
I
yy
i
I
yz
i
I
zz
i
m
i
x
i
m
i
y
i
m
i
z
i
m
i
]
T
(2.11)
Furthermore X becomes the (70 x 1) vector that contains all robot inertial parameters dened as:
X =
_
X
T
1
X
T
2
. . . X
T
n

T
(2.12)
Because of equation (2.3), except for friction, being linear in its inertial parameters it can be
linearly rewritten using the last denitions stated above as:
(t) = H(q(t), q(t), q(t))X +
f
(t) (2.13)
where H is the regression matrix.
8
2.5.2 Elimination of inertial parameters
Let E denote the kinetic energy and U denote the potential energy of the system. Because of the
linearity in inertial parameters these can be rewritten as:
E =
70

j=1
E
X
j
X
j
(2.14)
U =
70

j=1
U
X
j
X
j
(2.15)
Here X
j
refers to the jth element of X.
If
E
X
j
= 0 and
U
X
j
= constant then the corresponding X
j
has no eect on (t) and is therefore
not an element of the BPS. It can be eliminated by putting it to zero in equation (2.13).
In case of the PA-10 the following inertial elements satisfy these conditions and therefore can be
eliminated, i.e. put to zero:
I
xx
1
, I
xy
1
, I
xz
1
, I
yy
1
, I
yz
1
, mX
1
, mY
1
, mZ
1
, m
1
, mZ
2
, m
2
(2.16)
2.5.3 Regrouping
The next step to reduce the number of elements of the BPS is to regroup elements that are linearly
dependent. In case of the PA10 all joints are rotational and the inertial elements I
yy
i
, mZ
i
and
m
i
can be eliminated by regrouping them using the following relations, starting from the last
link down to the rst:
I
xxR
i
= I
xx
i
I
yy
i
I
xxR
i1
= I
xx
i1
+I
yy
i
+2d
i
mZ
i
+d
2
i
m
i
I
yyR
i1
= I
yy
i1
+cos
2
(
i
)I
yy
i
+2d
i
cos
2
(
i
)mZ
i
+(a
2
i
+d
2
i
cos
2
(
i
))m
i
I
zzR
i1
= I
zz
i1
+sin
2
(
i
)I
yy
i
+2d
i
sin
2
(
i
)mZ
i
+(a
2
i
+d
2
i
sin
2
(
i
))m
i
I
xyR
i1
= I
xy
i1
+a
i
sin(
i
)mZ
i
+a
i
d
i
sin(
i
)m
i
I
xzR
i1
= I
xz
i1
a
i
cos(
i
)mZ
i
a
i
d
i
cos(
i
)m
i
I
yzR
i1
= I
yz
i1
+cos(
i
)sin(
i
)I
yy
i
+2d
i
cos(
i
)sin(
i
)mZ
i
+d
2
i
cos(
i
)sin(
i
)m
i
mX
R
i1
= mX
i1
+a
i
m
i
mY
R
i1
= mY
i1
sin(
i
)mZ
i
d
i
sin(
i
)m
i
mZ
R
i1
= mZ
i1
+cos(
i
)mZ
i
+d
i
cos(
i
)m
i
m
R
i1
= m
i1
+m
i
_

_
(2.17)
So the subscript R denotes a replacement for the original inertial parameter and represents
a group of original inertial parameters. The original inertial elements that are regrouped to such
a replacement can be put zero in (2.13). With the MDH parameters as given in table 2.3 the
regrouping relations for each specic joint can be calculated as presented in appendix (A).
9
The remaining elements in the BPS for each link then become:
Link 7: I
xxR
7
, I
xy
7
, I
xz
7
, I
yz
7
, I
zz
7
, mX
7
, mY
7
Link 6: I
xxR
6
, I
xy
6
, I
xz
6
, I
yz
6
, I
zzR
6
, mX
6
, mY
R
6
Link 5: I
xxR
5
, I
xy
5
, I
xz
5
, I
yz
5
, I
zzR
5
, mX
5
, mY
R
5
Link 4: I
xxR
4
, I
xy
4
, I
xz
4
, I
yz
4
, I
zzR
4
, mX
4
, mY
R
4
Link 3: I
xxR
3
, I
xy
3
, I
xz
3
, I
yz
3
, I
zzR
3
, mX
3
, mY
R
3
Link 2: I
xxR
2
, I
xy
2
, I
xz
2
, I
yz
2
, I
zzR
2
, mX
2
, mY
R
2
Link 1: I
zzR
1
The last step is to detect the existence of supplementary combinations for link 1 to 7, where
link k denotes the rst link for which holds:
E
X
j
= 0 and
U
X
j
= constant for all the ten inertial
parameters corresponding to that specic link. No supplementary regrouping is possible though.
So by using elimination and regrouping the number of inertial elements that describe the
dynamic model is reduced from 60 to a BPS of 43. Expressing (2.13) using the BPS elements
gives:
(t) = H
bps
(q(t), q(t), q(t))X
bps
+
f
(t) (2.18)
Where X
bps
is a (43 x 1) vector of the BPS elements and H
bps
contains the collomns of the original
H that correspond to the original elements in X that have nowbecome part of the base parameter
set X
bps
.
To guarantee that our reduced model is correct and still represents the original system as
described in (2.3) it is veried using a robotic model made in the Robotic Toolbox [4]. From this
toolbox it is known that it gives reliable and accurate results.
2.6 Friction modeling.
The last element modeled in equation (2.3),
f
(t), is to include the inuence of friction forces.
Friction is a very complex issue on its own and the emphasis of this report is initially on the
identication of the inertial parameters. Later however, while describing experimental results
(Chapter 5), the inuence of friction becomes signicant. Several general friction models will be
used throughout this report. Despite the fact that the used models here are generally known in
common engineering practice there is a large variety in terminology. The following denitions
represent the terminology and corresponding mathematical formulas for the several friction mod-
els as used throughout this report
1
:
Coulomb friction: [gure 2.3] This friction component delivers a force when joint velocity is
present, q = 0. The applied force has opposite sign to the velocity direction and the magnitude
1
Note that the direction of the force in the gures representing friction models is opposite from the force
direction as really experienced by the joints. The reason for this is that in the model, equation (2.3), the
friction forces are already translated to the right-hand sight.
10
0
0
f
c
f
c
Velocity
F
r
i
c
t
i
o
n

Figure 2.3: Coulomb friction.

0
0
f
v
Velocity
F
r
i
c
t
i
o
n

Figure 2.4: Viscous friction.

of the force is independent of the velocity magnitude. The equation used to represent Coulomb
friction is as follows:

f
c
(t) = f
c
sign( q) (2.19)
Viscous friction: [gure 2.4] Similar to the Coulomb friction the viscous friction force is only
present in case of joint velocity, q = 0. The applied force is linear proportional to the magnitude
of the velocity. Again the force the joint experiences has opposite sign of the velocity. The equation
used to represent viscous friction is as follows:

f
v
(t) = f
v
q (2.20)
0
0
f
c
f
c
f
v
Velocity
F
r
i
c
t
i
o
n

Figure 2.5: Coulomb and viscous friction com-

bined.
0
0
f
s
f
s
Velocity
F
r
i
c
t
i
o
n

Figure 2.6: Coulomb and viscous friction with

added static friction.
Static friction: [gure 2.6] In contrast to the Coulomb and viscous friction, that are only present
in case of velocity being present, this friction delivers a force only at zero velocity. The static
friction represents a force that has to be overcome by an initial force to initiate velocity. It is
11
opposite to the force initially applied and disappears as the joint has overcome the static friction
value f
s
, that is the joint breaks away and possesses a velocity.
Stribeck friction: A friction phenomenon that occurs in the low velocity region. This phe-
nomenon is responsible for a nonlinear eect that decreases the friction force in the specic
region till viscous friction dominates. However, as mentioned before, the main emphasis of this
report is not on detailed friction identication and throughout this report the nonlinear phe-
nomena that occur in the friction force because of the Stribeck eect are neglected. This is a
reasonable assumption for the research done here because rst of all it only eects the low ve-
locity region and secondly detailed friction identications done at Mitsubishi PA-10 robots have
shown that the presence of the Stribeck eect is not really signicant [18,19].
In practice the friction models described above are not stand-alone but occur all together in the
robot joints. One of the friction combinations that is used throughout this report is the combina-
tion of Coulomb friction with viscous friction. So obvious this model becomes:

f
vc
(t) = f
v
q +f
c
sign( q) (2.21)
The corresponding friction graph is visible in gure (2.5)
12
Chapter 3
Literature review
3.1 Introduction
To implement model based control using an accurate dynamic model it is necessary to know all
the dynamic parameters. The dynamic model parameters as provided by the manufacturer are
insucient. To acquire the unknown dynamic parameter values, identication experiments are
necessary. SIMtech has already proposed such an experimental method as described in the next
Chapter. Although this paper focuses on the use of that proposed method this Chapter briey
mentions other identication procedures for completeness.
3.2 Some available identication methods
The most straightforward identication method, regarding mass and inertia properties, would
be dismantling the robotic manipulator. However, in most cases this is not a realistic option be-
cause of the manipulators complexity. Another obvious method that is mentioned several times
in literature [6,7,8] is the use of CAD drawings. This is not attractive though because robot man-
ufacturers do not provide these drawings for all parts of the robot and it will give approximate
values of which the accuracy is doubtful [7]. Other identication methods are known to make use
of modal analysis [8] or frequency spectrum analysis [9].
However, most literature encountered during browsing describe slight variations on the following
general procedure [7]:
The dynamic model of the robot manipulator is derived in symbolic form.
The minimal set of, to be identied, dynamic parameters that still fully describes the ma-
nipulator dynamics, i.e. the base parameter set, is determined.
Model simplications are done by reasonably neglecting the contribution, for example be-
cause of the manipulators structure, of certain parameters.
Specic trajectories are determined to gain a set of experimental data that optimizes the
identication result.
Finally the required values of the dynamic parameters are obtained by applying a suitable
estimation algorithm.
13
Dierences in the general procedure can be found in dening optimized trajectories and esti-
mation algorithms but are mainly found in the use of dierent dynamic models for use of iden-
tication. Three models are commonly mentioned. First the so called explicit dynamic model,
a dynamic model which depends on the joint acceleration and needs an explicit dierentiation
of the velocity. The main disadvantage of this is the introduction of noise in the signal by dif-
ferentiation of the velocity to end up with the acceleration. Secondly a model that avoids this
disadvantage is a model that dierentiates a function of the velocity, the so called implicit dy-
namic model. However, the computational burden of this implicit model is much higher because
of a complicated expression and no improvement in identication can be noticed [10]. The last
method that gets a lot of attention in research literature is the use of a power model, which expres-
sion is obtained from rewriting an energy model in a dierential form. The main advantage here
is not the necessity to use acceleration information. The use of a power model receives positive
criticism because of its simplicity and accuracy [10,11].
The above approaches suer however from an important drawback [12]; they depend on
torque information gained from actuator current measurements. The real torque applied to the
joint deviates from the used measured torque information because of friction inuences and
limited precision of the actuator torque constants. An alternative model combines this internal
measured information with externally gained information of the reaction forces and torques on
the robotic manipulators bedplate. This is measured by means of an external force/torque plat-
form and thus totally independent of internal torques such as joint friction, allowing an accurate
estimate of the inertial parameters. However, it must be noted that if the researcher is also inter-
ested in identifying joint friction parameters then a dynamic model that combines internal and
external information is obviously a necessity. Such a combined model results in only a slight
increase of the number of model parameters to be identied but a considerable increase of the
number of force/torque data involved [12]. A obvious disadvantage is the addition of an external
platform to the robotic manipulator to be able to measure reaction forces/torques. This is not a
likely option in an industrial environment.
Despite this range of research methods, experimental identication of dynamical robotic manip-
ulator parameters remains an active research topic for use in model-based control applications.
Although almost all papers show control improvement by experimental means when identied
parameters are used, no consensus is reached yet [13]. For exactly the same robotic manipula-
tor very dierent identied values exist([note; after expressing with use of the same denitions]).
This deviation in identied results leaves the possibility of trying other methods to identify dy-
namical parameters for use in model-based control.
SIMtech has, in cooperation with the National University of Singapore (NUS), developed an iden-
tication method that is simple and straightforward with the main advantage of not using accel-
eration information [14]. Here, for identifying system parameters of the Mitsubishi PA-10 robotic
manipulator, this proposed method is used. In the next Chapter the theory of the identication
method will be explained, later in Chapter 5 experimental results will be shown and commented
on.
14
Chapter 4
Experimental identication of the
inertial parameters
4.1 Introduction
To make the implementation of model based control possible it is necessary to know the values
of the inertial parameters. The values of the link masses and their center of masses are supplied
by the manufacturer Mitsubishi [5]. The inertias tensor elements are not known though. They
should be identied by experimental means. In section 4.2 a last model simplication is done to
reduce the number of unknowns that have to be experimentally determined. In section 4.3 the
experimental method to determine the unknown inertial parameters is described.
4.2 Model reduction
To reduce the number of unknowns in the BPS the assumption is made that the robotic links
are perfectly symmetric. In that case the inertia tensor expressed with respect to a frame that
coincides with the center of mass of link i becomes[2,3]:
I
cm
i
i
=
_
_
I
cm
xx
6
0 0
0 I
cm
yy
6
0
0 0 I
cm
zz
6
_
_
(4.1)
So the products of inertia become zero. If this result is translated to a inertia tensor which is with
respect to coordinate frame i of link i, instead of the frame in the center of mass, then the values
of those diagonal elements are known because of the Steiner translation formula:
I
i
i
=
_
_
I
xx
i
I
xy
i
I
xz
i
I
xy
i
I
yy
i
I
yz
i
I
xz
i
I
yz
i
I
zz
i
_
_
=
_
_
I
cm
xx
i
0 0
0 I
cm
yy
i
0
0 0 I
cm
zz
i
_
_
+m
i
S(r
i
)
T
S(r
i
) (4.2)
where r
i
denotes the position vector of the center of mass link i expressed in frame i and the
skew symmetric matrix S dened as:
S
i
=
_
_
0 r
z,i
r
y,i
r
z,i
0 r
x,i
r
y,i
r
x,i
0
_
_
(4.3)
15
With the assumption of all links being symmetric, all the values of the products of inertia are
known. They are equal to the contribution of the Steiner relation. Furthermore the masses and
center of masses are give by the manufacturer so the rst moments of inertia can also be regarded
as known parameters. These known inertia elements reduce the unknown parameters in the BPS
to 13.
4.3 Methodology
An specic experimental procedure is developed to determine the unknown inertial parameters
of the BPS [14]. With this method a dominant inertia I
d
is measured during experiments for each
particular robot joint. This dominant inertia is dened as the total inertia of a manipulator link,
or a group of links treated as a single link, at a certain manipulator conguration. To do this all
joints are locked except for the joint from which the dominant inertia is measured. The robot ma-
nipulator can now be considered as one rigid body at that specic conguration. The dominant
inertia measured at that particular joint consists of the total inertia contribution of all the links
chained to the joint at that conguration. Furthermore it is important to note that the method as-
sumes that viscous friction makes the dominant contribution to the total friction being present in
the joints. The procedure of measuring the dominant inertia by experimental means is as follows.
If all joints are locked except for the specic measured joint then the single joint dynamics can
be modeled as:
I
d
q +
f
( q) +g(q) =
u
(4.4)
where
f
( q) represents the friction present. At rst the identication procedure assumes viscous
friction to be dominantly present in the joints. So the following friction model is used:

f
( q) = f
v
q (4.5)
As control input
u
a normal PD control combined with gravity compensation is used:

u
= K
p
(q
d
q) +K
d
( q
d
q) + g(q) (4.6)
Where g(q) represents the approximated gravity compensation. If it is assumed that g(q) = g(q)
then the closed loop becomes:
I
d
q +f
v
q = K
p
(q
d
q) +K
d
( q
d
q) (4.7)
By setting the desired position and desired velocity to zero the resulting equation is:
q +
f
v
+K
d
I
d
q +
K
p
I
d
q = 0 (4.8)
If | K
d
| = f
v
an undamped second order system remains as:
q +
K
p
I
d
q = 0 (4.9)
For such a system the undamped natural frequency
n
can be calculated as follows:

n
=

K
p
I
d
T
n
=
2
_
K
p
I
d
(4.10)
16
Equation (4.10) can be rewritten as:
I
d
=
k
p
T
2
n
4
2
(4.11)
from which the dominant inertia for the specic conguration can be solved if T
n
is known from
experiments.
Thus a possible procedure to solve for the dominant inertia, based on the theory above, is as
follows:
put desired angle q
d
and desired velocity q
d
zero as reference input for the PD control
put K
p
zero and then manually tune K
d
to compensate for the viscous friction present so
that the joint acts undamped. This is done by giving it an initial push by hand to overcome
the static friction. Then, after releasing, the joint should maintain constant velocity.
give K
p
an positive value so that the system acts as an undamped mass-spring system.
give an initial angle and release.
the joint will now experience an undamped oscillation from where the undamped period
T
n
can be calculated by logging the angular feedback data.
Using this experiment several dominant inertias corresponding to several dierent robotic con-
gurations can be determined.
To solve for the unknown inertia elements in the BPS it must be noted that a dominant inertia
I
d
i
related to released joint i corresponds to element d
ii
from the inertia matrix D as in equation
(2.3). Otherwise stated I
d
i
= d
ii
. Such a element d
ii
is linear in some inertial elements of the
BPS. So element d
ii
can also be expressed as:
d
ii
() = f
1
()I
1
+f
2
()I
2
+ +f
m
()I
m
(4.12)
Where f
k
() denotes a specic formula dependent of , the robotic conguration, and I
k
denotes
the corresponding specic inertia element of the BPS.
To determine the m unknown inertial elements of the BPS the element d
ii
() should be
experimentally determined m times for linear independent congurations. Then the following
relations appear:
_

_
d
ii
(
1
)
d
ii
(
2
)
.
.
.
d
ii
(
m
)
_

_
=
_

_
f
1
(
1
) f
2
(
1
) . . . f
m
(
1
)
f
1
(
2
) f
2
(
2
) . . . f
m
(
2
)
.
.
.
.
.
.
f
1
(
m
) f
2
(
m
) . . . f
m
(
m
)
_

_
. .
B
_

_
I
1
I
2
.
.
.
I
m
_

_
(4.13)
From this the inertial parameters can be determined as they are the only unknowns. It should be
noted that if the m congurations are not linearly independent then matrix B becomes singular
and no solution is available.
17
In table 4.3 it is visible on which elements of the BPS the diagonal elements of D depend; one
if it is just dependent with that factor or () if there is relation that depends on the joint angles
and thus on the specic manipulator conguration. The bold elements are the only unknowns
though and thus the only inertial parameters that should be experimentally determined.
So, some of the inertial parameters where d
ii
() depends on are already known. If these
elements are expressed as I
k
then equation (4.13) can be rewritten as:
_

_
d
ii
(
1
)
d
ii
(
2
)
.
.
.
d
ii
(
k
1
)
.
.
.
d
ii
(
k
z
)
.
.
.
d
ii
(
m
)
_

_
=
_

_
f
1
(
1
) f
2
(
1
) . . . f
k
1
(
1
) . . . f
k
z
(
1
) . . . f
m
(
1
)
f
1
(
2
) f
2
(
2
) . . . f
k
1
(
2
) . . . f
k
z
(
2
) . . . f
m
(
2
)
.
.
.
.
.
.
f
1
(
k
1
) f
2
(
k
1
) . . . f
k
1
(
k
1
) . . . f
k
z
(
k
1
) . . . f
m
(
k
1
)
.
.
.
.
.
.
f
1
(
k
z
) f
2
(
k
z
) . . . f
k
1
(
k
z
) . . . f
k
z
(
k
z
) . . . f
m
(
k
z
)
.
.
.
.
.
.
f
1
(
m
) f
2
(
m
) . . . f
k
1
(
m
) . . . f
k
z
(
m
) . . . f
m
(
m
)
_

_
_

_
I
1
I
2
.
.
.
I
k
1
.
.
.
I
k
z
.
.
.
I
m
_

_
(4.14)
Because of this known inertial parameters the number of equations to solve for the unknown
inertial parameters can be reduced, and likewise the number of experiments necessary. Equation
(4.14) then becomes:
_

_
d
ii
(
1
)
d
ii
(
2
)
.
.
.
d
ii
(
m
)
_

_
f
k
1
(
1
) . . . f
k
z
(
1
)
f
k
1
(
2
) . . . f
k
z
(
2
)
.
.
.
f
k
1
(
m
) . . . f
k
z
(
m
)
_

_
_

_
I
k
1
.
.
.
I
k
z
_

_
=
_

_
f
1
(
1
) f
2
(
1
) . . . f
m
(
1
)
f
1
(
2
) f
2
(
2
) . . . f
m
(
2
)
.
.
.
.
.
.
f
1
(
m
) f
2
(
m
) . . . f
m
(
m
)
_

_
. .
Bm
_

_
I
1
I
2
.
.
.
I
m
_

_
(4.15)
Where the columns and rows corresponding to the known inertial parameters are removed at the
right side. Now the unknown inertial parameters can be solved. From table 4.1 it is obvious that
the experiments should start from the last link down to the rst because of the dependency on
the inertial parameters.
18
d
11
d
22
d
33
d
44
d
55
d
66
d
77
1. I
zzR
1
1
2. I
xxR
2
()
3. I
xy
2
()
4. I
xz
2
5. I
yz
2
6. I
zzR
2
7. mX
2
8. mY
R
2
9. I
xxR
3
() ()
10. I
xy
3
() ()
11. I
xz
3
()
12. I
yz
3
()
13. I
zzR
3
() 1
14. mX
3
()
15. mY
R
3
()
16. I
xxR
4
() () ()
17. I
xy
4
() () ()
18. I
xz
4
() ()
19. I
yz
4
() ()
20. I
zzR
4
() () 1
21. mX
4
() ()
22. mY
R
4
() ()
23. I
xxR
5
() () () ()
24. I
xy
5
() () () ()
25. I
xz
5
() () ()
26. I
yz
5
() () ()
27. I
zzR
5
() () () 1
28. mX
5
() () ()
29. mY
R
5
() () ()
30. I
xxR
6
() () () () ()
31. I
xy
6
() () () () ()
32. I
xz
6
() () () ()
33. I
yz
6
() () () ()
34. I
zzR
6
() () () () 1
35. mX
6
() () () ()
36. mY
6
() () () ()
37. I
xxR
7
() () () () () ()
38. I
xy
7
() () () () () ()
39. I
xz
7
() () () () ()
40. I
yz
7
() () () () ()
41. I
zz
7
() () () () () 1
42. mX
7
() () () () ()
43. mY
7
() () () () ()
Table 4.1: Inertial dependency of diagonal elements of D(q).
19
4.3.1 Experiment in practice
In theory (4.15) could be solved by using the exact inverse of B
n
and so by using only m experi-
ments. But, however in practice by using real experimental data this would give very inaccurate
results. To avoid this problem equation 4.15 should be rewritten as:
_

_
d
ii
(
1
)
d
ii
(
2
)
.
.
.
d
ii
(
w
)
_

_
f
k
1
(
1
) . . . f
k
z
(
1
)
f
k
1
(
2
) . . . f
k
z
(
2
)
.
.
.
f
k
1
(
w
) . . . f
k
z
(
w
)
_

_
_

_
I
k
1
.
.
.
I
k
z
_

_
=
_

_
f
1
(
1
) f
2
(
1
) . . . f
m
(
1
)
f
1
(
2
) f
2
(
2
) . . . f
m
(
2
)
.
.
.
.
.
.
f
1
(
w
) f
2
(
w
) . . . f
m
(
w
)
_

_
. .
Bw
_

_
I
1
I
2
.
.
.
I
m
_

_
(4.16)
with w >> m. So w dierent congurations, and thus w experiments, to solve for m unknowns.
Then the munknown inertia parameters can be solved with use of the pseudo inverse dened as:
(B
T
w
B
w
)
1
B
T
w
(4.17)
To guarantee the existence of such a pseudo inverse it must hold that: det(B
T
w
B
w
) > 0.
20
Chapter 5
Experimental results
5.1 Introduction
In Chapter 4 the theory on which the experimental determination of the inertial parameters is
based is addressed. In this Chapter experimental results of using that identication method will
be stated. First some practical problems encountered during experiments will be addressed in
section 5.2. The results of the original experiment as suggested in Chapter 4 are given in 5.3
There will be commented on and the section will be followed by two sections, section 5.4 and 5.5,
that reveal some new insights and critical aspects of the original method.
5.2 Practical problems
For a schematic representation of the experimental set-up and information about used hard-
ware/software and corresponding abbreviations the reader is referred to appendix (B).
Amplier gain settings
After some review it appears that the actual torque, T
act
, that results at the joint end diers from
the desired torque, T
des
, send by the source code. To search for the cause of this dierence
all steps in between are considered. In the original code there appears to be an unexplainable
correction gain present which is probably used before to correct for this dierence. This gain is
removed. Furthermore the ampliers contain a gain, A
g
, to convert their received voltage input
into a corresponding ampere output to the joint actuator. The gains originally present in the
ampliers seem to be randomly chosen.
To solve the problem the value of this gain, A
g
, is redetermined as follows. From the manu-
facturer the data as depicted in table 5.1 is known.
Table 5.1: Actuator specications supplied by manufacturer.
Joint 1, 2 3, 4 5, 6, 7
Rated torque T
r
[Nm] 5.3 2.0 0.36
Torque constant K
T
[
Nm
A
] 0.692 0.501 0.0831
21
The rated torque is chosen to be the maximum that can be demanded of the specic actuator. The
output range that the STG can deliver to the amplier is minus 10[V ] to plus 10[V ]. To utilize
this full range optimally 10[V ] received by the amplier should correspond to an ampere output
that will cause the actuator to deliver the rated torque. To ensure this the gain A
g
is set as:
A
g
= T
r
1
K
T
1
10
[
A
V
] (5.1)
With use of the values in table 5.1 equation (5.1) gives the following gains:
Joint 1, 2 3, 4 5, 6, 7
A
g
[
A
V
] 0.776 0.399 0.433
Table 5.2: Calculated amplier gains.
To ensure that T
des
= T
act
the source code is adapted correspondingly. After the desired torque
is calculated it is transformed to the corresponding voltage output, U
des
, and send by the STG to
the ampliers as:
U
des
= T
des
1
A
g
1
K
T
1
n
[V ] (5.2)
It is assumed that the manufacturer values are correct. Then it is easily veried that T
des
= T
act
by using the AMC software
1
and directly measure the current received by the actuator. From this
the output torque can easily be calculated as:
T
act
= i
measured
K
T
n [Nm] (5.3)
The complete scheme as implemented in the code and present in the physical parts is visible in
gure (5.1).
T
des
[Nm]
1
Ag
1
KT
1
n
U
des
[V ]
A
g
K
T
n
T
act
[Nm]
source code amplier actuator transmission
i
des
[A] T
mot
[Nm]
Figure 5.1: Schematic representation of gains.
Saturation settings
The output range of the STG as well as the allowable input range to the ampliers are 10[V ].
A common error in practice is the correct handling of desired control inputs that exceed the
allowable range. So here a |U
des
| > 10[V ]. It appears that the present algorithm in the source
1
Advanced Motion Controller: computer software provided by the manufacturer of the ampliers to adjust
and visualize amplier settings, measurements etcetera.
22
code does not suce correctly to handle exceeded limits. Though the algorithm should nicely
cap the exceeded control input to the maximum allowable input, the real torque output becomes
highly unexpected when limits are exceeded.
To cope with this problem a simple if loop is implemented in the source code after U
des
is
calculated . The loop just acts as a saturation before the limits are exceeded. An important thing
to notice however, is that for the saturation limits not the real limits of 10[V ] are used. The
limits are chosen more moderate. This has to be done because of an oset being present in the
ampliers. If the source code would exactly cap at the limits of 10[V ], so the STG can handle,
then the STG would nicely sends 10[V ] to the amplier. The receiving amplier can still not
handle this properly because of the oset being present that is added on to the receiving STG
signal. To be on the safe side and prevent exceeding the limits of the amplier, the source code
loop acts as a saturation with limits of 9.7[V ].
Saturated velocity
Another problem encountered during experiments is the joint velocity saturating before the max-
imum joint velocity, as supplied by the manufacturer, is reached. Originally the experimental
set-up was working with a power supply delivering 24[V ]. It appears that the following physical
balances limit the velocity before the maximum joint velocity is reached in case of using 24[V ] as
power supply:
U
pow
= K
T

m
+iR +L
d
dt
i (5.4)
The last term can be neglected because of L << R. Furthermore the current i is related to the
mechanical balance as:
K
T
i = (5.5)
So a certain required torque output requires a corresponding i related as in (5.5). Using this
required i (5.4) still has to be in balance so the velocity

will adapt accordingly to comply to
this. In case of U
pow
being 24[V ] this balances will result in

saturating before reaching the
maximum joint velocity as supplied by the manufacturer. To avoid this saturation problem the
simple solution is to use 48[V ] as power supply.
5.3 Results original procedure
To verify that the method as suggested in Chapter 4 indeed works, the experiment is rst applied
for only one conguration, namely the robot in upright position, and only on joint 1. So the
robotic manipulator is temporarily regarded as an 1 d.of. system. Then the experiment deter-
mines the inertia corresponding to that specic upright position for movement of joint 1, I
exp
.
Next this measured inertia can be used in trajectory control of joint 1 with the robot still in fully
upright position as:

u
= I
exp
q
d
+K
p
(q
d
q) +K
d
( q
d
q) (5.6)
Adding this feed-forward term should clearly improve the control performance if I
exp
is correct.
Like this the identication methodology can be easily veried by doing so. It must be noted
23
that the encountered phenomena as described below hold for all joints as is experienced during
experiments.
The rst step as described before in Chapter 4 is to manually tune K
d
in the PD-control to
compensate for the viscous friction present in the system. After slowly increasing K
d
(making
it more negative) the following phenomena appear to be the case:
| K
d
| < 12 results in a damped motion after releasing the joint with an initial velocity.
| K
d
| = 12 is the rst value that can keep the velocity going but only if the initial velocity
when released is close to the maximum velocity of the joint. If that is the case, then the
velocity will remain the maximum joint velocity. Else, in the case of a lower initial velocity,
the system still acts damped and the velocity settles at zero.
This last phenomenon is clearly not as expected by the theory. According to the model as in (4.8)
there are three possible situations. For clearness the model is recapitulated as:
q +
f
v
+K
d
I
d
q +
K
p
I
d
q = 0 (5.7)
Then the three possible situations are:
1. | K
d
| < f
v
, in which the system is undercompensated and all initial velocities will go to
zero after releasing.
2. | K
d
| = f
v
, in this case the viscous friction is exactly compensated and the system will
perfectly act as an undamped system. If that is the case then all initial velocities should
maintain constant after releasing.
3. | K
d
| > f
v
, in which the system is overcompensated and all initial velocities will keep
accelerating till the maximum joint velocity is reached.
Clearly none of this three situations correspond to the encountered phenomenon; some initial
velocities that completely damp to zero and initial velocities that result in maintaining maximum
joint velocity.
Hypothesis
The original procedure as suggested is fully based on the assumption that viscous friction domi-
nates the total friction being present in the joints. However, the rst results do not seem to sup-
port this assumption as clearly none of the three possible situations above ts the phenomenon.
A plausible explanation is that Coulomb friction seems to be signicantly present and modeling
only viscous friction compensation is insucient.
If this is indeed the case then the situation present is probably as in gure 5.2. So in the case of
Coulomb friction being signicantly present it is never possible to align the compensated friction
with the real friction if there is only compensated for viscous friction. Then there are two possible
situations that can occur. The compensated friction has no intersection point with the real friction
and is thus always smaller. This results in all initial velocities damp to zero. Another possibility
is that the compensated friction intersects with the real friction. This situation is visualized
in the gure by intersection point a. If this kind of situation is present then the encountered
phenomena can be explained as follows;
24
0
0
Velocity
F
r
i
c
t
i
o
n

a
real fric.
comp. fric.
Figure 5.2: Expected relation between real friction being present and compensated friction.
if the initial velocity after releasing is to the left side of intersection point a then the system
is undercompensated and the velocity will damp.
if the initial velocity after releasing is to the right side of intersection point a then the system
is overcompensated and will keep on acceleration till the maximumjoint velocity is reached.
Friction identication
To prove the stated hypothesis, i.e. Coulomb friction is signicantly present and therefore an
intersection point between the compensated and real friction is possible which causes the phe-
nomenon, a rough friction identication is done. The friction identication procedure used is as
follows; as a reference signal for the velocity a trapezoidal input V
des
is used, see gure 5.3. The
corresponding position reference signal is calculated. These two reference signals are used for
normal PD control.
On the horizontal parts of the trapezoidal the velocity is constant and thus the acceleration is
zero. For the parts where the acceleration is zero the following change of the original single joint
model
I q +
f
( q) +g(q) =
u
(5.8)
in to

f
( q) +g(q) =
u
(5.9)
holds. For the robot in upright position g(q) has no eect on following this desired trapezoidal
velocity pattern for joint 1, 3 and 5. For joint 2 and 4 however, the robot is brought in horizontal
position to cancel gravitational eects. Then the gravitational term can be canceled and (5.9)
becomes:

f
( q) =
u
(5.10)
From (5.10) it is clear that the control input force
u
is equal to the experienced friction force

f
. For one data point the mean of the input control torque, gure 5.4, on the horizontal parts
25
0 2 4 6 8 10 12 14 16 18 20
0.4
0.3
0.2
0.1
0
0.1
0.2
0.3
0.4
t[s]
V
e
l
o
c
i
t
y
[
r
a
d
/
s
]
V
des
V
Figure 5.3: Desired velocity with measured ve-
locity.
0 2 4 6 8 10 12 14 16 18 20
15
10
5
0
5
10
15
t[s]
T
o
r
q
u
e
i
n
p
u
t
[
N
m
]

u
Figure 5.4: Torque input
u
.een beetje plaat-
sopvulling
of the trapezoidal velocity pattern is taken. This force value represents the friction force for the
corresponding velocity as given by the trapezoidal input. With use of this method several data
points for dierent velocities are gathered to do a friction identication. The friction identication
consists of a least-square t through the gathered data points. The model used for this t is based
on Coulomb and viscous friction like described before in equation (2.21):

f
vc
(t) = f
v
q +f
c
sign( q) (5.11)
Here only the result of this friction identication for joint 1 is discussed. All other joint identi-
cations and corresponding identied values f
v
and f
c
are given in appendix (C).
0 0.2 0.4 0.6 0.8 1
0
2
4
6
8
10
12
14
Velocity [rad/s]
T
[
N
m
]
data point
real fric.
viscous |K
d
| = 12
Figure 5.5: Relation between identied real friction being present and compensated friction.
The identication result (dashed line) for joint 1 is visible in gure 5.5. To prove that the friction
as identied is indeed correct it is added to normal PD control to compensate for the friction. For
26
the same PD control values and the same reference signal the error reduction because of added
friction compensation is visible in gure 5.6
0 1 2 3 4 5 6
0.01
0.005
0
0.005
0.01
t[s]
E
r
r
o
r
[
r
a
d
]
PD
PD+fric.
Figure 5.6: Error reduction with use of identied friction compensation.
With a considerable error reduction being obviously present it is proven that the friction compen-
sation as added to the PD control is of good quality. The friction compensation is based on the
experimentally identied friction as visible in gure 5.5. Thus the experimentally procedure used
gives an identied friction model of good quality. Then, looking at the identied friction model
in gure 5.5 it is clear that the Coulomb friction adds a signicant part to the total friction. The
hypothesis stated earlier, that a Coulomb part being signicantly present is the cause to the phe-
nomenon, is correct because the behaviour matches with gure 5.5. The compensation as used
in the original procedure is also visible in the gure for K
d
= 12. This is the rst value that can
keep the velocity constant after releasing, though only when released at almost maximum joint
velocity (q
max
= 1.2[rad/s]). This can be explained with use of the graph. K
d
= 12 is the rst
friction compensation that has a slope steep enough to ensure that the intersection point between
the compensated friction and the real friction is within the possible velocity range. Indeed it is
correct that the velocity always acts undercompensated, so damped, if it is to the left side of this
intersection point and indeed it can only keep constant velocity at maximum velocity because
there the systems friction is overcompensated.
This result aects the original identication method as described in Chapter 4. To be able to get
an undamped system, and consequently being able to measure the undamped eigenfrequency,
the original model for the procedure as suggested in equation (4.4) has to be expanded with
Coulomb friction. Thus the new model becomes:
I
d
q +
f
vc
( q) +g(q) =
u
with
f
vc
( q) = f
v
q +f
c
sign( q) (5.12)
27
5.4 Results using identied friction
Now the friction present in the system is experimentally identied this result is used to make the
joint act as an undamped system using the following control:

u
=
f
vc
with
f
vc
=

f
v
q +

f
c
sign( q) (5.13)
Here
f
vc
represents the identied friction, with the identied values

f
v
and

f
c
that follow from
the least-square t through the data.
To verify that the system indeed acts undamped, joint 1 is released for dierent initial veloc-
ities initiated with use of trapezoidal control as used before to do friction identication. If the
velocity is on the constant part of the trapezoidal velocity trajectory then the controller used to
track the trapezoidal reference signal is stopped by switching to the controller of equation (5.13).
For an undamped system the velocity should remain nicely constant. The result is visible in g-
ure (5.7). It is clearly visible that the velocity does not remain constant as expected. The gure
gives not enough information to draw real conclusions about the systems behaviour.
0 2 4 6 8 10 12 14
0
0.1
0.2
0.3
0.4
0.5
0.6
0.7
0.8
0.9
1
t[s]
V
e
l
o
c
i
t
y
[
r
a
d
/
s
]
Figure 5.7: Velocity behaviour using identied friction compensation.
Because it is obvious that practical results always dier from the ideal case, and on rst sight
there does not seem to be a better alternative then the friction identication done before, that is
to guarantee an undamped system, the identication procedure is continued to nd the inertia
value for joint 1 in upright position. Therefore K
p
is added to equation (5.13). The new control
input thus becomes:

u
= K
p
q +
f
vc
with
f
vc
=

f
v
q +

f
c
sign( q) (5.14)
If the systeemis indeed undamped then a constant oscillation wil remain with an amplitude equal
to an initial angle oset. The result is visible in gure(5.8). Although the system indeed ends up
in a constant undamped oscillation it is remarkable that the amplitude for dierent initial osets
ends up to be the same. So the system behaves as if there is a stable limit cycle present. In the
next section this will be conrmed.
28
2 3 4 5 6 7 8 9 10 11 12 13
0.1
0.05
0
0.05
0.1
t[s]
A
n
g
l
e
[
r
a
d
]
Figure 5.8: Angular motion.
To be conclusive that the improved inertia identication method, i.e. with identied friction
instead of manually tuned friction compensation, does not give satisfactory results the inertia is
calculated for dierent K
p
. In case of an undamped system the period should change accordingly
to K
p
so that equation (4.11) always gives the same inertia. In gure (5.9) the calculated inertia
is plotted against the K
p
used. Each experiment with a dierent value for K
p
is performed twice
to verify its repeatability. A horizontal line is expected but the result appears to be that a dierent
K
p
does aect the calculated inertia outcome. It is clear that this experimental result is not
satisfactory and that the applied improved identication procedure is not suitable.
0 200 400 600 800
1
1.5
2
2.5
3
3.5
4
4.5
Kp [N/m]
I
[
K
g
m
2
]
Figure 5.9: Inuence of K
p
.
29
5.5 Friction inuence in practice
In the previous section it is shown that, despite the usage of an experimentally identied friction
model, the friction compensation still does not result in satisfying a fully undamped system. Even
if the friction is identied by experimental means it is never perfectly identical to the real friction.
The original inertia identication method as described in Chapter 4 is based on the assumption
that this dierence between compensated and real friction does not eect the systems behaviour
signicantly. It expects the system to approach undamped behaviour good enough to gain useful
experimental data. However in practice this seems not to be the case. The imperfection in friction
compensation results in a specic corresponding system behaviour.
0
0
Velocity
F
r
i
c
t
i
o
n

s
1
s
2
real fric.
comp. fric.
Figure 5.10: Friction compensation v.s. real
friction present; situation s
1
and s
2
.
0
0
Velocity
F
r
i
c
t
i
o
n

a
s
4
s
3
v
a
real fric.
comp. fric.
Figure 5.11: Friction compensation v.s. real
friction present; situation s
3
and s
4
.
Because the emphasis of this report is not on friction identication but on the usage of the
identication method of Chapter 4 for identication of system parameters, our analysis of the
inuence of real friction versus compensated friction is restricted to the widely excepted and
used friction model that consists of a viscous and Coulomb part as in 2.21. For a much more
thorough research to the eects of imperfect friction compensation in mechanical systems the
reader is referred to [20].
If the friction models are only based on viscous and Coulomb then there are four dierent
situations of imperfect friction compensation. The four situations are depicted in gures 5.10
and 5.11. Situation s
1
and s
2
are straightforward in their behaviour. In the case of situation s
1
the system is always undercompensated for the real friction and therefore the initial velocity will
always, regardless of its value, damp to zero. In the case of situation s
2
the system is always
overcompensated for the real friction and therefore the initial velocity, again regardless of its
value, will always keep increasing till the maximum possible joint velocity is reached.
More interesting are situation s
3
and s
4
. They dier from the other two in having an inter-
section point (a) with the real system friction. Situation s
3
is almost similar as experienced with
the original friction compensation; only now inclusive some small Coulomb part. The behaviour
for this two situations is expected to be as follows; in the case of situation s
3
the system will act
damped if the intial velocity is in the range to the left of intersection point a (v < v
a
). Other-
wise, if the initial velocity is in the right range of intersection point a (v > v
a
), the system will
be overcompensated and velocity will increase till the maximum joint velocity is reached. In the
30
case of situation s
4
the system will always settle at the velocity that corresponds to the intersec-
tion point a, v
a
. If the initial velocity is in the left range of intersection point a then the system
is overcompensated and the velocity will increase. This increasing stops if the velocity crosses
the velocity that corresponds to the intersection point. After the intersection point the system is
undercompensated and the velocity will experience damping. So the velocity will decrease again
till the intersection point. If the initial velocity is in the right range of the intersection point the
opposite happens, rst decreasing velocity and then increasing velocity till it is settled in the ve-
locity again that corresponds to the intersection point. In the previous section a limit cycle was
encountered. With use of simulations with a Matlab Simulink model this observed behaviour is
conrmed when a stiness K
p
is added; for the case of situation s
4
there appears to be a stable
limit cycle. Physically this can be explained as follows; for low velocities the system gains en-
ergy because of the control overcompensation and at higher velocities the system is subjected to
energy loss because of friction undercompensation, so damping.
The described behaviour above, the settled end velocity for dierent initial velocities, for the
dierent situations can be graphically represented as visible in gures (5.12) and (5.13).
0
0
s
1
s
2
q
max
v
a
Initial velocity
S
e
t
t
l
e
v
e
l
o
c
i
t
y

Figure 5.12: Settled velocity versus initial ve-

locity; situation s
1
and s
2
.
0
0
s
3
s
4
q
max
v
a
v
a
Initial velocity
S
e
t
t
l
e
v
e
l
o
c
i
t
y

Figure 5.13: Settled velocity versus initial ve-

locity; situation s
3
and s
4
.
That the behaviour of the system for situation s
3
and s
4
is indeed as described before is rst
veried using experiments. Based on the identied friction for joint 1 as visible in gure (5.5)
two dierent friction compensations are used that are sure to guarantee situation s
3
and s
4
. The
experimental results for situation s
3
and s
4
are visible in gure (5.14) and gure (5.15) respectively.
They conrm the expected behaviour. In gure (5.14) velocity v
a
, and thus the corresponding
intersection point a is somewhere between 0.4[rad/s] and 0.7[rad/s]. In gure (5.15) the settled
velocity v
a
is approximately 0.5[rad/s].
All the system behaviour encountered during experiments for system identication, based on
the theory of Chapter 4, t one of the described behaviours above corresponding to one of the
four dierent situations. This shows that even a slight dierence between compensated and real
friction gives a system behaviour that can not be regarded as fully undamped.
31
0 0.5 1 1.5 2 2.5 3 3.5 4 4.5 5
0.2
0
0.2
0.4
0.6
0.8
1
t[s]
V
e
l
o
c
i
t
y
[
r
a
d
/
s
]
Figure 5.14: Velocity behaviour for situation s
3
.
0 2 4 6 8 10 12 14
0
0.2
0.4
0.6
0.8
1
t[s]
V
e
l
o
c
i
t
y
[
r
a
d
/
s
]
Figure 5.15: Velocity behaviour for situation s
4
.
32
Chapter 6
Conclusions and Recommendations
6.1 Introduction
The research objective is recapitulated as; describe the complete dynamics of the Mitsubishi
PA-10 robotic manipulator and identify all unknown system parameters, used in the obtained
dynamic model, by experimental means. The applied identication method, as proposed by
SIMtech, is explained in Chapter 4. Experimental results using the identication method are
shown in Chapter 5. These results show that the identication method in practical applications
does not give the results that correspond to the theoretical foundations. In the next paragraph
research results and issues regarding the identication method will be discussed followed by
recommendations for future research.
6.2 Contribution
The complete symbolically dynamical model is derived for the Mitshubishi PA-10 manip-
ulator. This model is reduced in complexity to make it suitable to use in identication
procedures. This reduction is done in a well-structured way that gives the base parameter
set of 43 inertia parameters. Furthermore some other assumptions led to further reduc-
tion of the mathematical model (13 unknowns left) to diminish the burden of experimental
inertia identication.
It is shown that the original identication procedure, based on the assumption of real fric-
tion being dominated by viscous friction, is not suitable for the PA-10. Using friction iden-
tication by experimental means it is proven that Coulomb friction is signicantly present
in the joints.
Relating the previous observation that only viscous friction compensation does not suce,
the real friction present is experimentally identied for all joints. It is shown that the gained
identied friction model gives substantial improvement if used as friction compensation in
trajectory control.
Furthermore it is shown that even by using an experimentally identied friction model for
friction compensation, the proposed identication method is too sensitive to imperfections
in correct friction compensation. Four dierent situations of imperfect friction compensa-
tion are considered and their expected behaviour is conrmed by experimental means.
33
The general conclusion of the research is that the identication method based on nding an
undamped system is very nice in theory but not suitable for practical applications. The reason for
this is that in practice researchers can never succeed in achieving perfect compensation for the
real friction being present. So this unavoidable imperfection will always result in joint behaviour
not equal to a nicely undamped system. The sensitivity of the identication method for this im-
perfection is that large that identied results are worthless. Therefore the identication method
should be regarded impractical and it can not be used to identify system parameters as supported
by erroneous identication results.
6.3 Recommendations
The main obstacle present for succeeding in the implementation of model-based control is iden-
tifying the full dynamics of the robotic manipulator. It is shown that the originally proposed iden-
tication method is not suitable in practice to achieve this. If the general objective of SIMtech is
kept in mind, implementing force control using model-based control, then the following recom-
mendations for future research are done:
Worldwide there has been a lot of research on the topic of system parameter identica-
tion. In this report no attention is spent to investigating alternatives. In Chapter 3 some
of the common identication methods are mentioned though. Based on experimental re-
sults shown in literature, making use of an identication method using derived velocity
information in combination with external torque/force information measured on a bed-
plate supporting the robotic manipulator seems a promising alternative. Investing time
in nding a suitable alternative is recommended. For most identication methods the dy-
namic model, including the derivation of the base parameter set and other simplications,
as described in Chapter 2 is still suitable to use.
The inuence of friction compensation has been topic of many researches. Some re-
searches focus specically on the eects of under- and overcompensation on system be-
haviour [20]. Maybe there are possibilities to use all this available knowledge to adapt the
original identication procedure and making it less sensitive for the eects of frictions.
Another possibility to avoid the diculty of system identication is the use of adaptive
robot control. By using an adaptive model system parameters will converge during control
experiments and no system identication in advance of model-based control is required
[21].
34
References
[1] John J. Craig, "Introduction to Robotics: mechanics and control", 2nd edition, Addison-
Wesley, 1989
[2] L. Sciavicco and B. Siciliano, "Modeling and control of robot manipulators", 2nd edition,
Springer,
[3] Mark. W. Spong, Seth Hutchinson and M. Vidyasagar, "Robot Modeling and Control", Wiley,
2006
[4] P.I. Corke, "A Robotics Toolbox for Matlab", IEEE Robotics and Automation Magazine, vol.3,
no.1, pp.24-32, 1996
[5] Mitsubishi Heavy Industries LTD, "PA-10 portable general purpose intelligent arm, operating
manual", 9100025 rev1.
[6] Jan Swevers, Walter Verdonck and Joris de Schutter, "Dynamic Model Identication for In-
dustrial Robots", IEEE Control Systems Magazine, oct., pp.58-71, 2007
[7] G. Antonelli, F. Caccavale and P. Chiacchio, "A systematic procedure for the identication
of dynamic parameters of robot manipulators", Robotica, vol. 17, pp.427-435, Cambridge
University Press, 1999
[8] Fariborz Behi and Delbert Tesar, "Parametric Identication for Industrial Manipulators Using
Experimental Modal analysis", IEEE Transactions on Robotics and Automation, vol. 7, no.5,
pp.642-652, October 1991
[9] P.O. Vandanjon, M. Gautier and P. Desbats, "Identication of Robots Inertial Parameters by
Means of Spectrum Analysis", IEEE International Conference on Robotics and Automation,
pp. 3033-3038, 1995
[10] M. Gautier, "A Comparison of Filtered Models for Dynamic Identication of Robots", Pro-
ceedings of the 35th Conference on Decision and Control, pp.875-880,Kobe, Japan, 1996
[11] M. Gautier, "Dynamic Identication of Robots with Power Model", Proceedings of the 1997
IEEE International on Robotics and Automation, pp.1922-1927, Albuquerque, New Mexico,
April 1997
[12] J. Swevers et. all, "Experimental identication of robot dynamics for control", Proceedings
of the 2000 IEEE International Conference on Robotics and Automation, pp.241-246, San
Francisco, CA, April 2000
[13] Peter I. Corke and Brian Armstrong-Helouvry, "A Search for Consensus Among Model Pa-
rameters Reported for the PUMA 560 Robot", IEEE, pp.1608-1613, 1994
[14] Rodrigo S. Jamisola, JR., "Full Dynamics Identication and Control of PUMA 560 and Mit-
subishi PA-10 robots", thesis produced at National University Singapore, ME department,
2001
35
[15] M. Gautier and W. Khalil, "A direct determination of minimuminertial parameters of robots",
IEEE, 1988
[16] M. Gautier and W. Khalil, "Direct Calculation of Minimum Set of Inertial Parameters of Serial
Robots", IEEE Transactions on Robotics and Automation, vol. 6, no.3, pp.368-373, June 1990
[17] W. Khalil and F. Bennis, "Comments on Direct Calculation of Minimum Set of Inertial Pa-
rameters of Serial Robots", IEEE Transactions on Robotics and Automation, vol. 10, no.1,
pp.78-79, February 1994
[18] Nikolas A. Bompos et all., "Modeling full identication and Control of the Mitsubishi PA-10
Robot Arm", Control Systems Lab, School of Mechanical Eng. National Technical University
of Athens, publication date unknown
[19] Christopher W. Kennedy et all., "Modeling and Control of the Mitsubishi PA-10 Robot Arm
Harmonic Drive System", IEEE Transactions on Mechatronics, vol. 10, no.3, pp.263-274, June
2005
[20] D. Putra, H. Nijmeijer and N. van de Wouw, "Analysis of undercompensation and overcom-
pensation of friction in 1DOF mechanical systems", Automatica, no.43, pp.1387-1394, Elsevier,
June 2007
[21] H. Berghuis, H. Roebbers and H. Nijmeijer, "Experimental Comparison of Parameter Esti-
mation Methods in Adaptive Robot Control", Automatica, vol. 31, no.9, pp.1275-1285, Elsevier,
1995
36
Appendix A
Equations Base Parameter Set
Link 7:
I
xxR
7
= I
xx
7
I
yy
7
I
xxR
6
= I
xx
6
+I
yy
7
2d
7
mZ
7
+d
2
7
m
7
I
zzR
6
= I
zz
6
+I
yy
7
2d
7
mZ
7
+d
2
7
m
7
mY
R
6
= mY
6
mZ
7
d
7
m
7
m
R
6
= m
6
+m
7
_

_
(A.1)
So the base parameters of link 7 are: I
xxR
7
, I
xy
7
, I
xz
7
, I
yz
7
, I
zz
7
, mX
7
, mY
7
.
Link 6:
I
xxR
6
= I
xxR
6
I
yy
6
I
xxR
5
= I
xx
5
+I
yy
6
I
zzR
5
= I
zz
5
+I
yy
6
mY
R
5
= mY
5
+mZ
6
m
R
5
= m
5
+m
6
_

_
(A.2)
So the base parameters of link 6 are: I
xxR
6
, I
xy
6
, I
xz
6
, I
yz
6
, I
zzR
6
, mX
6
, mY
R
6
.
Link 5:
I
xxR
5
= I
xxR
5
I
yy
5
= I
xx
5
+I
yy
6
I
yy
5
I
xxR
4
= I
xx
4
+I
yy
5
+2d
5
mZ
5
+d
2
5
m
R
5
I
zzR
4
= I
zz
4
+I
yy
5
+2d
5
mZ
5
+d
2
5
m
R
5
mY
R
4
= mY
4
mZ
5
d
5
m
R
5
m
R
4
= m
4
+m
R
5
= m
4
+m
5
+m
6
_

_
(A.3)
So the base parameters of link 5 are: I
xxR
5
, I
xy
5
, I
xz
5
, I
yz
5
, I
zzR
5
, mX
5
, mY
R
5
.
37
Link 4:
I
xxR
4
= I
xxR
4
I
yy
4
= I
xx
4
+I
yy
5
+ 2d
5
mZ
5
+d
2
5
m
R
5
I
yy
4
I
xxR
3
= I
xx
3
+I
yy
4
I
zzR
3
= I
zz
3
+I
yy
4
mY
R
3
= mY
3
+mZ
4
m
R
3
= m
3
+m
R
4
= m
3
+m
4
+m
5
+m
6
_

_
(A.4)
So the base parameters of link 4 are: I
xxR
4
, I
xy
4
, I
xz
4
, I
yz
4
, I
zzR
4
, mX
4
, mY
R
4
.
Link 3:
I
xxR
3
= I
xxR
3
I
yy
3
= I
xx
3
+I
yy
4
I
yy
3
I
xxR
2
= I
xx
2
+I
yy
3
+2d
3
mZ
3
+d
2
3
m
R
3
I
zzR
2
= I
zz
2
+I
yy
3
+2d
3
mZ
3
+d
2
3
m
R
3
mY
R
2
= mY
2
mZ
3
d
3
m
R
3
_

_
(A.5)
So the base parameters of link 3 are: I
xxR
3
, I
xy
3
, I
xz
3
, I
yz
3
, I
zzR
3
, mX
3
, mY
R
3
.
Link 2:
I
xxR
2
= I
xxR
2
I
yy
2
= I
xx
2
+I
yy
3
+ 2d
3
mZ
3
+d
2
3
m
R
3
I
yy
2
I
zzR
1
= I
zz
1
+I
yy
2
_
(A.6)
So the base parameters of link 2 are: I
xxR
2
, I
xy
2
, I
xz
2
, I
yz
2
, I
zzR
2
, mX
2
, mY
R
2
.
Link 1:
Because of Theorem 3 the minimum parameter of link 1 is I
zzR
1
.
So this regrouping of inertial parameters reduces the number of inertial parameters from 60
to 43. (22 if you would not take the diagonal terms of the inertia matrices, so the products of
inertia (for example I
xz
i
etc.) in to account).
Theorem 7: minimum base parameters set n 7 4 = 45. Here 43 because: 2 because axis
z
1
is equal to the gravity vector.
38
Appendix B
Experimental setup
SOURCE CODE
2
STG
3
encoder pulses
reference signal 10[V ]
AMPLIFIER
4
ACTUATOR
5
drive signal [A] resolver feedback [V ]
COMPUTER SYSTEM
1
ROBOTIC MANIPULATOR
6
Figure B.1: Experimental set-up.
A schematic representation of the experimental set-up is visible in gure B.1. Technical spec-
ications of the computer system used for control and measurements and some other useful
technical data can be found below.
1. Computer system: Microsoft Windows XP
in combination with Realtime Extension for Windows (RTX)
sample rate: 1000[Khz]
2. Source code: c-code build with Microsoft Visual Studio 6.0,
GUI build with Microsoft Visual Basic 6.0
3. Server To Go (STG): model 2-8 axis ISA servo i/o card
4. Amplier: Advanced Motion Controller (AMC),
model Digiex; digital servo drive CDD005BB
5. Actuator: See technical specications reference [5], chapter 5 pg.13-14
6. Robotic manipulator: 7 d.o.f. robotic manipulator, see gures (1.1), (2.1).
See technical specications reference [5]
39
Appendix C
Friction identications
1
0.6 0.4 0.2 0 0.2 0.4 0.6
30
20
10
0
10
20
30
v [rad/s]
T
[
N
m
]
Figure C.1: Friction identication joint 2.
2 1.5 1 0.5 0 0.5 1 1.5 2
8
6
4
2
0
2
4
6
8
v [rad/s]
T
[
N
m
]
Figure C.2: Friction identication joint 3.
2 1.5 1 0.5 0 0.5 1 1.5 2
15
10
5
0
5
10
15
v [rad/s]
T
[
N
m
]
Figure C.3: Friction identication joint 4.
2 1.5 1 0.5 0 0.5 1 1.5 2
8
6
4
2
0
2
4
6
8
v [rad/s]
T
[
N
m
]
Figure C.4: Friction identication joint 5.
1
During the research period joint 6 was broke so no results of joint 6 are available.
40
2 1.5 1 0.5 0 0.5 1 1.5 2
2
1.5
1
0.5
0
0.5
1
1.5
v [rad/s]
T
[
N
m
]
Figure C.5: Friction identication joint 7.
The friction compensations are tted with the formula of equation (2.21) as:

f
vc
(t) = f
v
n
q +f
c
n
sign( q) q < 0 (C.1)

f
vc
(t) = f
v
p
q +f
c
p
sign( q) q > 0 (C.2)
The corresponding tted parameters are visible in table C.
joint f
v
p
f
c
p
f
v
n
f
c
n
2 3.57 14.26 7.64 18.25
3 1.64 3.27 1.70 4.15
4 2.58 6.93 2.89 4.95
5 1.61 3.75 1.64 3.64
6 - - - -
7 0.27 0.82 0.21 1.23
Table C.1: Friction compensation t parameters.
41
..
42