Projet COROUSSO Livrable n°1.4 Modèles des interactions dynamiques

INTRODUCTION

Currently, robots become more and more popular for a variety of technological processes, including highspeed precision machining. For this process, external loading caused by the machining force is applied on the robot tool. This force is generated by the interaction between the tool mounted on the robot end-effector and the workpiece during the material removal [START_REF] Dépincé | Active integration of tool deflection effects in end milling. Part 2. Compensation of tool deflection[END_REF]. It is a contact force and it is distributed along the affected area of the tool cutting part. To evaluate the influence and to analyze the robot behavior while machining, the cutting force should be defined either experimentally or using accurate mathematical model.

To evaluate the force caused by interaction between the tool and the workpiece, two approaches can be used. The static approach allows computing the average cutting force without any consideration of dynamic aspect in machining system. This force serves as an external loading of the robot. This approach is widely used in analysis of conventional machining processes using CNC machines [START_REF] Altintas | Manufacturing automation, metal cutting mechanics, machine tool vibrations and CNC design[END_REF], where the stiffness is high. In contrast, robots have relatively low structural stiffness. For this reason, in the case of robotic-based machining, an additional source of dynamic displacements of the end-effector with respect to the desired trajectory induced by robot compliance may arise. Such behavior leads to the variable contact between the machining tool and the workpiece. Thus, the generated contact force depends on the current position of the robot end-effector on the trajectory. Consequently, the cutting force cannot be evaluated correctly using the static approach. In this case, the dynamic approach, which will be used in the paper, is required. It is based on computing of the force at each instant of machining process that defines loading of the robot for the next instant of processing. As a result, the dynamic aspect of robot motion under such variable cutting force can be examined for whole process.

Usually, in the robot-based machining this force causes essential deflections that decrease the quality of the final product. The problem of the robot error compensation can be solved in two ways that differ in degree of modification of the robot control software: (a) by modification of the manipulator model, which better suits to the real manipulator and is used by the robot controller (in simple case, it can be limited by tuning of the nominal manipulator model, but may also involve essential model enhancement by introducing additional parameters, if it is allowed by a robot manufacturer); (b) by modification of the robot control program that defines the prescribed trajectory in Cartesian space (here, using relevant error model, the input trajectory is generated in such way that under the loading the output trajectory coincides with the desired one, while input trajectory differs from the target one).

Moreover, with regard to the robot-based machining, there is a solution that does not require force/torque measurements or computations [START_REF] Dépincé | Active integration of tool deflection effects in end milling. Part 2. Compensation of tool deflection[END_REF], where the target trajectory for the robot controller is modified by applying the "mirror" technique. An evident advantage of this technique is its applicability to the compensation of all types of the robot errors, including geometrical and compliance ones. However, this approach requires carrying out additional preliminary experiments which are quite expensive. So, it is suitable

COROUSSO

for the large-scale production only. Another compensation methodology has been proposed by Eastwood and Webb [START_REF] Eastwood | A gravitational deflection compensation strategy for HPKMs[END_REF] that was used for gravitational deflection compensation for hybrid parallel kinematic machines. This paper focuses on the modification of control program that is considered to be more realistic in practice. This approach requires also accurate stiffness model of the manipulator. From point of view of stiffness analysis, the external and forces directly influence on the manipulator equilibrium configuration and, accordingly, may modify the stiffness properties. So, they must be undoubtedly taken into account while developing the stiffness model. However, in most of the related works the Cartesian stiffness matrix has been computed for the nominal configuration [START_REF] Alici | Enhanced stiffness modeling, identification and characterization for robot manipulators[END_REF][START_REF] Kao | Conservative Congruence Transformation for Joint and Cartesian Stiffness Matrices of Robotic Hands and Fingers[END_REF]. Such approach is suitable for the case of small deflections only. For the opposite case, the most important results have been obtained in [START_REF] Kövecses | The stiffness matrix in elastically articulated rigid-body systems[END_REF][START_REF] Pashkevich | Enhanced stiffness modeling of manipulators with passive joints[END_REF][START_REF] Tyapin | Kinematic and elastostatic design optimization of the 3-DOF Gantry-Tau parallel kinamatic manipulator[END_REF], which deal with the stiffness analysis of manipulators under the end-point loading.

Thus, to compensate errors caused by the machining process, it is required to have an accurate stiffness model and precise cutting force model. In contrast to the previous works, the compliance error compensation technique presented in this work is based on the non-linear stiffness model of the manipulator [START_REF] Pashkevich | Enhanced stiffness modeling of manipulators with passive joints[END_REF] and dynamic model of technological process that generates the cutting force.

COROUSSO

PROBLEM STATEMENT

For the compliance errors, the compensation technique must rely on two components. The first of them describes distribution of the stiffness properties throughout the workspace and is defined by the stiffness matrix as a function of the joint coordinates. The second component describes the forces/torques acting on the end-effector while the manipulator is performing its machining task (manipulator loading).

The stiffness matrix required for the compliance errors compensation highly depends on the robot configuration and essentially varies throughout the workspace. From general point of view, full-scale compensation of the compliance errors requires essential revision of the manipulator model embedded in the robot controller. In fact, instead of conventional geometrical model that provides inverse/direct coordinate transformations from the joint to Cartesian spaces and vice versa, here it is necessary to employ the so-called kinetostatic model [START_REF] Su | A Polynomial Homotopy Formulation of the Inverse Static Analyses of Planar Compliant Mechanisms[END_REF]. It is essentially more complicated than the geometrical model and requires rather intensive computations that are presented in Section 3.

The dynamic behavior of the robot under the loading F caused by technological process can be described as In general, the cutting force F c has a nonlinear nature and depends on many factors such as cutting conditions, properties of workpiece material and tool cutting part, etc [START_REF] Ritou | A new versatile in-process monitoring system for milling[END_REF]. But, for given tool/workpiece combination, the force F c could be approximated as a function of an uncut chip thickness h, which represents the desired thickness to cut at each instant of machining.

C C C + + = M δt C δt K δt F && & (1) where
Hence, to reduce the errors caused by the cutting forces in the robotic-based machining it is required to obtain an accurate elastostatic model of the robot and elastodynamic model of the machining process. These problems are addressed in the following sections taking into account some particularities of the considered application (robotic-based milling).

Elastostatic model

Elastostatic model of a serial robot is usually defined by its Cartesian stiffness matrix, which should be computed in the neighborhood of loaded configuration. Let us propose numerical technique for computing static equilibrium configuration for a general type of serial manipulator. Such manipulator may be approximated as a set of rigid links and virtual joints, which take into account elastostatic properties (Fig. 1). Since the link weight of serial robots is not negligible, it is reasonable to decompose it into two parts (based on the link mass centre) and apply them to the both ends of the link. All this loadings will be aggregated in a vector [ ]

1 ... n = G G G
, where i G is the loading applied to the i-th node-point. Besides, it is assumed that the external loading F (caused by the interaction of the tool and the workpiece) is applied to the robot endeffector.

Following the principle of virtual work, the work of external forces , G F is equal to the work of internal forces θ τ caused by displacement of the virtual springs δθ ( )

T T T θ 1 δ δ δ n j j j= ⋅ + ⋅ = ⋅ ∑ G t F t τ θ (2)
where the virtual displacements δ j t can be computed from the linearized geometrical model derived from

( ) θ δ δ , 1.. j j j n = = t J θ
, which includes the Jacobian matrices ( )

( ) θ , j j
= ∂ ∂ J g q θ θ with respect to the virtual joint coordinates.

So, expression (2) can be rewritten as

( ) ( ) T ( ) T ( ) T θ θ θ 1 δ δ δ n j n j j= ⋅ ⋅ + ⋅ ⋅ = ⋅ ∑ G J θ F J θ τ θ (3)
which has to be satisfied for any variation of δθ . It means that the terms regrouping the variables δθ have the coefficients equal to zero. Hence the force balance equations can be written as

( )T ( )T θ θ θ 1 n j n j j = = ⋅ + ⋅ ∑ τ J G J F (4)
These equations can be re-written in block-matrix form as

(G)T (F)T θ θ θ = ⋅ + ⋅ τ J G J F (5)
where

(F) ( ) θ θ n = J J , (G) (1) ( ) θ θ θ T T T ... n ⎡ ⎤ = ⎣ ⎦ J J J , T T T 1 ... n ⎡ ⎤ = ⎣ ⎦ G G G
. Finally, taking into account the virtual spring

reaction θ θ = ⋅ τ K θ , where (
)

1 n θ θ θ ,..., diag = K K K
, the desired static equilibrium equations can be presented as

(G)T (F)T θ θ θ ⋅ + ⋅ = ⋅ J G J F K θ (6)
To obtain a relation between the external loading F and internal coordinates of the kinematic chain θ corresponding to the static equilibrium, equations ( 6) should be solved either for different given values of F or for different given values of t . Let us solve the static equilibrium equations with respect to the manipulator configuration θ and the external loading F for given end-effector position ( ) = t g θ and the function of auxiliary-loadings ( )

G θ ( ) ( ) (G)T (F)T θ θ θ ; ; ⋅ = + = = K θ J G J F t g θ G G θ (7)
where the unknown variables are ( )

, θ F .
Since usually this system has no analytical solution, iterative numerical technique can be applied. So, the kinematic equations may be linearized in the neighborhood of the current configuration i θ ( ) ( ) ( )

(F) θ 1 1
;

i i i i i + + = + ⋅ - t g θ J θ θ θ (8)
where the subscript 'i' indicates the iteration number and the changes in Jacobians (G) (F) θ θ , J J and the auxiliary loadings G are assumed to be negligible from iteration to iteration. Correspondingly, the static equilibrium equations in the neighborhood of i θ may be rewritten as

(G)T (F)T θ 1 θ 1 θ i i + + ⋅ + ⋅ = ⋅ J G J F K θ (9)
Thus, combining ( 8), ( 9) and expression for

1 (G)T (F)T θ θ θ ( ) - = ⋅ + ⋅ θ K J G J
F , the unknown variables F and θ can be computed using following iterative scheme

( ) ( ) ( ) ( ) θ θ θ 1 1 1 1 i i i i i i i i θ θ - - - - + + + + = ⋅ ⋅ - + - = ⋅ + ⋅ F J K J t g θ J θ J K J G θ K J G J F ( 10
)
The proposed algorithm allows us to compute the static equilibrium configuration for the serial robot under external loadings applied to any point of the manipulator and the loading from the technological process.

Stiffness matrix

In order to obtain the Cartesian stiffness matrix, let us linearize the force-deflection relation in the neighborhood of the equilibrium. Following this approach, two equilibriums that correspond to the manipulator state variables ( , , ) F θ t and ( δ ,

δ , δ ) + + +
F F θ θ t t should be considered simultaneously. Here, notations δF , δt define small increments of the external loading and relevant displacement of the end-point. Finally, the static equilibrium equations may be written as ( )

(G)T (F)T θ θ θ ; = ⋅ = ⋅ + ⋅ t g θ K θ J G J F (11) and ( ) ( ) ( ) ( ) ( ) (
)

T T (G) (G) (F) (F) θ θ θ θ θ δ δ δ δ δ δ δ + = + ⋅ + = + ⋅ + + + ⋅ + t t g θ θ K θ θ J J G G J J F F (12) where θ , , , , t F G
K θ are assumed to be known.

After linearization of the function ( ) g θ in the neighborhood of the loaded equilibrium, the system (11), ( 12) is reduced to equations

(F) θ (G) (G) (F) (F) θ θ θ θ θ δ δ δ δ δ δ δ = ⋅ = + + + t J θ K θ J G J G J F J F (13)
which defines the desired linear relations between δt and δF . In this system, small variations of Jacobians may be expressed via the second order derivatives (F)

(F) θ θ θ δ δ = ⋅ J H θ , (G) (G) θ θ θ δ δ = ⋅ J H θ , where (G) 2 (F) 2 θθ θθ 1 2 2 ; T T j n j j= = ∂ ∂ =∂ ∂ ∑ H g G θ H g F θ (14)
Also, the auxiliary loading G may be computed via the first order derivatives as

δ δ = ∂ ∂ ⋅ G G θ θ
Further, let us introduce additional notation

(F) (G) (G)T θθ θθ θθ θ = + + ⋅∂ ∂ H H H J G θ (15)
which allows us to present system (13) in the form

(F) θ (F)T θ θ θ θ δ δ δ ⎡ ⎤ ⎡ ⎤ ⎡ ⎤ = ⋅ ⎢ ⎥ ⎢ ⎥ ⎢ ⎥ -+ ⎣ ⎦ ⎣ ⎦ ⎣ ⎦ 0 J t F 0 θ J K H (16)
So, the desired Cartesian stiffness matrices C K can be computed as ( )

1 (F) 1 (F)T C θ θ θ θ θ ( ) - - = - K J K H J (17)
Below, this expression will be used for computing of the elastostatic deflections of the robotic manipulator.

Reduced mass matrix

To evaluate the dynamic behavior of the robot under the loading, in addition to the Cartesian stiffness matrix C K it is required to define the Cartesian mass matrix C M . This mass matrix has the same dimension as C

K

and can be obtained using some model reduction techniques. Comprehensive analysis and definition of this matrix have been proposed in [START_REF] Briot | Reduced elastodynamic modelling of parallel robots for the computation of their natural frequencies[END_REF]. Here, let us summarize the main results that will be used further.

To reduce the mass matrix dimension, model reduction techniques are applied for decreasing the size of the link mass matrices and also for the robot total mass matrix. Two main ways can be followed to reduce the size of the link mass matrices. The first one consists in discretizing the beam j into p j rigid links and springs and to express their displacements as a function of the beam extremity displacements. However, such numerical method must be repeated for each link and, thus, increases the size of the algorithm and decreases its efficiency. As a result, it is preferred to use the following procedure which allows analytical expressions to be obtained for the reduced link mass matrices.

Let us consider the link j, modeled as a beam (Fig. 2). At this beam is attached a local frame represented by the vectors x j , y j and z j . Before any deformation of the system, the beam j is linked to beams (j-1) and (j+1) at points O j and O j+1 , respectively (Fig. 2). After deformation of the robot, the beam extremity located at O j is displaced from

1 2 6 1 1 1 1 , , , δ δ δ -
- - - ⎡ ⎤ = ⎣ ⎦ K T j j j j t t t δt
and the one located at O j+1 is displaced from

1 2 6
, , ,

⎡ ⎤ = ⎣ ⎦ K T j j j j t t t δt δ δ δ
, where the three first components of each vector correspond to the translational displacements along local x j , y j and z j axes, respectively, and the three last components to the rotational displacements along the same axes.

The general formula for the kinetic energy of an elastic Bernoulli beam is equal to:

( ) 0 1 2 ; diag , , , , ρ = = ∫ & & j L T p y z j jj jj j j j jj j j T d x AAA I I I δ Q δ Q (18)
In this expression, & j δ represents the velocity of the beam cross-section located at position x from the local reference frame (Fig. 2), L j is the length of the beam j, ρ j the mass density at cross-section x, A j its area, p j I its torsional constant and y j I , z j I , the quadratic momentums along y j and z j , respectively.

COROUSSO

For the l-th natural mode l ω , the kinetic energy can be rewritten as:

( )

2 2 0 1 2 cos j L T jl l l l j j j j T t d x ω ω φ ρ = + ∫ δ Q δ (19)
δ j being the amplitude of the displacement of the beam cross-section located at position x from the local reference frame (Fig. 2).

In the Rayleigh-Ritz approximation, considering that the deformations due to the natural vibrations are similar to those obtained when an external load is applied at the robot end-effector only, each link of the structure will deform due to the stresses transmitted through the robot joints at points O j . As a result, the deformations ε j of the beam cross-section can be approximated by the deformations of a tip-loaded beam ( ) diag , , , , ,

j j j j j j j j f g g f h h = ε δ θ (20)
where ( )

j j j x L = = δθ ε
represents the deformation of the beam at its tip and (

, 0.5 3 , 2 0.5

j j j j j j j j f x x L g x x L x L h x x L x L = = - = - (21)
As a result, the global displacement δ j of the beam cross-section at x can be expressed as a sum of two terms:

( )

3 1 3 3 j j j × - ⎡ ⎤ = + ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ I D δ δ t ε 0 I , with ( ) 0 0 0 0 0 0 0 D x x × ⎡ ⎤ ⎢ ⎥ = - ⎢ ⎥ ⎢ ⎥ ⎣ ⎦ (22)
In this sum, the left terms corresponds to the displacement of the undeformed beam due to the displacement of the node located at O j .

Introducing (20) to ( 22) into (9) leads to the following equation:

( )

1 2 2 1 1 2 cos j T T r e d jl l l l j j j j T t ω ω ϕ - - ⎛ ⎞ ⎡ ⎤ ⎡ ⎤ = + ⎜ ⎟ ⎢ ⎥ ⎣ ⎦ ⎜ ⎟ ⎣ ⎦ ⎝ ⎠ δt δt δt M δt ( 23
)
where the expressions of each components of matrix red j M are given in [START_REF] Briot | Reduced elastodynamic modelling of parallel robots for the computation of their natural frequencies[END_REF].

Using these results, the total kinetic energy of the system for the l-th node is:

COROUSSO

the mass matrix can be recomputed into the Cartesian coordinates associated with the tool end-point using the Jacobian matrix θ J defined at expression (3) (which depend on the robot configuration q and computed with respect to virtual joint coordinates θ ) using following expression

θ C θ θ T = M J M J (25)
Thus, using expressions (25), it is possible to compute the reduced mass matrix

MACHINING PROCESS

Let us obtain the model of the cutting force which depends on the relative position of the tool with respect to the workpiece at each instant of machining. As follows from previous works [START_REF] Brissaud | The Fractional Model for the Determination of the Cutting Forces[END_REF], for the known chip thickness h, the cutting force F c can be expresses as

( ) ( ) 2 0 , 0 1 s s c s p h h r h h F h k a h h h + = ≥ + (26)
where p a is a depth of cut,

0 1 r k k ∞ =
< depends on the parameters k ∞ , k 0 that define the so called stiffness of the cutting process for large and small chip thickness h respectively (Fig. 3) and h s is a specific chip thickness, which depends on the current state of the tool cutting edge. The parameters k 0 , h s , r are evaluated experimentally for a given combination of tool/working material. To take into account the possible loss of contact between the tool and the workpiece, the above expression should be supplement by the case of

0 h < as ( ) 0, if 0 c F h h = < (27)
For the multi-edge tool the machining surface is formed by means of several edges simultaneously. The number of working edges varies during machining and depends on the width of cut. For this reason, the total force F c of such interaction is a superposition of forces F c,i generated by each tool edge i, which are currently in the contact with the workpiece. Besides, the contact force F c,i can be decomposed by its radial F r,i and tangential F t,i components (Fig. 4). In accordance with Merchant's model [START_REF] Merchant | Mechanics of metal cutting process. I-Orthogonal cutting and type 2 chip[END_REF], the t-component of cutting force F t,i can be computed with the equation ( 26). The r-component F r,i is related with F t,i by following expression [START_REF] Laporte | Phenomenological model for drilling operation[END_REF] , , r i r t i

F k F = ( 28
)
where the ratio factor k r depends on the given tool/workpiece characteristics.

COROUSSO

It should be mentioned that in robotic machining it is more suitable to operate with forces expressed in the robot tool frame {x,y,z}. Then, the corresponding components F x , F y (Fig. 4) of the cutting force F c can be expressed as follows

, , 1 1
, , 1 1 cos sin sin cos

z z z z n n x r i i t i i i i n n y r i i t i i i i F F F F F F ϕ ϕ ϕ ϕ = = = = = - + = + ∑ ∑ ∑ ∑ (29)
where n z is the number of currently working cutting edges, φ i is the angular position of the i-th cutting edge (the cutting force in z direction F z is negligible here). So, the vector of external loading of the robot due to the machining process can be composed in the frame {x,y,z} using the defined components F x , F y as F=[F x ,F y ,0,0,0,0] T .

Fig. 4. Forces of tool/workpiece interaction

It should be stressed that the cutting force components F r,i , F t,i mentioned in equation ( 26), (28) are computed for the given chip thickness h i , which should be also evaluated. Let us define model for h i using mechanical approach. Then the chip thickness h i removed by i-th tooth depends on the angular position φ i of this tooth and it can be evaluated using to the geometrical distance between the position of the given tooth i and the current machining profile (Fig. 4). It should be mentioned, that the main issue here is to follow the current relative position between the i-th tooth and the working material or to define whether the i-th tooth is involved in cutting for given instant of process. Because of the robot dynamic behavior and the regenerative mechanism of surface formation [START_REF] Tlusty | Basic non-linearity in machining chatter[END_REF] this problem cannot be solved directly using kinematic relations. In this case it is reasonable to introduce a special rectangular grid, which decomposes the workpiece area into segments and allows tracking the tool/workpiece interaction and the formation of the machining profile (Fig. 5). , N w is the number of nodes) of the grid can be marked as "1" or "0": "1" corresponds to nodes situated in the workpiece area with material (rose nodes in Fig. 6), "0" corresponds to nodes situated in workpiece area that was cut away (white nodes in Fig. 6).

Fig. 6. Evaluating the tool/workpiece intersection A i and computing the corresponding chip thickness h i

In order to define the number of currently cut nodes by the i-th tooth, the previous instant of machining process should be considered. Let us define A i as an amount of working material that is currently cut away by the i-th tooth (Fig. 6). So, if node j marked as "1" is located inside the marked sector (green nodes in Fig. 6), it changes to "0" and A i is increasing by

COROUSSO

thickness h i , removed at given instant of the process by the i-th tooth, can be estimated by

, i i i h A R α = Δ 1, z i N =
. The angle Δφ i determines the current angular position of the i-th tooth regarding to its position at the instant τ-Δτ and referred to the position of TCP at τ-Δτ.

Fig. 7. Algorithm for numerical simulation of robotic machining process dynamics

Described mechanism of chip formation and the machining force model (26) allow computing the dynamic behavior of the robotic machining process where models of robot inertia and stiffness are discussed in the section 3 of the paper. The detailed algorithm that is used in numerical analysis is presented in Fig. 7, where the analysis of the robot dynamics is performed in the tool frame with respect to the dynamic displacement of the tool δt dyn fixed on the robot end-effector around its position on the trajectory.

COMPLIANCE ERROR COMPENSATION TECHNIQUE

In industrial robotic controllers, the manipulator motions are usually generated using the inverse kinematic model that allows us to compute the input signals for actuators 0 ρ corresponding to the desired end-effector location 0 t , which is assigned assuming that the compliance errors are negligible. However, if the external loading F is essential, the kinematic control becomes non-applicable because of changes in the end-effector location. It can be computed from the non-linear compliance model as ( )

1 F 0 | f - = t F t (30)
where the subscripts 'F' and '0' refer to the loaded and unloaded modes respectively, and ' | ' separates arguments and parameters of the function ( )

f
. Some details concerning this function are given in our previous publication [START_REF] Pashkevich | Enhanced stiffness modeling of manipulators with passive joints[END_REF].

To compensate this undeterred end-effector displacement from 0 t to F t , the target point should be modified in such a way that, under the loading F , the end-platform is located in the desired point 0 t . This requirement can be expressed using the stiffness model in the following way ( ) To solve equation (31) for (F) 0 t , similar numerical technique can be applied. It yields the following iterative scheme ( )

(F) 0 0 | f = F t t (31)
(F) (F) 1 (F) 0 0 0 0 ( | ) • f α - = + - ′ t t t F t (32)
where the prime corresponds to the next iteration,

(0,1)
α ∈ is the scalar parameter ensuring the convergence.

More detailed presentation of the developed iterative routines is given in Fig. 8.

Hence, using the proposed computational techniques, it is possible to compensate a main part compliance errors by proper adjusting the reference trajectory that is used as an input for robotic controller. In this case, the control is based on the inverse kinetostatic model (instead of kinematic one) that takes into account both the manipulator geometry and elastic properties of its links and joints. Implementation of developed compliance error compensation technique presented in Fig. 9.

EXPERIMENTAL VERIFICATION

The developed compliance error compensation technique has been verified experimentally for robotic milling with the KUKA KR270 robot along a simple trajectory in aluminum workpiece. It is assumed that at the beginning of the technological process the robot is in the configuration q (see Table 1, Fig. 10). The parameters of the stiffness model for the considered robot have been identified in [START_REF] Dumas | Joint stiffness identification of six-revolute industrial serial robots[END_REF] and are presented in Table 1. Link masses required for the mass matrix of the robot are presented also in Table 1.

For the milling, the cutter with the external diameter D=20 mm and four teeth (N z =4) distributed uniformly over the tool is used. For the given combination of the tool and the workpiece material the following parameters correspond to the cutting force model defined in (26): k 0 =

COROUSSO

Taking into account that the workpiece has a straight borders let us assume that at the instant t=0 one of the teeth of the tool is in contact with the workpiece material as it is shown in the Fig. 11. It is also assumed that the machining process is performing with the constant feed rate v f =4 m/min (applied in x-direction of the robot tool frame) and the constant spindle rotation Ω=8000 rpm along the straight line of 80 mm. Experimental verification and numerical simulation of the described case of the milling process with KUKA KR-270 robot using the algorithm shown in Fig. 7 allows us to trace the evolution of machining force x,ycomponents for the whole process (Fig. 12). The corresponding dynamic displacement of the tool around its current position on the trajectory is shown in Fig. 13. In accordance with the obtained results the system robot/machining process realize complex vibratory motion. The high frequency component of this motion (about 700 Hz, Fig. 12) is related to the spindle rotation and the number of tool teeth N z . In certain cases such behavior can excites the dynamics of the robot (natural modes) but this study remains out the frame of the presented paper. On the contrary, the low frequency component of robot/tool motion (about 7 Hz, Fig. 13), especially in the y-direction (that is perpendicular to the applied feed) influences directly the quality of final product. Such motion is related to the robot compliance and it can be compensated using the error compensation technique described in the paper. Hence, let us form the modified trajectory based on the dynamic displacement of the robot end-effector in the y-direction (Fig. 14).

COROUSSO

It should be stressed that the time step between referenced points of this modified trajectory is limited with the characteristics of the controller used in the robot (in the presented case this step is chosen 0.05 sec). The corresponding feed rate v fy for the modified trajectory has been computed. So, this new data (feed f y and feed rate v fy ) with the data defined in the beginning of this section allow us to compensate the trajectory error during machining caused by the robot compliance. The resulted compensated trajectory in the y-direction (in time domain) is presented in Fig. 15.

Fig. 15. Evolution of the dynamic displacement obtained after involving the error compensation technique into the analysis of robotic milling process

It should be noted that the part of the trajectory while machining tool is engaging into the workpiece does not have effect on the quality of final product (surface). During this stage the contact area between the tool and the workpiece is increasing progressively. Hence, at each instant of processing the cutter corrects the machining profile and eliminates trajectory errors produced during all previous instants. On the contrary, during the stage of machining with the fully engaged tool the trajectory in x,y-directions define directly the final machining profile and this part of trajectory is analyzed here (Fig. 15). Comparison results presented in Fig. 13 and Fig. 15 are summarized in Table 2. So after applying error compensation technique the static deviation in y direction has been reduced from 0.058 mm to 0.00014 mm (99.8%). Maximum defilation in the machining profile has been reduced from 0.063 mm to 0.0047 mm (92.6%). Low frequency remained the same for both cases.

Table 2 Milling trajectory accuracy before and after compliance error compensation

CONCLUSION

In robotic-based machining, an interaction between the workpiece and technological tool causes essential deflections that significantly decrease the manufacturing accuracy. Relevant compliance errors highly depend on the manipulator configuration and essentially differ throughout the workspace. Their influence is especially important for heavy serial robots. To overcome this difficulty this paper presents a new technique for compensation of the compliance errors caused by technological process. In contrast to previous works, this technique is based on the non-linear stiffness model and the reduced elasto-dynamic model of the robotic based milling process.

The advantages and practical significance of the proposed approach are illustrated by milling with of KUKA KR270. It is shown that after error compensation technique significantly increase the accuracy of milling. In future the proposed technique will be integrated in a software toolbox.

ANR

ANR

ANR- 10

 10

Fig. 1 . 4 25

 14 Fig. 1. VJM model of industrial robot with end-point and auxiliary loading

COROUSSO

ANR- 10

 10

COROUSSO

ANR- 10

 10

COROUSSO

Fig. 2 . 4 25

 24 Fig. 2. Displacements and elastic deformations of a beam.

4 25

 4 Then, assuming that the first natural modes of vibrations, i.e. the modes that have the most energy, lead to deformations that are close to the static deformations of the robot under a load applied on the end-effector, ANR-10-SEGI-003-LI1.

CM 4 25

 4 for a given robot configuration q . The performances of this model reduction are shown in[START_REF] Briot | Reduced elastodynamic modelling of parallel robots for the computation of their natural frequencies[END_REF].ANR-10-SEGI-003-LI1.

COROUSSO

Fig. 3 . 4 25

 34 Fig. 3. Fractional cutting force model F c (h)

ANR- 10 Fig. 5 .

 105 Fig. 5. Meshing of the workpiece area

4 25

 4 Δ . Analyzing all potential nodes and computing A i , the chip ANR-10-SEGI-003-LI1.

ANR

COROUSSO

where (F) 0 t

 0 denotes the modified target location. Hence, the problem is reduced to the solution of the nonlinear equation (31) for (F) 0 t , while F and 0 t are assumed to be given. It is worth mentioning that this equation completely differs from the equation where the unknown variable is t . It means that here the compliance model does not allow us to compute the modified target point (F) 0 t straightforwardly, while the linear compensation technique directly operates with Cartesian compliance matrix[START_REF] Gong | Nongeometric error identification and compensation for robotic system by inverse calibration[END_REF].

Fig. 8 .Fig. 9 . 4 25

 894 Fig. 8. Procedure for compensation of compliance errors

COROUSSO

Fig. 10 .

 10 Fig. 10. Starting pose of the KUKA KR270 robot to perform the operation of milling

Fig. 11 .

 11 Fig. 11. Starting relative position of the tool with respect to the workpiece

Fig. 12 . 4 25Fig. 13 .

 12413 Fig. 12. Variation of machining force components F x (a) and F y (b) for whole milling process

Fig. 14 . 4 25

 144 Fig. 14. Modified trajectory f y and corresponding feed rate v fy in y-direction, computed based on the original dynamic displacement of the tool δt dyn

4 25

 4 deviation y MAX , [mm] 63.2e-3 4.70e-3 Hence, obtained results show that the developed compliance error compensation allows us significantly increase the accuracy of the robotic-based machining. ANR-10-SEGI-003-LI1.

COROUSSO

ANR

Table 1 .

 1 Initial data for robotic-based milling

	Joint coordinates, [deg]


Bilal Amghar
email: b.amghar@epmi.fr

A-Moumen Darcherif
email: m.darcherif@epmi.fr

Jean-Pierre Barbot
email: barbot@ensea.fr

Djamel Boukhetala
email: djamel.boukhetala@enp.edu.dz

Modeling and control of series resonant converter for high voltage applications

Keywords: Power converter, Petri nets, hybrid dynamic systems, resonant converter

Most industrial processes that use electricity as a source of voltage or current, require power circuit for the control of its physical quantities. Power converters are a multi-site privileged in the control of high power systems.In this article we show that the series resonant DC/DC converter, which is a hybrid system. The SRC achieves an output DC voltage equal to n times the input voltage of the converter, the major disadvantage of this type of converter is the the stresses on the power components and Capacitor Charge Time. Furthermore, a control strategy for minimizing the no-load conduction losses is proposed and the transient behavior in case of load steps including output shortcircuit is discussed based on digital simulations.

I. INTRODUCTION

High voltage dc power supplies have wide application range. There are several challenging and unique issues to be addressed while designing these power supplies. The Series Resonant Converter (SRC) is a preferred converter topology for high power, high voltage power supplies. The advantages commanded by the SRC are a good compromise over its disadvantages. The classical SRC is a frequency controlled converter. Control can also be achieved by keeping the operating frequency fixed but varying the duty cycle. This article studies the SRC for use in high voltage applications. Complete steady state analysis of the converter for operation above the resonant frequency is presented. Dynamic analysis of the converter is done through numerical simulation. Interest in reprate operation of pulsed power systems has led to an increase in power requirements of high voltage (HV) capacitor charging power supplies (CCPS). The requirements for a power supply charging the capacitor are severe The main requirements for these applications should be as follow: -Charging voltage regulated up to few kV.

-Charging time as low as few ms.

-The capability of withstanding short circuit conditions. High Voltage power of the required level is generally not available in the raw form. A suitable HV generating equipment has to be used to convert the low voltage available from the utility supply, on board generators or battery pack to the required voltage level. The HV dc was obtained by rectifying the stepped up voltage using either vacuum tubes or other electromechanical devices [START_REF] Garcia | Electronic driver without electrolynic capacitor for dimming high brightness LEDs[END_REF]. These were not only large but also difficult to control, inefficient and unreliable. Modern HV power generators are solid state having high frequency designs and provide superior performance and control. The first section of this paper is devoted to the modeling of the converter with a study of cell switching. A reconfiguration of the control switch is proposed in response to a variation between the current branches, close control by Petri nets is one of the solutions best suited to this type of system [START_REF] Hong | A High Efficiency Controller IC for LLC Resonant Converter in 0.35 µm BCD[END_REF]. Finally the results simulation are presented to validate the theoretical results.

II. THE RESONANT DC/DC CONVERTER MODEL

In Figure 1 we show a simplified nonlinear circuit representing the series resonant DC/DC power converter.

The SRC (Fig. 1) is a constant frequency, duty ratio controlled converter. The fixed switching frequency of the SRC is an advantage over the variable switching frequency operation of the classical SRC. The duty ratio is varied through phase modulation [START_REF] St Microelectronics | Reference design: wide range 200W L6599-based HB LLC resonant converter for LCD TV and flat panels[END_REF]. Conduction of switches on the leading leg of the inverter (T 11 , T 12 ) is phase shifted with respect to the conduction of switches on the lagging leg (T 21 , T 22 ) resulting in a quasi-square excitation A direct computation shows that the controlled nonlinear differential equations modelling the circuit are given by [START_REF] Garcia | Electronic driver without electrolynic capacitor for dimming high brightness LEDs[END_REF], [START_REF] St Microelectronics | Reference design: wide range 200W L6599-based HB LLC resonant converter for LCD TV and flat panels[END_REF], [START_REF] Amghar | TN )-Observability and control of parallel multicell chopper using Petri nets Author[END_REF].

               L r di dt = -R s i -v -mv o sign(i) + E(t) C r dv dt = i C o dvo dt = mabs(i) -vo R -I o (1)
where v and i are, respectively, the series capacitor voltage and the inductor current in the resonant series tank, while v o is the output capacitor voltage feeding both the load R and the sink current I o which, for simplicity, we assume to be of value zero m is transformation ratio. E(t) is a control input signal (usually taken as a periodic square wave of variable amplitude and frequency), even though for design purposes the amplitude is considered fixed leaving the frequency as the only control variable. Which is usually restricted to take values in the discrete set -V dc , V dc where V dc is a fixed given constant [START_REF] Garcia | Electronic driver without electrolynic capacitor for dimming high brightness LEDs[END_REF], [START_REF] St Microelectronics | Reference design: wide range 200W L6599-based HB LLC resonant converter for LCD TV and flat panels[END_REF].

The SRC is, therefore, equivalent to an RsLrCr circuit supplied by two alternating square voltages v AB (generated by the action of the switches Q1,Q2, Q3 and Q4) with an amplitude equal to the input voltage, and voltage V dc , as presented in Fig. 2. The resonant inductor (L r ) and resonant capacitor (C r ) are in series. The resonant tank and the load act as a voltage divider DC gain is always lower than 1 (maximum gain happens at the resonant frequency). The impedance of resonant tank can be changed by varying the frequency of driving voltage (E(t)). Let X ∈ R n be the continuous space of the inverter states and let Q ∈ {q 1 , q 2 , , q m } be a finite set of discrete states. In case of the inverter from Fig. 1 and the control strategy presented below, n = 3 and m = 5. The continuous state space specifies the possible values of the continuous states for all switching configurations q ∈ Q of the converter. For each q ∈ Q the continuous dynamics can be modeled by differential equations of the form [START_REF] Figge | Paralleling of LLC resonant converter using frequency controlled current balancing[END_REF], [START_REF] Amghar | Observability analysis for parallel muticell chopper[END_REF]:

ẋ = f (x, q, t) y = h(x, q, t) (2)
The states of the circuit from Fig. 1 are: x

(t) = [ i, v, v o ] T
. The discrete states of the system are defined by the states of the two-quadrant switches T 11 , T 12 , T 21 and T 22 . There are several possible strategies for the control of the inverter from Fig. 1.

The converter has four modes of operation defined by the four combinations of the values of the two power supplies Fig. 4.

Table I identifies the operation modes and the correspondent conducting switches and generated voltages. Systems including both continuous characteristics and discrete are called hybrid dynamic systems. In a very simplified one SDH has two sub sets, a continuous block, a block discrete: -The continuous block symbolizes the dynamic evolution of the state continue in our case the RLC resonant circuit and HV transformer.

Config 1 2 v AB +V dc -V dc u LC +V dc -vp -Rsi -V dc -vp -Rsi Switch Q 1 , Q 2 D 3 , D 4 Config 3 4 v AB -V dc +V dc u LC -V dc + vp -Rsi +V dc + vp -Rsi Switch Q 3 , Q 4 D 1 , D 2
-The block has the discrete system is discrete event receives internal events, external conditions, the converter is the state of the switches of switching cells.

In this work, we were interested in modeling method and control systems hybrid event-based dominant the use of Petri nets. The method is illustrated in Fig. 5. Figure 6 represents the Petri net control of the switches, the places P0, P1 and P2 respectively are modeling the initial state of converter and state of the switches of cells CELL1 and CELL2.

This algorithm is developed in order to control the output voltage of system. The transition from one place to another is dependent on the resonant current state of converter. The closure of the switch of the cell (Cell i ) depends on the validation of the transition t i0 . The Petri net contol is composed of tow arcs inhibitors, their role is to prevent the presence of more than one token in places P1 and P2, [START_REF] Amghar | TN )-Observability and control of parallel multicell chopper using Petri nets Author[END_REF], [START_REF] Niu Huapeng | Design of High Power Density DC Bus Converter Based on LLC Resonant Converter with Synchronous Rectifier[END_REF], [START_REF]LLC resonant half-bridge converter design guideline[END_REF].

Depending on the current i(t) and error e i = i ref -i the converter is allowed to set in 4 configurations.

The significance of all places and transition is shown in Table 1 and2 The switch of the first cell P 2

The switch of the second cell with : -δ v = Undulation of voltage output

-δ i = Undulation of current
The time response of system is one of the major problems of this type of converter . Pollution of the power system harmonics reactions is one of the other consequences of this problem, we will show in the simulation result the contribution of our approach on the input current and clean up the power grid.

B. Second case

In this part of work we studied two resonant converters in parallel. The regulation consists of two parts: -Two PI regulators for regulating of output voltage.

-A controller based on Petri nets for the regulation of currents i 1 and i 2 . the outputs of the networks Petri nets are the states of the switches cells.

The same as the first case the command consists of two parts, a continuous and a discrete part. The first is based on a control loop PI controller for controlling the both outputs voltages v s1 and v s2 (Fig. 8) . This loop has as inputs the error e s1 = v sref 1 -v s1 , e s2 = v sref 2 -v s2 and I ref 1 , I ref 2 as outputs current references. The second control loop is modeled by a Petri net whose mission is the currents regulation i 1 and i 2 respectively to value i ref 1 and i ref 2 calculated by the PI.

Figure 9 represents the Petri net control of the switches, the places P0, P1, P2, P3 and P4 respectively are modeling the initial state of converter and state of the switches of cells CELL1, CELL2, CELL3 and CELL4. This algorithm is developed in order to control the output voltage of system. The transition from one place to another is dependent on the resonant current state of converter. The closure of the switch of the cell (Cell i ) depends on the validation of the transition t i0 . The Petri net control is composed of arcs inhibitors, their role is to prevent the presence of more than one token in places Pi [START_REF] Amghar | TN )-Observability and control of parallel multicell chopper using Petri nets Author[END_REF], [START_REF] Tian | Modelling of asymmetrical pulse width modulation with frequency tracking control using phasor transformation for half-bridge series resonant induction cookers[END_REF], [START_REF] Choi | Analysis and Design of LLC Resonant Converter with Integrated Transformer[END_REF]. The significance of all places and transition is shown in Table 3 and4 The switch of cell1 P 2

The switch of cell2 P 3

The switch of cell3 P 4

The switch of cell 4 The Fig. 10 and Fig. 11 show respectively a simulation of the evolution of the output V s voltage and output current i s in both control. The figure Fig. 12 represents the results obtained by using th proposed control law on the converter. The result shows the hybrid control using Petri net control regulate the output voltage V s1 , V s2 and V s in a finite time. The figure Fig. 13 shows the evolution in Output current . Note that if the system has an imbalance current phase, the control by a Petri net consumes less current wavy and this is an important factor for the electrical systems. The dynamics of the output voltage is the image of the total current branch. The response time of the system is in the range of 3.10 -3 s. Simulation results have shown that the hybrid controller is efficient with respect to PI controller. The difference is in the instantaneous reaction to changes in physical parameters of the system.

V. CONCLUSION

This work consists of a modeling and control of Series Resonant Converter, an algorithm for controlling the switches of the converter using Petri nets has been proposed to solve the problem related to the time response and robustness of control. Petri nets are among the powerful tools for modeling and control of such systems which have discontinuities in their mathematical models. The algorithm is based on the statements of resonant current , the resonant current reference calculated by the PI and the configurations of switch cells. Finally simulation results show the convergence of output voltage V o to a neighborhood of the value of the reference output voltage V o in finite time. The evolution of the resonant current and primary voltage of transformer shows a remarkable performance of our approach.

Fig. 1 .

 1 Fig. 1. Principle of series resonance converter (SRC).

Fig. 2 .

 2 Fig. 2. Series Resonant DC-DC Converter equivalent circuit.

Fig. 3 .

 3 Fig. 3. DC Gain Caracteristic of Series Resonant DC-DC Converter.

Fig. 4 .

 4 Fig. 4. The five configurations status of the inverter.

Fig. 5 .

 5 Fig. 5. The schematic of close-loop control of DC-DC resonant converter.

Fig. 6 .

 6 Fig. 6. Hybrid control of DC-DC resonant converter.

Fig. 7 .

 7 Fig. 7. The Petri net control of system.

Fig. 8 .

 8 Fig. 8. DC-DC resonant converter.

Fig. 9 .

 9 Fig. 9. The Petri net control of system.

  with : -δ i1 = Undulation current through the first resonant circuit -δ i2 = Undulation current through the second resonant circuit IV. SIMULATION RESULTS The simulation results are obtained using the parameters of converter as follows:: R r = 2 mΩ, L s = 24 10 -6 H, L p = 150 10 -6 H, C s = 0.75 10 -6 F , V sref = 50 10 +3 V , R load = 10 Ω, C o = 100 10 -9 F , F o = 15 kHz and V dc = 560 V. After application of the close control of the converter, simulation results are shown in the following figures:

Fig. 10 .

 10 Fig. 10. Output voltage Vs variations in both control.

Fig. 11 .

 11 Fig. 11. Output current Is variations in both control.

Fig. 12 .

 12 Fig. 12. Output voltages V s1 , V s2 and Vs

Fig. 13 .

 13 Fig. 13. Output currents I s1 , I s2 and Is .

TABLE II THE

 II

TRANSITIONS Transitions Designations t 01 i ref > 0 and e i > δ i t 02 i ref < 0 and e i < -δ i t 10 i ref < 0 or (i ref > 0 ande i < δ i ) t 20 i ref > 0 or(i ref < 0 and e i > -δ i )

  > 0 and e i1 > δ i1t 02 i ref 1 < 0 and e i1 < -δ i1 t 03 i ref 2 > 0 and e i2 > δ i2 t 04 i ref 2 < 0 and e i2 < -δ i2 t 10 i ref 1 < 0 or (i ref 1 > 0 ande i1 < δ i1 ) t 20 i ref 1 > 0 or(i ref 1 < 0 and e i1 > -δ i1 ) t 30 i ref 2 < 0 or (i ref 2 > 0 ande i2 < δ i2 ) t 40 i ref 2 > 0 or(i ref 2 < 0 and e i2 > -δ i2 )

		TABLE IV
		THE TRANSITIONS
	Transitions	Designations
	t 01	i ref 1


Zhongkai Chen
email: zhongkai.chen@ensam.eu

Naf Ssa

Lakbakbi El

Gabriel Abba
email: abba@enim.fr

Optimal Walking of an Underactuated Planar Biped with Segmented Torso

Keywords: Biped Robot, Optimal Walk, Segmented Torso, Zero Dynamics, Energy Eff ciency

Recently, underactuated bipeds with pointed feet have been studied to achieve dynamic and energy eff cient robot walking patterns. However, these studies usually simplify a robot torso as one link, which is different from a human torsos containing 33 vertebrae. In this paper, therefore, we study the optimal walking of a 6-link planar biped with a segmented torso derived from its 5-link counterpart while ensuring that two models are equivalent when the additional torso joint is locked. For the walking, we suppose that each step is composed of a single support phase and an instantaneous double support phase, and two phases are connected by a plastic impact mapping. In addition, the controlled outputs named symmetry outputs capable of generating exponentially stable orbits using hybrid zero dynamics, are adopted to improve physical interpretation. The desired outputs are parameterized by Bézier functions, with 5-link robot having 16 parameters to optimize and 6-link robot having 19 parameters. According to our energy criterion, the segmented torso structure may reduce energy consumption up to 8% in bipedal walking, and the maximum energy saving is achieved at high walking speeds, while leaving the criteria at low walking speeds remain similar for both robots.

INTRODUCTION

One of the fundamental tasks in designing biped robots is to f nd robot structures and control strategies capable of providing energy-eff cient and stable walking motions. As indicated in Figure 1, many researchers have tried to study the energy saving benef ts of bionic robot designs by adopting certain human resembling parts such as curved feet [START_REF] Martin | The effects of curved foot design parameters on planar biped walking[END_REF], rolling knees [START_REF] Hobon | Gait optimization of a rolling knee biped at low walking speeds[END_REF], locked knees [START_REF] Sinnet | A human-inspired hybrid control approach to bipedal robotic walking[END_REF], swinging arms [START_REF] Kaddar | Arms swing effects on a walking planar biped[END_REF] or springs (Bauer et al., 2014). Nevertheless, the theoretical advantages of a human-like robot torso are still rarely studied among the vast literature on humanoid robots, despite the studies on realizing human-like biped torsos from a practical point of view [START_REF] Mizuuchi | Development of musculoskeletal humanoid kotaro[END_REF].

Unveiling the theoretical advantages of humanlike robot structures requires numerical optimization. By far, the most common method is through designing reference trajectories with a limited number of parameters, whether the trajectories are functions of time or functions of robot state, and then optimizing these parameters using a certain cost function. Many researchers have used optimization in the analysis

Curved feet Underactuation

Rolling knees Knee-locking

Swinging arms

Springs and dampers

Segmented torso

Compliant torso and design of biped walking motions with trajectories def ned as functions of time [START_REF] Cabodevila | Quasi optimal gait for a biped robot using genetic algorithm[END_REF][START_REF] Chevallereau | Optimal reference trajectories for walking and running of a biped robot[END_REF]Beletskii et al., 1982), and these early studies usually adopted the zero moment point (ZMP) principle to ensure the feasibility of various walking motions. However, a periodic and feasible walking motion might not be stable because biped robots are hybrid systems with impacts during each step, which might serve as a discrete positive feedback of error. Some researchers also attempted to eliminate the impact completely and thus creating walking motions that consume no mechanical energy, but these studies are only valid for passive biped robots walking in a highly controlled environment [START_REF] Rouhollah | Energy-conserving gaits for point-foot planar bipeds: A f ve-dof case study[END_REF]. For activated robots with motors, it's practically impossible to eliminate the energy consumption caused by Joule heating.

Recently, with more and more researchers striving to achieve more dynamic biped walking patterns, a planar biped robot with pointed feet and the name "RABBIT" was created, which witnessed the birth of the concept called "virtual constraints" [START_REF] Chevallereau | Rabit: A testbed for advanced control theory[END_REF]. By using a monotonic scalar quantity representing the progression of a step, [START_REF] Westervelt | Design of asymptotically stable walking for a 5-link planar biped walker via optimization[END_REF] proposed a method to def ne trajectory functions with a closed-form stability indicator. Armed with this method, Ames (2012) proposed human outputs to relate the walking motions of biped robots to the human walking data collected in a laboratory.

The work presented here is inspired by both the previously studied human-like parts and the optimization methods with provable stability. Our intention is to analyze the energy-saving potentials of the segmented torso. The optimal walking of a 5-link robot and a 6-link robot with underactuated pointed feet is studied to quantify the difference between these two robot conf gurations. The segmented torso of the 6-link robot is assumed to be straight at the beginning and the end of each step. Based on the stability analysis proposed by [START_REF] Westervelt | Design of asymptotically stable walking for a 5-link planar biped walker via optimization[END_REF], we def ne symmetry outputs to observe the mirror-like behavior of a virtual stance leg and a virtual swing leg. This interpretation of bipedal walking motions is related to a previous research on asymptotic stable walking with simplif ed posture control and swing leg advancement strategies [START_REF] Grizzle | Asymptotically stable walking for biped robots: Analysis via systems with impulse effects[END_REF]. By running optimization with symmetry outputs at different walking speeds, we have found obvious energy eff ciency improvement of the 6-link robot compared with its 5-link rival, especially at high walking speeds.

The paper is structured as follows. Section 2 delineates the biped models of both 5-link and 6-link robots. Section 3 introduces the symmetry outputs together with the feedback design using optimized HZD. Section 4 deals with the optimization problem using Bézier trajectories with provable stability. Section 5 presents the optimized walking motions of the two biped robots at different speeds. Finally, Section 6 concludes this paper by discussing the advantage of the additional torso joints and possible subsequent studies.

ROBOT MODELS

Robot Parameters

Both 5-link robot and 6-link robot have pointed feet which are not actuated during each step, and the rotor inertias are the same for all the actuators. As shown in Figure 2, the 5-link robot analyzed here comprises one torso, two femurs and two tibiae. We assume that all links are considered rigid and the motion is constrained to the sagittal plane. The 6-link robot, depicted in in Figure 3, is derived from the 5-link robot, and the only difference is that the additional torso joint divides the torso into an upper torso and a lower torso. In order to compare the performance of both robots, the 6-link robot should be equivalent to its 5-link counterpart once the additional joint is locked. This is realized by the following equations

M T I T L T p T M L f p f M M f I f M t I t p t M L t (x,y) q 5 q 1 q 2 q 4 q 3 x y u 1 u 2 u 3 u 4 x y
M T = M T 1 + M T 2 (1) M T p M T = M T 1 p M T 1 + M T 2 (L T 1 + p M T 2 ) (2)
I T =I T 1 + M T 1 (p M T -p M T 1 ) 2 + I T 2 + M T 2 (p M T -p M T 2 ) 2 (3)
where the notations used above is explained in Table 1. For the actuators of the femurs and tibiae, the rotor inertias equal 0.83 kgm 2 , and the motor rotor inertia of the additional torso joint is 0.2 kgm 2 .

Dynamic Models

The generalized coordinate vector describing both biped robots in the sagittal plane is def ned as q e = [x, y, q 1 , q 2 , ..., q n ] ′ where n = 5 or 6 equaling the number of links. It contains two coordinates (x and y) describing the position of pelvis and n angles q = [q 1 , q 2 , ..., q n ] ′ for the orientation of legs and torsos. We assume that the positive angles are anticlockwise.

( ) + ( , ) + ( ) = ( × ) ( × ) = [ , ] + = Δ( -) = Δ - Δ ( -) - Δ Δ ( -) × + = [ + , + ] -= [ -, -] -( ) = τ→ -(τ) + ( ) = τ→ + (τ) Σ = ( ) + ( ) ∈ + = Δ( -) -∈ = {( , )| ( ) = , ( ) < } ( ) ( ) ( ) = --( )( ( , ) + ( )) ( ) = × - ( -) = ( ) = ( ) -• θ( ) ( ) • θ( ) θ( ) θ( ) θ θ( ) = = --/ - θ ( ) = = [θ , , , ] θ = / + / + / + / + + π θ θ ( ) [θ , , , , ] s ¡¡ θ( ) ( ) ξ = θ( ) ξ = γ( , ) θ( ) γ( , ) ( ) ξ = θ - ζ = (ξ ) ζ ξ = κ (ξ ) κ (ξ ) κ (ξ ) κ (ξ ) κ (ξ ) = ∂θ ∂ ∂ ∂ γ - κ (ξ ) = - ∂ ∂ ( ) γ = ∂ /∂ = ( ) ξ - ξ + ξ + = δ ξ - δ δ = γ ( + )Δ ( -)σ ( -) σ ( -) = ∂ ∂ ( -) γ ( -) - (ξ ) = - ξ θ + κ (ξ) κ (ξ) ξ δ -δ (θ -) + < , (δ = ) = (ξ ) θ + ≤ ξ ≤ θ - < δ < ξ * = -| Ψ -, ζ * = - (θ -) -δ ξ * ζ * ξ (ξ , ζ * ) ξ (ξ , ζ * ) = -(δ ζ * - (ξ )) ( ) ( ) = - (ξ ) ξ = ∂ ∂ γ - ξ ( ) ( ) = -( ( , )) - ( , ) (ξ )

Figure 1 :

 1 Figure 1:Possible methods for reducing energy consumption.

Figure 2 :

 2 Figure 2: 5-link planar biped robot.