A Zavala-R Ío

On the Control of a One Degree-of-Freedom Juggling Robot

Keywords: nonsmooth mechanical systems, juggling, hybrid control, robustness, viability, non-prehensile manipulation

This paper is devoted to the feedback control of a one degree-of-freedom (dof) juggling robot, considered as a subclass of mechanical systems subject to a unilateral constraint. The proposed approach takes into account the whole dynamics of the system, and focuses on the design of a force input. It consists of a family of hybrid feedback control laws, that allow to stabilize the object around some desired (periodic or not) trajectory. The closed-loop behavior in presence of various disturbances is studied. Despite good robustness properties, the importance of good knowledge of the system parameters, like the restitution coefficient, is highlighted. Besides its theoretical interest concerning the control of a class of mechanical systems subject to unilateral constraints, this study has potential applications in non-prehensile manipulation, extending pushing robotic tasks to striking-and-pushing tasks.

Introduction

The feedback control of mechanical systems subject to unilateral constraints has received some attention. However, it remains largely unexplored. Recent research concerning analysis and control of such nonsmooth systems has focused on stabilization of manipulators [START_REF] Brogliato | On the control of finite dimensional mechanical systems with unilateral constraints[END_REF], [START_REF] Huang | Time optimal control for a robotic contour following problem[END_REF], wellposedness and system theoretical issues [START_REF] Ten Dam | The contact problem for linear continuous-time dynamical systems: A system theoretical approach[END_REF], [START_REF] Van Der Schaft | The complementary-slackness class of hybrid systems[END_REF], bipedal locomotion [START_REF] Hemami | Modelling and control of constrained dynamic systems with application to biped locomotion in the frontal plane[END_REF], [START_REF] Hurmuzlu | Dynamics of bipedal gait. Part 1: Objective functions and the contact event of a planar 5-link biped. Part 2: Stability analysis of a planar 5-link biped[END_REF], control of juggling and catching robots [START_REF] Aboaf | Task-level robot learning: Juggling a tennis ball more accurately[END_REF], [START_REF] Buehler | Robotic tasks with intermittent dynamics[END_REF], [START_REF] Buehler | A family of robot control strategies for intermittent dynamical environments[END_REF], [START_REF] Buehler | Planning and control of robotic juggling and catching tasks[END_REF], [START_REF] Rizzi | Preliminary experiments in robot juggling[END_REF], [START_REF] Rizzi | Progress in spatial robot juggling[END_REF], [START_REF] Rizzi | Further progress in robot juggling: the spatial two-juggle[END_REF], [START_REF] Rizzi | Further progress in robot juggling: solvable mirror laws[END_REF], [START_REF] Rizzi | An active visual estimator for dexterous manipulation[END_REF], [START_REF] Rizzi | Distributed real-time control of a spatial robot juggler[END_REF], [START_REF] Schaal | Open loop stable control strategies for robot juggling[END_REF], [START_REF] Schaal | Robot juggling: implementation of memory based learning[END_REF], [START_REF] Vincent | Controlling a ball to bounce at a fixed height[END_REF], [START_REF] Vincent | Controllable targets on or near a chaotic attractor[END_REF], [START_REF] Wang | Dynamics and planning of collisions in robotic manipulation[END_REF], [START_REF] Wang | Dynamic modeling and stability analysis of mechanical systems with time-varying topologies[END_REF], stabilization of polyhedral objects in some manipulation tasks [START_REF] Zumel | Balancing of a planar bouncing object[END_REF], non-prehensile manipulation [START_REF] Higuchi | Application of electro-magnetic impulsive force to precise positioning tools in robot systems[END_REF], [START_REF] Huang | Impulsive manipulation[END_REF], [START_REF] Lynch | Stable pushing: mechanics, controllability, and planning[END_REF], [START_REF] Peshkin | Planning robotic manipulation strategies for workpieces that slide[END_REF]. One common feature of these studies is that the open-loop models used are basically rigid body dynamics with a set of unilateral constraints on the generalized position. Such systems may be represented as follows:

M(q) q + C(q, q) q + g(q) = T u + ∇ q f (q, t)λ [START_REF] Aboaf | Task-level robot learning: Juggling a tennis ball more accurately[END_REF] f (q, t) ≥ 0 (2)

λ ≥ 0 , λ T f (q, t) = 0 (3)
where the classical dynamics of Lagrangian systems is in [START_REF] Aboaf | Task-level robot learning: Juggling a tennis ball more accurately[END_REF], the set of unilateral constraints is in [START_REF] Brach | Mechanical Impact Dynamics, Rigid Body Collisions[END_REF] with f (q, t) ∈ R , and the so-called complementarity conditions [START_REF] Van Der Schaft | The complementary-slackness class of hybrid systems[END_REF] are in [START_REF] Branicky | A unified framework for hybrid control: Background, model and theory[END_REF] where λ ∈ R is the Lagrange multiplier vector. We assume frictionless constraints, although this assumption is worthless concerning the one dof juggler (impacts are always central). q ∈ = {q ∈ R n : f (q, t) ≥ 0}, M(q) ∈ R n×n , C(q, q) q ∈ R n , g(q) ∈ R n and u ∈ R m are the generalized coordinate vector, the inertia matrix, the Coriolis and centrifugal terms, the gravity forces and torques, and the control input respectively, with T ∈ R n×m . In order to render the dynamical system complete, one must add to (1)-(3) a so-called restitution law that relates post and preimpact velocities. Such physical rules are necessary to render the domain invariant with respect to the system's dynamics. The most widely used restitution rule is known as Newton's conjecture [START_REF] Brach | Mechanical Impact Dynamics, Rigid Body Collisions[END_REF]. It is based on the knowledge of restitution coefficients e i , and is represented in its generalized form as follows:

q(t + k ) T ∇ q f i (q(t k ), t k ) = -e i q(t - k ) T ∇ q f i (q(t k ), t k ) -(1 + e i )∇ t f i (q(t k ), t k ) ( 4
)
where ∇ q f i ∈ R n , 1 ≤ i ≤ , and t k generically denotes the impact times (the superindices + andstand respectively for the instants just after and just before the collisions), with q(t - k ) T ∇ q f i (q(t k ), t k ) < 0. In case of a codimension 1 constraint ( = 1), e ∈ [0, 1] from energetical arguments. The system in (1)-( 4) is complete, in the sense that given preimpact velocities q(t - k ), one is able to calculate the postimpact velocities and continue the integration after the collision has occurred. Let us notice that the system in (1)-( 4) is a complex hybrid dynamical system [START_REF] Branicky | A unified framework for hybrid control: Background, model and theory[END_REF]. This class of dynamical systems can be divided further into subclasses. In particular the case when the free-motion dynamics are controllable has recently received attention [START_REF] Brogliato | On the control of finite dimensional mechanical systems with unilateral constraints[END_REF], [START_REF] Ten Dam | The contact problem for linear continuous-time dynamical systems: A system theoretical approach[END_REF]. However this class does not cover some impacting robotic systems like bipedal and juggling robots. Indeed let us write down the dynamical equations of a one dof juggler, i.e., a system composed of an object (a ball) subject to gravity, which rebounds on a controlled mass (a one dof robot) as shown in figure 1:

m 1 ÿ1 = -m 1 g (5) m 2 ÿ2 = -m 2 g + u (6) f (y 1 , y 2 ) = y 1 -y 2 ≥ 0 (7) ẏ1 (t + k ) -ẏ2 (t + k ) = -e ẏ1 (t - k ) -ẏ2 (t - k ) (8)
The dynamics in ( 5)-( 8) may be seen as a simple convex conical system [START_REF] Ten Dam | Representations of dynamical systems described by behavioral inequalities[END_REF] in the sense that if φ 1 and φ 2 are continuous-time trajectories of ( 5)- [START_REF] François | Contribution à la locomotion articulée dynamiquement stable[END_REF] then

βφ 1 + (1 -β)φ 2 , ∀β ∈ [0, 1]
, is also a continuous-time trajectory of the system. Most importantly, notice that the ball dynamics in [START_REF] Buehler | Robotic tasks with intermittent dynamics[END_REF] are not controlled since only gravity acts on the ball. Hence, the only way to influence the trajectory (y 1 , ẏ1 ) is through impacts, i.e., at times t k such that y 1 (t k )y 2 (t k ) = 0. This is a strong motivation for considering the feedback control of a one dof juggler as in ( 5)- [START_REF] François | Contribution à la locomotion articulée dynamiquement stable[END_REF] since it represents a simplified model of manipulation of objects through "controlled" collisions with a robotic device [START_REF] Zumel | Balancing of a planar bouncing object[END_REF]. Remark 1. More generally, juggling robots may also be considered as a particular case of a class of nonlinear systems:

ż1 = f 1 (z 1 , t, λ) (9) ż2 = f 2 (z 2 , u, λ) (10) h(q 1 , q 2 , t) ≥ 0 , λ T h(q 1 , q 2 , t) = 0 , λ ≥ 0 (11)
where

z i = q T i , qT i T ∈ R 2n i , i = 1, 2 (n 1
and n 2 may have the same or different values), u ∈ R m , λ ∈ R , and a restitution law must be added to complete the model. Our goal is to stabilize z 1 (t) around a desired trajectory. The controllability and stabilizability properties of such systems, which depend on the vector fields in [START_REF] Guckenheimer | Nonlinear Oscillations, Dynamical Systems and Bifurcations of Vector Fields[END_REF], [START_REF] Hemami | Modelling and control of constrained dynamic systems with application to biped locomotion in the frontal plane[END_REF] and on h in [START_REF] Higuchi | Application of electro-magnetic impulsive force to precise positioning tools in robot systems[END_REF], have not yet been fully understood. However, as we shall see later, some basic properties are quite useful. Among those: i) the explicit knowledge of the trajectories of (9) between impacts, ii) the controllability of [START_REF] Hemami | Modelling and control of constrained dynamic systems with application to biped locomotion in the frontal plane[END_REF], and iii) the controllability of (9) through collisions.

Property iii) will be made clear in remark 5. Interestingly enough, let us note that those assumptions are quite similar to the assumptions made in [START_REF] Kolmanovsky | Hybrid feedback laws for a class of cascade nonlinear control systems[END_REF], although the class of nonlinear systems they considered and the control objectives are quite different (they deal with asymptotic stabilization to the origin of a class of smooth cascaded systems modeling nonholonomic mechanisms). The interest for considering systems as in ( 9)- [START_REF] Higuchi | Application of electro-magnetic impulsive force to precise positioning tools in robot systems[END_REF] is that their study finds potential applications in all types of juggling robots, catching tasks, nonprehensile manipulation (pushing-and-striking tasks), and stabilization of manipulators on passive dynamical environments. It is clear that in practice, disturbances and parameter uncertainties will be an obstacle to the use of assumption i). Hence, robustness studies will be very important for such systems. Finally, it is worth noting the similarities between systems as in ( 9)- [START_REF] Higuchi | Application of electro-magnetic impulsive force to precise positioning tools in robot systems[END_REF] and walking or hopping machines. Indeed in both cases there is a part of the system that is not controllable during flight times (i.e., the center of mass motion for walking robots). We therefore conjecture that the dynamics in ( 9)-( 11) may represent a very large class of systems subject to unilateral constraints. The role of the vector fields f 1 , f 2 and of h in the stabilizability properties of such systems remains however an open field, as well as that of the restitution rule 1 (recall that a restitution law must be added to complete the model). The study of simple juggling robots is expected to be a nice manner to investigate their basic properties.

Review of Previous Work

The structure of mechanical systems as in ( 5)-( 8) naturally leads to the study of their dynamics and the control design at impact times only through their so-called impact Poincaré map P (the Poincaré section is chosen as = {(q, q, t) ∈ × R n × R + : f (q, t) = 0 , qT ∇ q f (q, t) > 0}). One way to obtain the closed-loop impact Poincaré map is to integrate [START_REF] Buehler | Robotic tasks with intermittent dynamics[END_REF] and ( 6) between impacts to get

y 1 (t k+1 ) = y 1 (t k ) + ẏ1 (t + k ) k - g 2 2 k (12) ẏ1 (t - k+1 ) = ẏ1 (t + k ) -g k ( 13
)
y 2 (t k+1 ) = y 2 (t k ) + ẏ2 (t + k ) k + t k+1 t k t t k u(s) ds dt (14) ẏ2 (t - k+1 ) = ẏ2 (t + k ) + t k+1 t k u(t) dt ( 15
)
t k+1 = t k + k ( 16
)
and then to seek a controller u such that y

1 (t k+1 ) = y 2 (t k+1 ) = y * , ẏ1 (t + k+1 ) = ẏ * 1 and ẏ2 (t + k+1 ) = ẏ * 2 (
using also the restitution law). However, it is not easy in general to find such a u. This philosophy has been employed in [START_REF] Buehler | Robotic tasks with intermittent dynamics[END_REF], [START_REF] Buehler | A family of robot control strategies for intermittent dynamical environments[END_REF], [START_REF] Buehler | Planning and control of robotic juggling and catching tasks[END_REF], [START_REF] Vincent | Controlling a ball to bounce at a fixed height[END_REF], [START_REF] Vincent | Controllable targets on or near a chaotic attractor[END_REF] to control juggling robots. Basically, these authors used the assumption that m 1 m 2 = 0, i.e., the robot's velocity ẏ2 is not affected by the shocks. Then the open-loop system is reduced to a 2-dimensional system (the ball dynamics in (5)), together with a time-varying unilateral constraint y 1 (t) ≥ y 2 (t). It is then possible to obtain an explicit form of P . Buehler, et al. [START_REF] Buehler | Robotic tasks with intermittent dynamics[END_REF], [START_REF] Buehler | A family of robot control strategies for intermittent dynamical environments[END_REF], [START_REF] Buehler | Planning and control of robotic juggling and catching tasks[END_REF] formulate the problem as designing a feedback control for P , where the inputs are first considered as being a sequence of flight-times k = t k+1t k and robot velocities ẏ2 (t - k ). Then, they propose a heuristic strategy which they call the mirror law (because it aims at making y 1 and y 2 more or less symmetrical with respect to the shock position y 1 = y 2 = 0). This work was extended in [START_REF] Rizzi | Preliminary experiments in robot juggling[END_REF], [START_REF] Rizzi | Progress in spatial robot juggling[END_REF], [START_REF] Rizzi | Further progress in robot juggling: the spatial two-juggle[END_REF], [START_REF] Rizzi | Further progress in robot juggling: solvable mirror laws[END_REF], [START_REF] Rizzi | An active visual estimator for dexterous manipulation[END_REF], [START_REF] Rizzi | Distributed real-time control of a spatial robot juggler[END_REF] to more complex jugglers. Experimental results have also been presented by these authors. The control strategies are based on continuous feedback of the object state, measured either via cameras [START_REF] Rizzi | An active visual estimator for dexterous manipulation[END_REF], [START_REF] Rizzi | Distributed real-time control of a spatial robot juggler[END_REF] or via a digitizing table [START_REF] Buehler | A family of robot control strategies for intermittent dynamical environments[END_REF], [START_REF] Buehler | Planning and control of robotic juggling and catching tasks[END_REF]. Vincent [START_REF] Vincent | Controlling a ball to bounce at a fixed height[END_REF], [START_REF] Vincent | Controllable targets on or near a chaotic attractor[END_REF] directly uses the 2-dimensional map P derived in [START_REF] Holmes | The dynamics of repeated impacts with a sinusoidally vibrating table[END_REF] and nicely combines results on the bouncing ball dynamics [START_REF] Guckenheimer | Nonlinear Oscillations, Dynamical Systems and Bifurcations of Vector Fields[END_REF] to obtain a hybrid controller. First y 2 (t) is set to create chaotic motion of the ball. Then the robot's motion is switched to stabilize (y 1 , ẏ1 ) around a periodic trajectory when its basin of attraction has been attained. The control strategy uses the measurements of flight times obtained via a microphone attached to the robot. However, the main drawbacks of this strategy are that the transient period length is not guaranteed, moreover its extension to more complex systems like two dof jugglers is not trivial. Finally, learning techniques have been applied to various juggling systems in [START_REF] Aboaf | Task-level robot learning: Juggling a tennis ball more accurately[END_REF], [START_REF] Schaal | Open loop stable control strategies for robot juggling[END_REF], [START_REF] Schaal | Robot juggling: implementation of memory based learning[END_REF].

In this work we consider the full dynamics ( 5)- [START_REF] François | Contribution à la locomotion articulée dynamiquement stable[END_REF], i.e. we allow for discontinuities in ẏ2 . Indeed although the assumption m 1 m 2 = 0 may be satisfied in some practical cases (like juggling with a ping-pong ball), this is not a universal property, and it may not be satisfied in the general case of pushing-and-striking manipulation tasks. Moreover, its relaxation makes it difficult to employ the above cited strategies, since P has no longer dimension 2, but 4. Most importantly, let us notice that the design of the force input u in (6) has never been presented in the previous works on the control of juggling robots. We are therefore interested in designing directly the u that is to be implemented on the robot. This note is organized as follows: In section 3 we present the controller and the closed-loop analysis for the one dof juggler. Section 4 is devoted to examine the closed-loop behavior when various types of disturbances (bad knowledge of the restitution coefficient e, bad velocity measurements) are present. Conclusions are given in section 5. Several numerical simulations illustrate the theoretical investigations.

Main Results

Before developing the proposed controller u in (6), let us define the objectives of the one dof juggling task: given the dynamics in ( 5)-( 8), make the ball attain a periodic motion from any initial condition x(0), where x T = (y 1 , y 2 , ẏ1 , ẏ2 ), and for any m = m 1 m 2 and e ∈ [0, 1]. In other words, since the trajectory (y 1 (t), ẏ1 (t)) can be controlled only through impacts with the robot at times t k , the goal is to design u in [START_REF] Buehler | A family of robot control strategies for intermittent dynamical environments[END_REF] such that the collisions occur at (k) and ẏd (k) to be tracked. Let us note that this definition is slightly different from some other works: for instance in [START_REF] Vincent | Controlling a ball to bounce at a fixed height[END_REF], [START_REF] Vincent | Controllable targets on or near a chaotic attractor[END_REF], only the apex is regulated while the impact position is not considered. In the next subsection we a priori introduce an input u. We shall see further (see remark 5) that it can be deduced from a more general control design methodology.

y 1 (t k ) = y 2 (t k ) = y d , ẏ1 (t + k ) = ẏd > 0,

The Force Input Controller u

CLAIM 1 Consider the dynamical system in ( 5)- [START_REF] François | Contribution à la locomotion articulée dynamiquement stable[END_REF]. Suppose that x(0) and u are such that it exists an impact or a contact time t 0 . Let us define the following control input:

u = m 2 g + m 2 v (17) v = A k (t -t k ) + B k (18) k ≥ 0, with A k = 6 d 2 k ẏ * 2 (k + 1) + ẏ2 (k) - 12 d 3 k y * (k + 1) -y(k) ( 19
)
B k = - 2 d k ẏ * 2 (k + 1) + 2 ẏ2 (k) + 6 d 2 k y * (k + 1) -y(k) ( 20
)
d k = ẏ1 (k) + ẏ2 1 (k) -2g (y * (k + 1) -y(k)) g (21) ẏ * 2 (k + 1) = 1 + m 1 + e ẏ * 1 (k + 1) + m -e 1 + e ẏ2 1 (k) -2g (y * (k + 1) -y(k)) ( 22
)
where

ẏ1 (k) = ẏ1 (t + k ), ẏ2 (k) = ẏ2 (t + k ), y(k) = y 1 (t k ) = y 2 (t k ), and y * (k +1) and ẏ * 1 (k +1) are chosen such that y * (k + 1) < h k = ⎧ ⎨ ⎩ y(k) + ẏ2 1 (k) 2g if ẏ1 (k) > 0 y(k) if ẏ1 (k) ≤ 0 ( 23
)
and

ẏ1 (t k + d k ) -ẏ * 2 (k + 1) < 0 (24)
Then, for all k ≥ 0:

1. y 1 (t) -y 2 (t) > 0, ∀ t ∈ (t k , t k + d k ) 2. t k+1 = t k + d k 3. y(k + 1) = y * (k + 1) 4. ẏ1 (k + 1) = ẏ * 1 (k + 1)
Proof: The proof can be found in [START_REF] Zavala-Río | Hybrid feedback strategies for the control of juggling robots[END_REF], [START_REF] Zavala-Río | Commande de robots jongleurs[END_REF] and is given in appendix A.

Remark 2. Controllers as in (52), ( 53) are basically open-loop. But one notices from ( 17)-( 24) that the controller is computed with x(t + k ). Consequently u is a state feedback for the system considered as a discrete-time operator at the shock times t k . More generally, the presented methodology for stabilization of the system in ( 9)-( 11) at some desired trajectory at impact times may be based on a continuous-time feedback of z 2 and a discrete-time feedback of z 1 . This is another common point with the assumptions in [START_REF] Kolmanovsky | Hybrid feedback laws for a class of cascade nonlinear control systems[END_REF] (see remark 1). 52), (53) that one can modify the matrix F to obtain different controllers in [START_REF] Kreisselmeier | Adaptive control of continuous-time overmodeled plants[END_REF] by merely modifying u in (17) (and provided that the obtained system in ( 6) is controllable with v as the new input). For instance let us choose

Remark 3. It clearly appears from (

u = m 2 g - m 2 ω 2 k y 2 + m 2 v with ω k = 2π d k .
Then, as shown in [START_REF] Zavala-Río | Commande de robots jongleurs[END_REF], [START_REF] Zavala-Río | On the feedback control of a simple juggling robot[END_REF], the equation f (y

1 (t), y 2 (t)) = y 1 (t) -y 2 (t) = f ( ) can be written as: f ( ) = y(k) + ẏ1 (k) - g 2 2 -y(k) - A k 2ω k cos(ω k ) - ẏ2 2 (k) ω k + A k 2ω 2 k + B k 2ω k sin(ω k ) ( 25
)
where

= t -t k > 0.
In order to prove that f ( ) > 0, ∀ ∈ (0, d k ), one must first demonstrate that f ( ) has no roots on (0, d k ) which, from [START_REF] Rizzi | Further progress in robot juggling: solvable mirror laws[END_REF], is not easy. Since point 1 in claim 1 is crucial for the overall scheme to work, we choose to disregard those controllers for which we cannot prove that it is satisfied. Such conditions of non-existence of "accidental" shocks are general in the analysis of vibro-impact systems. In particular for results concerning existence of periodic trajectories see [START_REF] Hindmarsh | On the motions of the offset impact oscillator[END_REF]. They may be called viability conditions. They are fundamental conditions in the study of vibro-impact systems which have sometimes been forgotten [START_REF] Masri | On the stability of the impact damper[END_REF], [START_REF] Shaw | The transition to chaos in a simple mechanical system[END_REF].

Remark 4. Let us define the vector x (k) T = (y(k), ẏ1 (k), ẏ2 (k), t k ). Then, the dynamics of the impact Poincaré map

P : → x (k) → x (k + 1) (26) with = {(x, t) ∈ R 4 × R + : y 1 -y 2 = 0, ẏ1 -ẏ2 > 0}
, are guaranteed in closed-loop via the chosen input force u in ( 17)-( 24). As we pointed out in the introduction, a possible path for the control design is to use equations ( 12)-( 16). For instance, introducing ( 17), ( 18) into ( 14), ( 15) one obtains:

y 2 (k + 1) = y 2 (k) + ẏ2 (k) k + B k 2 2 k + A k 6 3 k (27) ẏ2 (t - k+1 ) = ẏ2 (k) + B k k + A k 2 2 k (28)
Considering ( 12), ( 13), ( 27), ( 28), [START_REF] Hurmuzlu | Dynamics of bipedal gait. Part 1: Objective functions and the contact event of a planar 5-link biped. Part 2: Stability analysis of a planar 5-link biped[END_REF], and (59), (60), we get a 4-dimensional mapping

x (k + 1) = P (x (k), A k , B k , k ), given by x (k + 1) = ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ y(k) + ẏ1 (k) k -g 2 2 k m-e 1+m ( ẏ1 (k) -g k ) + 1+e 1+m ẏ2 (k) + k B k + 2 k 2 A k m(1+e) 1+m ( ẏ1 (k) -g k ) + 1-em 1+m ẏ2 (k) + k B k + 2 k 2 A k t k + k ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ (29)
with A k , B k , and k as control parameters, and an algebraic equation

k A k 6 2 k + B k + g 2 k + ( ẏ2 (k) -ẏ1 (k)) = 0 ( 30
)
It is noteworthy that although k may be considered an input signal for the application in [START_REF] Schaal | Robot juggling: implementation of memory based learning[END_REF], in reality k is a consequence of the dynamics and of the real input (A k , B k ). In fact, from [START_REF] Shaw | The transition to chaos in a simple mechanical system[END_REF], k can be expressed in terms of A k and B k and introduced in (29) to get an explicit system's Poincaré mapping

x (k + 1) = P (x (k), A k , B k ) ( 31
)
Notice that an expression for k cannot always be obtained explicitly. Besides, the resultant mapping ( 31) is usually non-linear in its inputs. This is what renders the design of the control input (A k , B k ) for tracking purposes directly from (31) a difficult task. One could analyze the stabilization of ( 31) by linearizing it around a desired trajectory. Then, a local linear feedback can be designed provided the linearization is controllable. Buehler, et al. [START_REF] Buehler | Robotic tasks with intermittent dynamics[END_REF], [START_REF] Buehler | A family of robot control strategies for intermittent dynamical environments[END_REF], [START_REF] Buehler | Planning and control of robotic juggling and catching tasks[END_REF] did similarly for their 2-dimensional mapping, but found that such a control strategy does not provide a satisfactory behavior in practice. We have overcome such difficulty using the Dead-beat controllers which allow us to directly assign a value to k .

Remark 5 (3-step recursive control design).

From the above developments it follows that considering the impact Poincaré map P in [START_REF] Sontag | Mathematical Control Theory, Deterministic Finite Dimensional Systems[END_REF] for control design purposes is not a good solution. As we announced in the introduction of section 3, we now describe a recursive control design method, which will enable us to recover the Dead-beat strategy. Let us define

v(t) = A k + B k , = t -t k , with A k and B k constants on (t k , t k+1
). From ( 12)-( 16) and (59) we get

⎧ ⎨ ⎩ y 1 (t k+1 ) = y 1 (t k ) + ẏ1 (t + k ) k -g 2 2 k ẏ1 (k + 1) = m-e 1+m ẏ1 (t - k+1 ) + 1+e 1+m ẏ2 (t - k+1 ) (32) and ⎧ ⎨ ⎩ y 2 (t k+1 ) = y 2 (t k ) + ẏ2 (t + k ) k + A k 3 k 6 + B k 2 k 2 ẏ2 (t - k+1 ) = ẏ2 (t + k ) + A k 2 k 2 + B k k (33)
Step 1: We choose k and ẏ2 (t - k+1 ) as the inputs of the system in [START_REF] Ten Dam | Representations of dynamical systems described by behavioral inequalities[END_REF], such that y

1 (t k+1 ) = y * (k + 1), ẏ1 (t + k+1 ) = ẏ * 1 (k + 1). This gives (32) ⇔ F 1 (x 1 (k), x * 1 (k + 1), k , ẏ2 (t - k+1 )) = 0 ⇒ k = d k ẏ2 (t - k+1 ) = ẏ * 2 (k + 1) (34)
where x T 1 = (y 1 , ẏ1 ); d k and ẏ * 2 (k + 1) are given in claim 1 (eqs. ( 21) and ( 22) respectively), and for some function F 1 .

Step 2: Introducing the values in (34) into (33), we get

F 2 (x(k), d k , ẏ * 2 (k+1), A k , B k ) = 0,
for some function F 2 , from which we deduce A k and B k , which are equal to those in claim 1 (eqs. [START_REF] Lynch | Stable pushing: mechanics, controllability, and planning[END_REF] and (20) respectively).

Step 3: Check viability, i.e., the sign of the function f

( ) = y 1 (t) -y 2 (t), on (t k , t k+1 ).
It is noteworthy that the success of the first step relies on the "controllability" (or invertibility) properties of the first subsystem with k and ẏ2 (t - k+1 ) as inputs. The second step mainly hinges on the satisfaction of the viability conditions. In other words one can choose another controller structure v(t). For instance

v(t) = n i=0 A k,i i (35)
The viability conditions can be studied from the sign of f ( ) as in (57), which is a polynomial in . Clearly for n ≥ 2 in [START_REF] Vincent | Controlling a ball to bounce at a fixed height[END_REF], the sign of f ( ) is not easy to obtain. This is why the solution n = 1 has been chosen. Let us finally notice that it is not proven yet that to each (polynomial) controller as in [START_REF] Vincent | Controlling a ball to bounce at a fixed height[END_REF] there corresponds a Dead-beat controller (possibly designed after a dynamic pre-feedback in order to get the right matrix F in (51)-( 53)). Extending the above procedure to systems as in ( 9)-( 11) is a topic of future research.

Definition of the Desired Trajectories

The next step is to define the signals y * (k + 1) and ẏ * 1 (k + 1) which can be regarded as the desired trajectory of the ball. We still consider the system at a generic impact time t k . COROLLARY 1 Let y * (k + 1) and ẏ * 1 (k + 1) be given as:

y * (k + 1) = y d if h k > y d y(k) + r if h k ≤ y d (36) ẏ * 1 (k + 1) = ẏd if h k > y d ẏ2 d + 2g (y d -y(k) -r ) if h k ≤ y d ( 37
)
where ẏd > 0, and r is a real value chosen such that condition [START_REF] Rizzi | Progress in spatial robot juggling[END_REF] 29), it is easy to verify that whenever the control input ( 17)-( 24) is applied, the system's impact Poincaré mapping (26) is given by

x (k + 1) = ⎡ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎢ ⎣ y * (k + 1) ẏ * 1 (k + 1) 1-em 1+e ẏ * 1 (k + 1) -e(1+m) 1+e ẏ2 1 (k) -2g (y * (k + 1) -y(k)) t k + ẏ1 (k)+ √ ẏ2 1 (k)-2g(y * (k+1)-y(k)) g ⎤ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎥ ⎦ (38)
From the iterative substitution of ( 36) and ( 37) in [START_REF] Wang | Dynamic modeling and stability analysis of mechanical systems with time-varying topologies[END_REF], the proof follows.

Notice that this consitutes mainly a theoretical result, since in practice more than two shocks will often be necessary before the convergence to the desired value, mainly due to input saturations. In the next subsection a strategy is proposed to cope with this problem. Example 1. Let us consider the application of the control strategy ( 17)-( 22), ( 36) and ( 37

Control Magnitude Reduction via Desired Motion Modification

The impact trajectories defined in [START_REF] Vincent | Controllable targets on or near a chaotic attractor[END_REF] and [START_REF] Wang | Dynamics and planning of collisions in robotic manipulation[END_REF] are designed in such a way that the ball reaches the desired apex s d = y d + ẏ2 d 2g just after the first control collision. This may require continuous-time input forces of big magnitudes during the first flight interval (t 0 , t 0 + d 0 ) (see figure 2b). Let us consider the realistic case where the final actuator of the juggling robot provides forces bounded in magnitude such that |u| ≤ u max . From [START_REF] Kreisselmeier | Adaptive control of continuous-time overmodeled plants[END_REF], this implies v 1 ≤ v ≤ v 2 . One sees from ( 36) and ( 37) that the signals y * (k + 1) and ẏ * 1 (k + 1) can be chosen time-varying in order to cope with the input magnitude. Intuitively, one easily conceives that the desired apex s d can be attained by imparting several small strikes to the ball rather than one big hit (this is for instance what one does with a racket to make a tennis ball rebound higher and higher). In this subsection we present such a strategy. COROLLARY 2 Let us consider that v 1 < -3+e-2m+2em 1+e g and v 2 > 3-2e+4m-em 1+e g, and assume that:

1. The impact states just after the first shock at t 0 are such that ẏ1 (0) > 0 and

- v 2 g ẏ1 (0) ≤ ẏ2 (0) ≤ min{M 1 , M 2 , 1} ẏ1 (0)
where

M 1 = -1 2 v 1 g + m-e 1+e + e M 3 1+m 1+e , M 2 = v 2 g -2(m-e) 1+e -2e M 3 1+m 1+e , M 3 = min{M 4 , M 5 }, M 4 = 1+e 1+m v 2 2g + 3-e+2m-2em 2(1+m)
, and

M 5 = -1+e 1+m v 1 g + 3-2e+4m-em 1+m .
2. The desired trajectories are designed, at each impact, as: y * (k + 1) = y(k) and

ẏ * 1 (k + 1) = ⎧ ⎪ ⎪ ⎪ ⎨ ⎪ ⎪ ⎪ ⎩ N 1 (k) if ẏd < N 1 (k) ẏd if N 1 (k) ≤ ẏd ≤ N 2 (k) N 2 (k) if ẏd > N 2 (k) ( 39
)
∀k ≥ 0, for any positive value of ẏd , where N

1 (k) = max N 3 (k), N 4 (k), e M 3 ẏ1 (k) , N 2 (k) = min {N 5 (k), N 6 (k)}, N 3 (k) = -1+e 1+m v 2 g + m-e 1+m ẏ1 (k)-2(1+e) 1+m ẏ2 (k), N 4 (k) = 1+e 1+m v 1 2g -m-e 1+m ẏ1 (k) -1+e 2(1+m) ẏ2 (k), N 5 (k) = -1+e 1+m v 1 g + m-e 1+m ẏ1 (k) -2(1+e) 1+m ẏ2 (k), and N 6 (k) = 1+e 1+m v 2 2g -m-e 1+m ẏ1 (k) -1+e 2(1+m) ẏ2 (k).
Then, the impact velocity of the ball converges to ẏd without saturating the input, i.e.,

v 1 < v(t) < v 2 , ∀t ≥ t 0 .
Proof: The proof can be found in [START_REF] Zavala-Río | Commande de robots jongleurs[END_REF], [START_REF] Zavala-Río | Hybrid feedback strategies for the control of juggling robots[END_REF] and is not given in this note for the sake of briefness.

Let us explain the underlying philosophy of the strategy presented in corollary 2. Assume that v 1 and v 2 are as pointed out in the corollary, and that an impact such that ẏ1 (k) > 0 takes place. Let us define the desired impact position as y * (k + 1) = y(k). Then, in order to have ẏ1 (k + 1) > 0 without saturating the control input, ẏ2 (k) must be within

1 (k) = ẏ2 (k) ∈ R : - v 2 g ẏ1 (k) ≤ ẏ2 (k) ≤ min(M 1 , M 2 , 1) ẏ1 (k)
(M 1 and M 2 are defined in point 1 of corollary 2), and the set of nonnegative values within which ẏ * 1 (k + 1) must be selected is

2 (k) = ẏ * 1 (k + 1) ∈ R : N 1 (k) ≤ ẏ * 1 (k + 1) ≤ N 2 (k)
(N 1 (k) and N 2 (k) are defined in point 2 of corollary 2). Besides, in order for 2 (k) to include ẏ1 (k), so that ẏ * 1 (k + 1) can be chosen greater or smaller than, or even equal to, ẏ1 (k) (having all these possibilities in the same time), ẏ2 (k) must be within

3 (k) = { ẏ2 (k) ∈ R : max (P 1 , P 2 ) ẏ1 (k) ≤ ẏ2 (k) ≤ min(P 3 , P 4 , 1) ẏ1 (k)} ⊂ 1 (k)
where

P 1 = v 1 g -2(1-e+2m) 1+e , P 2 = -1 2 v 2 g + 1-e+2m 1+e , P 3 = -1 2 v 1 g + 1-e+2m 1+e
, and

P 4 = v 2 g -2(1-e+2m) 1+e
. Moreover, the selection of ẏ * 1 (k+1) ∈ 2 (k) ensures that ẏ2 (k+1) will be within 3 (k + 1), for any ẏ2 (k) ∈ 1 (k). Now, notice that the conditions established in point 1 and the design criterion proposed in point 2 of corollary 2 are such that ẏ2 (0) ∈ 1 (0) and ẏ * 1 (k + 1) ∈ 2 (k), ∀k ≥ 0. Moreover, (39) always assigns to ẏ * 1 (k + 1) the closest value to ẏd > 0 inside 2 (k). Then, the convergence of ẏ1 (k) to ẏd and the non-saturation of the input are guaranteed simultaneously. This criterion is repeated at each impact until a collision such that ẏd ∈ 2 (k) takes place. From this moment, according to [START_REF] Zavala-Río | Commande de robots jongleurs[END_REF], the desired value is selected, i.e. ẏ * 1 (k + 1) = ẏd . Notice that the subsets i (k), i = 1, 2, 3, need not be calculated. Just the value of N 1 (k) or N 2 (k) in ( 39) must be computed at each impact.

Example 2. Let us consider the same juggler and the same desired trajectory as in the previous example, but with initial impact conditions given by y(0) = 0, ẏ1 (0) = 1 m/s, and ẏ2 (0) = -0.2 m/s. Assume, this time, that the force is bounded such that |u| ≤ 50 N. In the left-hand side of figure 3, one can see the result of the application of the strategy in corollary 2 (the dashed line indicates the force limit u max = 50 N). Observe that the ball approaches progressively the desired motion, and that the force remains always within its limits. In the right-hand side of figure 3, one can see the evolution of the subspaces

¯ 1 (k) = 1 (k) × 2 (k) (bounded by full lines) and ¯ 2 (k) = 3 (k) × 2 (k) (

bounded by full and dashed lines). The squares represent the values of ẏ *

1 (k + 1) chosen at every collision (according to ẏ2 (k)). The dark dot represents ẏd . The digits stand for the impact number k. Notice that ẏ2 (k) ∈ 3 (k), ∀k ≥ 1.

Figure 3 shows that the application of the foregoing algorithm not only helps maintaining the force within certain limits but also keeps the robot trajectory with small amplitudes. A numerical study showed that, in general, increasing the number of impacts during the transient reduces the limits within which y 2 evolves. In figure 4, the maximum and minimum values of y 2 at each flight time (corresponding to ẏ2 = 0) are plotted versus time for the same juggler in example 2, with initial collision velocities given by ẏ1 (0) = 0.5 m/s < ẏd = 4.43 m/s and ẏ2 (0) = -0.7 m/s. The bounded input algorithm was applied four times, each time increasing the number of collisions during the transient. The different cases were marked, in increasing order of transient shocks, by •, ×, +, * . Notice that the number of impacts in which the desired motion is reached is not predefined by the bounded input control algorithm (or by the designer). However, this one ensures as least as possible transient collisions since [START_REF] Zavala-Río | Commande de robots jongleurs[END_REF] assigns to ẏ * 1 (k + 1) the closest value to ẏd inside 2 (k) at each impact. Hence, in order to increase the number of transient shocks, the following criterion was established: A coefficient β < 1 was defined (recall that ẏ1 (0) < ẏd ). Then, at each collision, if N 2 (k) = ẏd , the actual value of ẏ * 1 (k + 1) used by the control law was β times the value chosen by the algorithm, i.e. ẏ *

1 (k + 1) = β N 2 (k).
Since smaller values of ẏ * 1 (k + 1) are used, more transient shocks take place. The coefficients used to obtain figure 4 were β • = 1, β × = 0.9, β + = 0.8, and β * = 0.7. Observe that increasing the number of transient impacts causes y 2 to remain within a smaller interval. This strategy can be used to keep not only u but also y 2 between two bounds.

Finally, let us notice that we could have defined ẏ2 (0) ∈ 3 (0) as initial condition in point 1 of corollary 2 (instead of ẏ2 (0) ∈ 1 (0)). This would have ensured sequences of the ball impact velocity of the form: 0) is a more restrictive condition compared with ẏ2 (0) ∈ 1 (0) since 3 (0) ⊂ 1 (0). That is why ẏ2 (0) is permitted to be within 1 (0) rather than within 3 (0). In any case, recall that the design criterion in point 2 of corollary 2 ensures ẏ2 (k) ∈ 3 (k) for all k ≥ 1 and the selection of ẏ * 1 (k + 1) > 0 (always positive) for all k ≥ 0. The difference is that if at impact 0 the shock velocities are such that: min{P 3 , P 4 } ẏ1 (0) < ẏ2 (0) ≤ min{M 1 , M 2 , 1} ẏ1 (0) and ẏ1 (0) ≤ ẏd , then the sequence of ball impact velocities is of the form { ẏ1 (k) : that even if ẏ1 (0) = ẏd , the sequence must first go away from its goal, ẏd , and come back to it later. An illustration of this phenomenon is given in figure 5 for a juggler with m 1 = 0.1 kg, m 2 = 1 kg (then m = 0.1), ẏd = 4.43 m/s, y(0) = 0, ẏ1 (0) = 3 m/s < ẏd , ẏ2 (0) = -4.5 m/s, and u max = 22 N.

{ ẏ1 (k) : 0 < ẏ1 (0) < ẏ1 (1) < ẏ1 (2) < • • • < ẏd } if ẏ1 (0) < ẏd , { ẏ1 (k) : ẏ1 (0) > ẏ1 (1) > ẏ1 (2) > • • • > ẏd > 0} if ẏ1 (0) > ẏd , or { ẏ1 (k) : ẏ1 (k) = ẏd > 0, ∀k ≥ 0} if ẏ1 (0) = ẏd . Nevertheless, ẏ2 (0) ∈ 3 (
ẏ1 (0) > ẏ1 (1) > 0, ẏ1 (1) < ẏ1 (2) < • • • < ẏd }; or if: -v 2 g ≤ ẏ2 (0) < min{P 1 , P 2 } ẏ1 (0) and ẏ1 (0) ≥ ẏd , then { ẏ1 (k)} has the form { ẏ1 (k) : 0 < ẏ1 (0) < ẏ1 (1) > ẏ1 (2) > • • • > ẏd > 0}. Notice

Relationship with Buehler-Koditschek's Mapping and Mirror Law

As recalled in [START_REF] Zavala-Río | Hybrid feedback strategies for the control of juggling robots[END_REF], [START_REF] Zavala-Río | Robust control of one degree-of-freedom jugglers[END_REF] the studies in [START_REF] Buehler | Robotic tasks with intermittent dynamics[END_REF], [START_REF] Buehler | A family of robot control strategies for intermittent dynamical environments[END_REF], [START_REF] Buehler | Planning and control of robotic juggling and catching tasks[END_REF] essentially focused on the derivation of a one-dimensional non-linear mapping for the ball velocity. However, only a heuristic control strategy (mirror law) has been proposed in [START_REF] Buehler | Robotic tasks with intermittent dynamics[END_REF], [START_REF] Buehler | A family of robot control strategies for intermittent dynamical environments[END_REF], [START_REF] Buehler | Planning and control of robotic juggling and catching tasks[END_REF]. In this section we show how their results may be used in the framework of the hybrid strategy proposed in ( 17)-( 22), ( 36) and [START_REF] Wang | Dynamics and planning of collisions in robotic manipulation[END_REF]. We also emphasize the relationships between their mirror law and our algorithm.

As we previously pointed out, the signals y * (k + 1) and ẏ * 1 (k + 1) in ( 19)-( 24) can be chosen differently. For instance, let us consider that [START_REF] Vincent | Controllable targets on or near a chaotic attractor[END_REF] holds, but that we choose ẏ *

1 (k +1) as ẏ * 1 (k + 1) = h 1 (w 1 (k)) if h k > y d h 2 • h 1 (w 2 (k)) if h k ≤ y d ( 40
)
where

h 1 (w) = 1 + γ ẏ2 d -w 2 w, h 2 (w) = w 2 + 2g(y d -y(k) -r ), w 1 (k) = ẏ2 1 (k) -2g (y d -y(k)), w 2 (k) = ẏ2 1 (k) -2gr, and 0 < γ < 1 ẏ2 d .
Then, after at most two collisions, the ball impact states will be given by y(k + 1) = y(k) = y d and

ẏ1 (k + 1) = 1 + γ ẏ2 d -ẏ2 1 (k) ẏ1 (k) (41)
Buehler, et al. [START_REF] Buehler | Robotic tasks with intermittent dynamics[END_REF], [START_REF] Buehler | A family of robot control strategies for intermittent dynamical environments[END_REF], [START_REF] Buehler | Planning and control of robotic juggling and catching tasks[END_REF] proposed a different way (the mirror law) to make the ball impact velocity follow the behavior stated by mapping [START_REF] Zavala-Río | On the feedback control of a simple juggling robot[END_REF] for a specific case of the one dof juggler (m = 0, y d = 0), and showed that ( 41) is asymptotically stable [5, p. 72]. Let us now suppose that ẏ * 1 (k + 1) is as in [START_REF] Zavala-Río | Hybrid feedback strategies for the control of juggling robots[END_REF] and that m > 0, but that the scheme is designed with the assumption m = 0, as done in [START_REF] Buehler | Planning and control of robotic juggling and catching tasks[END_REF], [START_REF] Buehler | A family of robot control strategies for intermittent dynamical environments[END_REF], [START_REF] Buehler | Robotic tasks with intermittent dynamics[END_REF]. Then from [START_REF] Rizzi | Preliminary experiments in robot juggling[END_REF] (with m = 0) and ( 59), one obtains

ẏ1 (k + 1) = 1 -m 1 + m + γ 1 + m ẏ2 d -ẏ2 1 (k) ẏ1 (k) (42)
∀k ≥ i ∈ {0, 1, 2}. It follows that the fixed point (y d , ẏd ) is modified to the one of the mapping in [START_REF] Zavala-Río | Robust control of one degree-of-freedom jugglers[END_REF], i.e., y d , ẏ2 d -2m γ . Although the mapping in ( 42) is always stable for a choice of γ ∈ 0, 1 ẏ2 d [5, p. 72] and for all m ∈ [0, 1) such that the fixed point of the real mapping exists (i.e., ẏ2 d -2m γ > 0), the assumption m = 0 may introduce convergence problems and significant tracking errors on the desired velocity and apex if in reality m > 0. Now, as shown previously, the formulation of the robotic juggling problem is divided in two phases: the design of the continuous-time input force and the definition of the impact trajectories. Buehler, et al. [START_REF] Buehler | Robotic tasks with intermittent dynamics[END_REF], [START_REF] Buehler | A family of robot control strategies for intermittent dynamical environments[END_REF], [START_REF] Buehler | Planning and control of robotic juggling and catching tasks[END_REF] proposed a continuous-time surface trajectory consisting of a mirror-like "reflex" of the ball's trajectory with a "distortion" gain (α k ) that changes every flight-time in such a way that (41) takes place, i.e.,

y 2 = -α k y 1 ( 43
)
where

α k = 1-e 1+e + γ 1+e ẏ2 d -ẏ2 1 (k) .
The controller given by ( 17)-( 22), ( 36) and ( 37) follows a similar behavior, but the ball's trajectory is rather reflected by the surface velocity. Indeed from [START_REF] Holmes | The dynamics of repeated impacts with a sinusoidally vibrating table[END_REF], it is not difficult to realize that (18) can be written as

v = -A k g ẏ1 (t) + C k , where C k = ẏ1 (k) g A k + B k .
After integrating the surface dynamics, one gets

ẏ2 (t) = - A k g y 1 (t) + D k + C k (t -t k ) ( 44
)
where

D k = ẏ2 (k)+ y(k) g A k .
Considering the desired trajectories ( 36) and ( 37), one realizes that after at most two impacts, the last term in the right-hand side of (44) becomes dependent just on m. This means that it helps the surface to recover from the perturbations due to impacts. If m = 0, it disappears. The second term is an offset that helps the surface to adjust the impact positions and velocities to the desired ones. If it did not exist and, for example, the desired impact position were zero, the surface collision speed would also be zero which would not make any sense for the generic case e ∈ [0, 1). Actually, choosing (43) constrained Buehler, et al. to define the origin as the only possible desired impact position, while (44) permits the definition of any value for y d or even a sequence of values. Finally, the first term expresses the mirror-like reflex of the ball's trajectory by the surface velocity.

Robustness Analysis

The class of feedback controllers derived in the foregoing section relies on three essential facts: First, it is supposed that the object's postimpact velocities ẏ1 (k) are exactly known.

Second, we assume that there is no disturbance on the ball during the flight-times, so that its motion is perfectly known. Third, the system's parameters are well-known. In practice, it may happen that one or several of these assumptions are not satisfied. Let us consider the following cases: i) ẏ1 (k) and the ball's motion are well-known but the coefficient e is estimated by ê = e. Then the signal ẏ * 2 (k + 1) in ( 22) will be different from the ideal one, but the collisions will occur at the desired position and at t k+1 = t k + d k . ii) ẏ1 (k) or the ball's motion or both are not well-known. Then the collision will not in general occur at the desired height and t k+1 = t k + d k . The main discrepancy between these two situations is that in the first case, the dynamics on the impact Poincaré section is known. Whereas in the second case, we will have to derive the expression of the impact Poincaré map P as a result of the closed-loop dynamics with disturbances. In other words the control law will no longer guarantee the location of the impacts (hence of the section ). In general the impact Poincaré map will be impossible to calculate explicitly.

Unknown Restitution Coefficient

Let us suppose that we compute ( 17)-( 22) using an estimate of the restitution coefficient ê which value is different from the real e, i.e., ê = e, but we still assume that the postimpact velocities are measured correctly. Then, after integrating the dynamics ( 5), ( 6) considering ( 17)- [START_REF] Rizzi | Further progress in robot juggling: the spatial two-juggle[END_REF], it is easy to realize that the next impact position does not change, i.e., y 1 (k +1) = y 2 (k + 1) = y * (k + 1). Moreover, the velocity of the surface just before the next impact will be the calculated one, i.e., ẏ2 (t - k+1 ) = ẏ * 2 (k + 1), which we would compute as (see [START_REF] Rizzi | Preliminary experiments in robot juggling[END_REF])

ẏ * 2 (k + 1) = 1 + m 1 + ê ẏ * 1 (k + 1) + m - ê 1 + ê ẏ2 1 (k) -2g (y * (k + 1) -y(k)) (45)
On the contrary, the system's impact velocities are the states that are affected by this discrepancy ( ê = e). Indeed, let us analyze the behavior of the ball impact velocities. From (59), (45), and the integrated dynamics of the ball ( 12), ( 13), one obtains

ẏ1 (k + 1) = 1 + e 1 + ê ẏ * 1 (k + 1) + e - ê 1 + ê ẏ2 1 (k) -2g(y * (k + 1) -y(k)) ( 46
)
Let us for the moment consider the case in which h k > y d . Then, from ( 36), [START_REF] Wang | Dynamics and planning of collisions in robotic manipulation[END_REF] in (46), we get

ẏ1 (k + 1) = 1 + e 1 + ê ẏd + e - ê 1 + ê ẏ2 1 (k) -2g(y d -y(k)) (47)
Assuming that the difference e -ê is such that h k+ j > y d , ∀ j ≥ 1, we have for all j ≥ 2

ẏ1 (k + j) = 1 + e 1 + ê ẏd + e - ê 1 + ê ẏ1 (k + j -1)
Then by induction, one can verify that

ẏ1 (k + j) = 1 - e - ê 1 + ê j ẏs + e - ê 1 + ê j ẏ2 1 (k) -2g(y d -y(k)) (48)
for all j ≥ 1, where ẏs = 1+e 1+2 ê-e ẏd . In the same way, for the case in which h k ≤ y d , we get for all j ≥ 2

ẏ1 (k + j) = 1 - e - ê 1 + ê j-1 ẏs + e - ê 1 + ê j-1 ẏ2 1 (k + 1) -2g(y d -y(k) -r ) (49)
It is easy to verify that e-ê 1+ê < 1, ∀(e, ê) ∈ [0, 1] × [0, 1] : (e, = (1, 0) ( 2 ). From (48) and (49), one realizes that whatever the relationship between h k and y d is, the impact velocity of the ball converges to ẏs , i.e., lim j→∞ ẏ1 (k + j) = ẏs . Following the same analysis for the impact velocity of the surface, we find that the system's collision states converge to x s = (y d , y d , ẏs , 1-e-2em 1+2 ê-e ẏd ). It is noteworthy that the convergence is ensured assuming that h k+ j > y d , ∀ j ≥ 1. From a detailed analysis, one can realize that a sufficient condition for this to be true is the subestimation of the restitution coefficient, i.e., ê < e. Indeed, let us again consider the case in which h k > y d . Since y(k + 1) = y * (k + 1), from (36) y(k + 1) = y d . Hence, ẏ1 (k + 1) must be positive for h k+1 > y d to be satisfied. From (47), it is easy to see that if e -ê > 0 then ẏ1 (k + 1) > 0. Continuing this analysis for j ≥ 2, one realizes from (48) that if e -ê > 0 then ẏ1 (k + j) > 0, ∀ j ≥ 2, and hence h k+ j > y d , ∀ j ≥ 1. For the case in which h k ≤ y d , if h k+1 > y d then it is easy to verify that the same analysis can be made for j ≥ 2 (see ( 49)) and that it leads to the same result. It can also be proved that for h k+1 > y d to be true, e > ê is also a sufficient condition. Finally, for the case in which ê > e, a small enough difference êe is desirable.

Errors in the Ball's Velocity Measurements

In order to show what happens when the impact velocities of the ball are not well measured, let us consider the same one dof juggler of Example 1, with the same desired trajectory (y d , ẏd ) (and initial conditions of case a). Several simulations were executed introducing a noise signal to the measured impact speed according to the following expression: ν1 (k) = (1-α +2αρ) ẏ1 (k), where ν1 (k) is the measured collision velocity, ρ is a random signal with values in [0, 1], and α ∈ [0, 1) is an indicator of the noise level, such that ν1 (k) is randomly within [(1 -α) ẏ1 (k), (1 + α) ẏ1 (k)], for all k ≥ 0. α was varied from 0 to 0.2 every 0.01 units. Each simulation includes around 30 impacts and they were repeated 10 times for each noise level. An average radius x(k)x d 2 was obtained for each simulation. The result is shown in figure 6. We conclude that the error distance is in general bounded and proportional to the noise level (despite the existent dispersion). The impact states remain in a neighborhood around the desired fixed point and a destabilisation of the closed-loop system does not in general occur.

Computed Impact Velocities and Restitution Coefficient Estimation

In practice, it is difficult to have exact impact velocity measurements and they are usually noisy [START_REF] Buehler | Robotic tasks with intermittent dynamics[END_REF], [START_REF] Buehler | Planning and control of robotic juggling and catching tasks[END_REF]. On the contrary, detecting the positions where the ball's velocity is zero (apex) s(k) is in general much easier. From this and the fact that the free-motion and impact dynamics are well-known, we propose a noiseless method for the knowledge of the ball's velocities (we suppose that the rest of the impact states are exactly measurable). It consists of the computation of this impact state using apex measurements. Indeed, from (59) we get an expression for the calculation of ẏ1 (k + 1) given by

ẏ1 (k + 1) = m - ê 1 + m ẏ1 (t - k+1 ) + 1 + ê 1 + m ẏ2 (t - k+1 )
where we have replaced the real Newton restitution coefficient e by an estimation ê. The term ẏ1 (t - k+1 ) can be calculated from the free-motion dynamics of the ball, the desired impact trajectories and the knowledge of s(k). A disadvantage of this method is the possibility of having ê = e. This case is analyzed via simulation for the juggler previously described in Example 1, with the same desired trajectory (y d , ẏd ) (and initial conditions of case a). Figure 7 shows the results obtained for: a) (e, ê) = (0.7, 0.5), and b) (e, ê) = (0.7, 0.9). We can observe that a periodic motion is always attained, i.e., lim k→∞ x(k) = x s . For the case ê > e, the ball impact states are characterized by y s < y d and ẏs < ẏd . On the contrary, when ê < e, it can be seen that y s > y d and ẏs > ẏd . As expected, an underestimation of e makes the robot hit the ball with a too large velocity. An overestimation implies the opposite effect. Notice (compare figures 7 and 2) that an error of 30% on e results in an error of 500% of the apex for the underestimation (remember that the desired apex was defined above as s d = 1 m) and 70% for the overestimation. One therefore concludes that the accurate knowledge of e is very important.

In order for the impact states to attain asymptotically the desired fixed point (x d , d ), a heuristic adaptive law was implemented for ê in terms of the apex error:

ê(k + 1) = ê(k) + a s(k)-s d s d
, where a ∈ R is the adaption gain. Figure 8 shows (considering the same juggler as above) the evolution of the errors y(k)y d , s(k)s d , and ê(k)e, for a = 0.025, e = 0.7 and: a) ê(0) = 0.5, and b) ê(0) = 0.9. These and several other numerical results have shown that for |ê(0) -e| and an adaption gain a both small enough, convergence of (x(k), ê(k)) → (x d , e) occurs. Object dynamics with damping. Finally, the robustness of the control scheme with respect to unknown damping acting on the object during the flight times has been studied in [START_REF] Zavala-Río | Commande de robots jongleurs[END_REF], [START_REF] Zavala-Río | Hybrid feedback strategies for the control of juggling robots[END_REF]. Since the impact Poincaré map is impossible to obtain explicitely in that case, numerical simulations have been presented. They show that for small enough (but strictly larger than zero) damping and the same desired trajectories as in example 1, the ball closed-loop trajectories still converge towards a periodic motion with one impact per period.

Conclusions

A family of feedback controllers is presented for one degree-of-freedom juggling robots. Due to the potential applications of such study in non-prehensile manipulation, in which the mass of the object is not necessarily negligible with respect to that of the robot, a 4-dimensional model of the system has been considered. The hybrid control strategy mainly consists of a discrete-time state feedback. It is flexible enough to cope with various juggling tasks and can be designed so that the impact position and velocity are free (and is not restricted to periodic trajectories). A strategy is presented for the case when the input is bounded. Robustness for this class of hybrid strategies is also investigated both numerically and analytically for different types of disturbances and uncertainties (measurement noise, bad knowledge of the restitution coefficient, unknown damping on the object during the flight-times). The possible extension of the presented method to a larger class of mechanical systems with unilateral constraints is also outlined.

A. Proof of Claim 1

The proof is divided in two parts: In Part 1 we treat the generic case in which e ∈ (0, 1] and an initial impact takes place at t 0 ≥ 0, i.e., y 1 (t 0 )y 2 (t 0 ) = 0 and ẏ1 (t + 0 ) -ẏ2 (t + 0 ) > 0. In Part 2 we suppose that the two bodies are in contact at a generic instant t k (in which the control input ( 17)-( 24) is applied), i.e., y 1 (t k )y 2 (t k ) = ẏ1 (t k ) -ẏ2 (t k ) = 0. The objective of Part 2 is to prove that points 1-4 of claim 1 are true also for the perfectly plastic case, i.e., e = 0, and for the case in which the bodies are in contact at t 0 (for any value of e ∈ [0, 1]).

Part 1. Let us for the moment consider the isolated dynamics of the ball (5) assuming initial conditions (y 1 (t k ), ẏ1 (t k )) = (y(k), ẏ1 (k)). Then, it is easy to verify that the position and velocity after a time interval d k as defined in [START_REF] Peshkin | Planning robotic manipulation strategies for workpieces that slide[END_REF] for any y * (k + 1) satisfying [START_REF] Rizzi | Progress in spatial robot juggling[END_REF]. Notice that [START_REF] Rizzi | Progress in spatial robot juggling[END_REF] determines the set of reachable positions of the ball given specific initial states. This condition not only ensures that the term in the square roots of ( 21) and (50) (see also [START_REF] Rizzi | Preliminary experiments in robot juggling[END_REF]) is non-negative but also ensures positive values of d k . Complex or negative time intervals would not make any sense.

f ( ) > 0, ∀ ∈ (0, d k ), which demonstrates point 1 of the claim. This and (55) prove points 2 and 3. Finally, from ( 5)- [START_REF] François | Contribution à la locomotion articulée dynamiquement stable[END_REF] Notice that this may be the result of two possible situations: On the one hand, when the restitution coefficient is zero, i.e., e = 0, (62) takes place at impact instants t k (see [START_REF] François | Contribution à la locomotion articulée dynamiquement stable[END_REF]). In the other hand, when the two bodies are in contact at the initial instant t 0 , (62) is verified with k = 0. Notice that this last situation does not depend on the value of e. Considering this assumption, and from ( 1)-(3), the system dynamics becomes where λ > 0 is the interaction force between the bodies. This dynamics holds until the moment t in which the interaction force vanishes and the relative acceleration becomes positive, i.e., λ(t ) = 0 and ÿ1 (t ) -ÿ2 (t ) > 0. From (63) and (64), it is easy to verify that this condition is translated to: ÿ1 (t ) -ÿ2 (t ) = 1 m 1 + 1 m 2 λ(t ) -u(t ) m 2 = -u(t ) m 2 > 0, i.e., u(t ) < 0. Considering the input force ( 17)-( 24) taking into account (62), we can prove that u(t + k ) = m 2 (g + B k ) < 0 (see [START_REF] Zavala-Río | Commande de robots jongleurs[END_REF], app. B, sect. B.2). We conclude that t = t + k , that is, the bodies loose contact at the immediate instant in which the control input ( 17)-( 24) is applied (just after the shocks for the perfectly plastic case). The dynamics in ( 5)-( 8) is recovered at t = t = t + k . Consequently, Part 1 of the proof holds for this case too.

Notes

1. For instance the work in [START_REF] François | Contribution à la locomotion articulée dynamiquement stable[END_REF] for stabilization of hopping robots strongly relies on the choice e = 0 which allows them to apply a control during the constrained motion phase. It is not clear how the results can be extended to the case e > 0. 2. It would be indeed difficult to think of the case in which the collisions are perfectly elastic while one tries to model them as perfectly plastic.

Figure 1 .

 1 Figure 1. One dof juggler.

Figure 2 .

 2 Figure 2. Results of example 1: a) h 0 > y d ; b) h 0 ≤ y d .

  ) to a one dof juggler with e = 0.7, m 1 = 0.3 kg, m 2 = 2 kg (which gives m = 0.15), y d = 0, ẏd = 4.43 m/s (which is equivalent to s d = 1 m), y(0) = -0.4 m, ẏ2 (0) = -0.2 m/s. In figure 2, one can observe the trajectories followed by the ball and the surface in the upper graphs and the control input force in the lower ones (where the dashed line indicates the surface weight, m 2 g = 19.62 N). Two cases are considered: a) ẏ1 (0) = 3.4 m/s which gives h 0 = 0.2 m > 0, and b) ẏ1 (0) = 2 m/s which gives h 0 = -0.2 m < 0. This shows the two possible initial cases: h 0 > y d and h 0 ≤ y d respectively.

Figure 3 .

 3 Figure 3. Bounded input force strategy.

Figure 4 .

 4 Figure 4. Decreasing the robot trajectory amplitudes.

Figure 5 .

 5 Figure 5. Bounded input strategy: { ẏ1 (k) : ẏ1 (0) > ẏ1 (1) > 0, ẏ1 (1) < ẏ1 (2) < • • • < ẏd }.

Figure 6 .

 6 Figure 6. Average of x(k)x d 2 for different noise levels.

Figure 7 .

 7 Figure 7. Computed impact velocities: a) ê < e; b) ê > e.

Figure 8 .

 8 Figure 8. Adaption errors: a) ê(0) < e; b) ê(0) > e.

  are

(y 1

 1 (t k + d k ), ẏ1 (t k + d k )) = y * (k + 1), -ẏ2 1 (k) -2g (y * (k + 1)y(k)) (50)

m 1

 1 ÿ1 + m 1 g = λ (63) m 2 ÿ2 + m 2 g = u -λ (64) f (y 1 , y 2 ) = y 1y 2 = 0 (65)

  after a possible transient period. Notice that this is equivalent to impose y 1 (t k ) = y d and t k+1 = t k + d , with d = 2 g ẏd , or y 1 (t k ) = y d and the trajectory apex s(k) = s d = y d +

	ẏ2 d 2g . Nevertheless, as we shall see later, one is
	also free to choose varying signals y d

  it follows that at a generic impact k + 1 Part 2. Let us now suppose that at a generic instant t k , both bodies are in contact and tend to move together, i.e.,y 1 (t k )y 2 (t k ) = ẏ1 (t k ) -ẏ2 (t k ) = 0 (62)

	ẏ1 (k + 1) =	m -e 1 + m	ẏ1 (t -k+1 ) +	1 + e 1 + m	ẏ2 (t -k+1 )	(59)
	ẏ2 (k + 1) =	m(1 + e) 1 + m	ẏ1 (t -k+1 ) +	1 -em 1 + m	ẏ2 (t -k+1 )	(60)
	Considering (50) and (54) in (59), one gets
	ẏ1 (k + 1) = -	m -e 1 + m		ẏ2 1 (k) -2g (y * (k + 1) -y(k)) +	1 + e 1 + m	ẏ * 2 (k + 1)	(61)
	Introducing (22) into (61), point 4 follows.

Acknowledgments

A. Zavala Río was supported by CONACYT-Mexico. Part of this work was performed while he held a post-doctoral fellowship at the Mechanical Engineering Laboratory, AIST, Japan.

Now, let us consider the isolated dynamics of the surface. From ( 6) and [START_REF] Kolmanovsky | Hybrid feedback laws for a class of cascade nonlinear control systems[END_REF], one obtains ÿ2 = v. This system is easily written as

for suitable (F, G) and z. Then, the control input

where

and (t) = exp(Ft) are respectively the Grammian and the state-transition matrix of (51), drives the state z from z k to z k+1 in the time interval d k (see [START_REF] Sontag | Mathematical Control Theory, Deterministic Finite Dimensional Systems[END_REF]). Such controllers are sometimes called Dead-beat algorithms [START_REF] Kreisselmeier | Adaptive control of continuous-time overmodeled plants[END_REF]. 52), we get ( 18)-( 20) as the input that leads the surface states from (y

Notice that the d k used in [START_REF] Lynch | Stable pushing: mechanics, controllability, and planning[END_REF] and ( 20) is designed using [START_REF] Peshkin | Planning robotic manipulation strategies for workpieces that slide[END_REF] which makes the final position of the isolated trajectory of the surface coincide with the one of the ball. In other words, 50) and (54) satisfy

Let us now consider the whole dynamics ( 5)-( 8) assuming initial conditions such that y 1 (k) ≡ y 2 (k) ≡ y(k) and

Notice that these are in general the conditions that result of an impact when e ∈ (0, 1]. Introducing ( 17), ( 18) in ( 6) and integrating two times ( 5) and ( 6) for t ∈ (t k , t k + d k ), we get

where = tt k . It is not difficult to verify that (57) can be rewritten as

where Y ( )

From [START_REF] Lynch | Stable pushing: mechanics, controllability, and planning[END_REF] and [START_REF] Masri | On the stability of the impact damper[END_REF], and considering conditions (56) and [START_REF] Rizzi | Further progress in robot juggling: the spatial two-juggle[END_REF], one can verify that Y ( ) > 0, ∀ ∈ (0, d k ). From this fact, one easily sees that


A Virazel
email: virazel@lirmm.fr

R David
email: rene.david@lag.ensieg.inpg.fr

P Girard
email: girard@lirmm.fr

C Landrault
email: landraul@lirmm.frpravo@lirmm.fr

S Pravossoudovitch

H J Wunderlich

Delay Fault Testing: Choosing Between Random SIC and Random MIC Test Sequences

Keywords: delay testing, random testing, robust test, non-robust test, BIST

The combination of higher quality requirements and sensitivity of high performance circuits to delay defects has led to an increasing emphasis on delay testing of VLSI circuits. In this context, it has been proven that Single Input Change (SIC) test sequences are more effective than classical Multiple Input Change (MIC) test sequences when a high robust delay fault coverage is targeted. In this paper, we show that random SIC (RSIC) test sequences achieve a higher fault coverage than random MIC (RMIC) test sequences when both robust and non-robust tests are under consideration. Experimental results given in this paper are based on a software generation of RSIC test sequences that can be easily generated in this case. For a built-in self-test (BIST) purpose, hardware generated RSIC sequences have to be used. This kind of generation will be shortly discussed at the end of the paper.

Introduction

Delay fault testing allows to test for delay faults. A delay fault occurs in a circuit when one or more paths in the circuit fail to propagate a signal within the time interval specified by the clock period. Detection of delay faults requires two-pattern tests. An initialization vector is applied and the circuit is allowed to stabilize. Then, the test vector is applied and the circuit outputs are sampled at clock speed. The response is then compared to that of the fault-free circuit to determine the presence or the absence of a delay fault. Hence, correct operation of a circuit at the intended speed can only be guaranteed if there is no delay fault in the circuit.

With the continuous increase in the operating speed of VLSI circuits, delay fault testing is likely to become industrially accepted in the near future [START_REF] Majhi | Tutorial: Delay Fault Models and Coverage[END_REF]. With delay fault testing (i.e. application of two-pattern tests to ensure temporal correctness of the design), at-speed testing (i.e. test at the intended operating speed of the circuit) is also becoming an essential part of the verification process of today's VLSI circuits since it allows to optimize the test time and provides the means to test for delay faults. From a general point of view, at-speed testing does not necessarily lead to delay fault testing, but delay fault testing needs to be performed at-speed to detect timing defects and test the performance of the circuit.

A problem that occurs with currently-used external testers, however, is that they are several times slower than the designs they have to test. Purchasing high speed testers that meet the performances of new designs requires a huge investment. Moreover, even with those high speed testers, it is not always possible to have a timing accuracy comparable to the IC internal speed [START_REF] Zorian | Testing the Monster Chip[END_REF]. In this context, BIST represents an attractive test solution since it allows at-speed testing of the circuit under test, thus solving timing accuracy and test time related problems encountered with traditional external testers. In addition, BIST drastically reduces the amount of test data exchanged with the tester, thus reducing the need for complex external testing equipment. BIST can hence be run on a very low cost tester. Finally, BIST solves the problem of tester capacity (very often, external testers do not have enough memory to store the entire test set to cover stuck-at, transition and path delay faults [START_REF] Hetherington | Logic BIST for Large Industrial Designs: Real Issues and Case Studies[END_REF]) and the problem of low accessibility of internal nodes of the design, that increases the test complexity [START_REF] Wunderlich | Built-In Self Test (BIST): Synthesis of Self-Testable Systems[END_REF].

Although the work presented in this paper comes within the general frameworks of delay testing and random testing, the targeted application is logic BIST for delay faults. BIST is generally based on pseudorandom testing [START_REF] Hetherington | Logic BIST for Large Industrial Designs: Real Issues and Case Studies[END_REF]. Pseudo-random testing refers to the application of test patterns that exhibit randomness, but which are generated using special-purpose hardware (LFSR or Cellular Automata), and are thus repeatable. In terms of delay fault coverage, conventional pseudo-random test patterns in which more than one bit change between two consecutive patterns are not efficient to robustly test combinational circuits in a reasonable test time [START_REF] Wang | Weighted Random Robust Path Delay Testing of Synthesized Multilevel Circuits[END_REF]. This fact is what has motivated the development of pseudo-random BIST techniques in which Single Input Change (SIC) test pairs are generated for testing delay faults [START_REF] Girard | An Optimized BIST Test Pattern Generator for Delay Testing[END_REF][START_REF] Wang | Weighted Random Robust Path Delay Testing of Synthesized Multilevel Circuits[END_REF]. SIC test pairs are sufficient to detect all robustly detectable path delay faults [START_REF] Smith | Model for Delay Faults Based upon Paths[END_REF], with a test length shorter than that required with Multiple Input Change (MIC) test pairs [START_REF] Wang | Weighted Random Robust Path Delay Testing of Synthesized Multilevel Circuits[END_REF]. Note that robustness of delay tests is important to guarantee timing correctness of the CUT since it allows to detect a delay fault even in the presence of other delay faults in the circuit. Additionally, as not all delay faults have a robust test, the effectiveness of SIC test pairs (compared with MIC test pairs) has also been observed in the case of a fault detected by a validatable non-robust test (a non-robust test validated by another test in the sequence) [START_REF] Crepaux-Motte | An Algebraic Method for Delay Fault Testing[END_REF].

Several authors have studied the performance of test sequences which are SIC but not truly random. In [START_REF] Craig | Pseudo-Exhaustive Adjacency Testing: A BIST Approach for Stuck-Open Faults[END_REF][START_REF] Girard | An Optimized BIST Test Pattern Generator for Delay Testing[END_REF][START_REF] Wang | Weighted Random Robust Path Delay Testing of Synthesized Multilevel Circuits[END_REF], after each generation of a classical pseudorandom vector, a string of SIC vectors is applied to the CUT. In these approaches, a high robust fault coverage is obtained thanks either to a weighted random generation [START_REF] Wang | Weighted Random Robust Path Delay Testing of Synthesized Multilevel Circuits[END_REF] or to an efficient shifting mechanism [START_REF] Girard | An Optimized BIST Test Pattern Generator for Delay Testing[END_REF]. Another approach, proposed in [START_REF] Glover | A Deterministic Approach to Adjacency Testing for Delay Faults[END_REF], consists in using a deterministic sequence for the stuckat faults. From this set of deterministic patterns, the authors produce a SIC sequence which is very effective for the path delay fault testing. From a general point of view, these papers concentrate their analysis on the coverage of faults having at least one robust test.

Our paper differs from those previous studies into two points. Firstly, we are interested in measuring the performance of a truly Random SIC (RSIC) sequence to test path delay faults, in comparison with the performance of a truly Random MIC (RMIC) sequence. Secondly, we aim at taking into account all the delay faults, even those that do not have a robust test. A fault with no robust test is worth being tested too! (a circuit can be synthesized such that all path faults are robustly testable [START_REF] Pramanick | On Multiple Path Propagating Tests for Path Delay Faults[END_REF], but most of the existing circuits do not have this property).

In practice, a truly RSIC sequence is difficult to produce since any generation (software or hardware) is pseudo-random in nature. However, a "careful" generation allows to obtain test sequences whose randomness is very good (we mean that they have roughly the same properties than a pure random generation). Since a good generation is easier by software than by hardware [START_REF] David | Random Testing of Digital Circuits: Theory and Applications[END_REF], the results in this paper are based on a software generation briefly explain in Section 2. Hardware generation (required for BIST) will be shortly discussed by the authors at the end of the paper.

The rest of the paper is organized as follows. Section 2 gives some preliminary definitions on delay fault testing. Section 3 analyses the delay fault coverage of RSIC test sequences on a case study, and Section 4 generalizes the results obtained to other circuits and other fault models. Concluding remarks and future work are discussed in Section 5.

Preliminaries

Correct operation of a circuit at the intended speed requires that any path delay exceeds the value determined by the clock period. This is usually verified by delay testing, using the path delay fault model [START_REF] Smith | Model for Delay Faults Based upon Paths[END_REF]. In this model, it is assumed that the presence of a delay fault increases the delay along the path. This model represents distributed delays in the circuit, often caused by device parameter variation [START_REF] Sivaraman | Diagnosis of Parametric Path Delay Faults[END_REF], as well as single isolated failures. This is the main advantage of this model over other existing delay fault models, namely the gate delay fault model and the transition delay fault model [START_REF] Abraham | Testing of MOS VLSI Circuits[END_REF][START_REF] Koeppe | Modeling and Simulation of Delay Faults in CMOS Logic Circuits[END_REF]. However, an important feature of the path delay fault model is that the single fault assumption is not realistic since a single defect usually affects a large number of paths. For this reason, a robust test is preferred to detect a path delay fault. A robust test is a test that detects a delay fault regardless of all other delays and delay faults in the circuit [START_REF] Smith | Model for Delay Faults Based upon Paths[END_REF]. In contrast, a test that detects a fault with the assumption that no other delay fault can exist in the circuit is called a non-robust test.

A test for a path delay fault consists of propagating a transition along the target path P = (g 0 , g 1 , . . . , g n ), where each g i is a gate except g 0 and g n which are the source and the destination of the path respectively (usually an input and an output of the circuit). The inputs of g i other than the output of g i -1 are called side-inputs of g i (or side-inputs of P). Each connection between g i -1 and g i is called an on-input of P. Sensitization conditions for non-robust and robust testability of a path delay fault are as follows. A non-robust test for a rising (or a falling) transition on a path P sets every side-input of P to a final non-controlling value. In a robust test, side-inputs must be at stable non-controlling values (with no static hazards) when on-path inputs have transitions to final controlling values. These conditions are illustrated in Fig. 1 for an AND-type logic gate.

As the purpose in this work was not to develop a delay fault simulator to evaluate the fault coverage of the experimented delay test sequences (random SIC and MIC test sequences), we used an existing industrial test evaluation package, TestGen of Synopsys [START_REF] Testgen | version Tg3.0.2 User Guide[END_REF], to perform test validations. In the delay fault tool suite of this package, several models are supported for testing delay paths. The conventional non-robust delay fault model is referred to as the weak non-robust delay model in TestGen, and delay fault simulation or test generation from this model is possible. Conversely, the conventional robust delay fault model is also supported by TestGen, but one of the main sensitization constraints in this model (no glitches on the side-inputs when onpath inputs have transitions to final controlling values) is not verified during delay fault simulation or test generation. As a consequence, fault-free results cannot be ensured when using such delay fault model with TestGen.

The closest delay fault model handled by this tool is the strong non-robust delay fault model, for which the sensitization conditions are the following: in a strong non-robust test, side-inputs must be at initial and final non-controlling value when on-path inputs have transitions to final controlling values (the only difference with the conventional robust delay model states in the acceptance of static hazards on side-inputs). We used the strong non-robust delay model (in addition to the non-robust delay model) as a metric to evaluate the quality of the test sequences generated for testing path delay faults. As the strong non-robust delay model is very close to the conventional robust delay model, this model will be referred to as the pseudo-robust delay model in the sequel.

Apart from robust and non-robust testable paths, two other classes of paths have been defined in the literature: functional sensitizable paths and functional redundant paths [START_REF] Lam | Delay Fault Coverage, Test Set Size, and Performance Trade-Offs[END_REF]. Functional redundant paths can never determine the performance of the circuit and do not have to be tested [START_REF] Cheng | Delay Testing for Non-Robust Untestable Circuits[END_REF]. On the other hand, defects on functional sensitizable paths may degrade the circuit performance when several path delay faults occur simultaneously. Although the number of functional sensitizable paths in a circuit may be not negligible, this model of path delay faults is not handled by TestGen. For this reason, only pseudo-robust and non-robust tests are considered in the rest of this study.

In general, two-pattern tests may vary in multiple bit positions. In this case, they are called multiple input change (MIC) pattern pairs. Test pattern pairs that differ in exactly one bit are called adjacent or single input change (SIC) pattern pairs. Let us now define what a RMIC and a RSIC sequences should be from a theoretical point of view (we assume implicitly the case of equal likelihood of all vectors). Let S = V (1) V (2) . . . V (l) . . . V (L) be a test sequence composed of L successive n-bit vectors V (l). Each vector takes a value from the set V = {V 0 , V 1 , . . . , V j , . . . , V 2 n -1 }, where V j corresponds to the n-bit vector associated with the decimal value j. For example, for n = 5, V 9 = 01001, i.e., x 5 = x 3 = x 2 = 0 and x 4 = x 1 = 1. In a RMIC sequence, the probability Pr[V (l) = V j ] = 1/2 n for any l and any j, and the probability

Pr[V (l) = V j ] is independent of the values V (i), i = 1, . . . , l -1. In a RSIC sequence, Pr[V (l) = V j |V (l -1) = V k ] = 1/n iff | j -k| = 2 a ,
where a is a non-negative integer (in other words, V (l) differs from V (l -1) by exactly one bit randomly drawn, V (l + 1) differs from V (l) by exactly one bit, etc).

The RSIC sequence is software generated as follows. We use an instruction "random", drawing a (pseudo) random number Y uniformly distributed in the range [0,1). Let V (l) = x 1 (l)x 2 (l) . . . x n (l) the vector at time l. If the random value Y (l + 1) is drawn, then the input variable x i such that i = Y (l + 1) • n is changed, i.e., x i (l + 1) = x i (l) and x j (l + 1) = x j (l) for any j = i (Section 10.2 in [START_REF] David | Random Testing of Digital Circuits: Theory and Applications[END_REF]). The principle of the hardware generation is similar, the vector Y is obtained from a LFSR.

A delay test consisting of adjacent or SIC vectors is called an adjacency or an asynchronous test. As a single transition is applied at the primary inputs of the CUT in an adjacency test, the probability of delay test invalidation due to hazards or multiple delay faults is greatly reduced. This is one of the main reasons for using such kind of tests. Moreover, SIC test pairs are sufficient to detect all robustly detectable path delay faults [START_REF] Smith | Model for Delay Faults Based upon Paths[END_REF]. Finally, the universe of pattern pairs considered for SIC test generation (O(n • 2 n )) is significantly smaller than that for MIC test generation (O(2 2n )). Hence, SIC fault coverage can be higher than MIC fault coverage for the same test length [START_REF] Crepaux-Motte | An Algebraic Method for Delay Fault Testing[END_REF][START_REF] Wang | Weighted Random Robust Path Delay Testing of Synthesized Multilevel Circuits[END_REF].

Delay Fault Coverage-A Case Study

This case study corresponds to experiments performed on the combinational part of circuit s382 of the ISCAS'89 benchmark set [START_REF] Brglez | Combinational Profiles of Sequential Benchmark Circuits[END_REF]. The performance of the random SIC and MIC test sequences is evaluated based on the results obtained by a deterministic ATPG (Test-Gen Tg3.0.2 of Synopsys [START_REF] Testgen | version Tg3.0.2 User Guide[END_REF]). For further comparison purpose, the results obtained from this ATPG are given in Table 1 for some of the ISCAS89 circuits (combinational part). Notations in Table 1 are as follows. Fc SA, Fc PR and Fc NR are the stuck-at fault coverage, the pseudo robust delay fault coverage and the non-robust delay fault coverage respectively. Ls and Ld are the test length of stuck-at fault and delay fault test sequences respectively. These lengths are the test lengths required to achieve the maximum fault coverage.

Table 1 should be interpreted as follows. Circuit s382 (its combinational part) has 24 inputs. All the stuck-at faults are detected with a test sequence of length 36. The delay test sequence is composed of 1398 input vectors. When this sequence is applied to the CUT, 88% of the delay faults are pseudo-robustly (PR) tested (12% of the faults have no robust test), 91.75% are non-robustly (NR) tested, and 8.25% are functionally sensitizable or redundant [START_REF] Cheng | Delay Testing for Non-Robust Untestable Circuits[END_REF]. As the set of faults covered by NR tests includes the set of faults covered by PR tests, there are exactly (91.75 -88.00) = 3.75% of delay faults that have at least one NR test but no PR test.

In order to measure the performance of a test sequence, let us denote by Eff (Efficiency) the ratio of faults detected by the sequence over the maximum number of testable faults. For example, the pseudorobust efficiency of a sequence S that detects 44% of the PR testable faults in circuit s382 is 44/88 = 0.5 (i.e. 50%). If it detects 88% of the PR testable faults, its PR efficiency is Eff = 100%. Table 2 given above presents simulation results obtained from one RMIC and one RSIC sequences. Note that results would be quite similar with others RMIC or RSIC sequences. In Table 2, the results are given for various test lengths. For example, for Ld = 1398 vectors, the following results are obtained: 1) RMIC sequence: a PR test has been found for 33.24% of the PR testable faults and a NR test has been found for 91.01% of the NR testable faults. 2) RSIC sequence: a PR test has been found for 52.84% of the PR testable faults and a NR test has been found for 57.22% of the NR testable faults.

The graph in Fig. 2 is a graphical representation of the results given in Table 2.

The following observations can be made:

1) The PR efficiency of the RSIC sequence is close to the NR efficiency, while there is a huge gap between these efficiencies for the RMIC sequence. The reason is that, due to the adjacency of successive vectors, most of the tests found in the RSIC sequence are robust. 2) For a short test length (not interesting in practice), the fault efficiencies achieved by the RSIC test sequence are lower than those achieved by the RMIC sequence. For average and long test lengths, the NR efficiency of the RMIC sequence is much higher than that of the RSIC sequence, while the RSIC sequence is more efficient in terms of PR efficiency. 3) For very long test lengths, the RSIC test sequence may be as efficient as the RMIC sequence regarding the NR efficiency (last line in Table 2). It is much more efficient regarding the PR efficiency.

To summarize, we can say that the RMIC sequence achieves a better NR efficiency while the RSIC sequence leads to a better PR efficiency (for average and long test lengths). Consequently, a question should be asked: what is the best test sequence to be used when both PR efficiency and NR efficiency are targeted? Let us then try to find a single criterion to help selecting the best test sequence.

Let us remind that, for circuit s382, 88% of the faults have at least one PR test and 91.75% of the faults have at least one NR test (Table 1). In other words, 91.75% of the delay faults are NR testable and q = 88/91.75 = 0.959, i.e. 95.9% of the NR testable faults have at least one PR test. Let us try now to answer the following question: what is the fraction of NR testable faults which is effectively detected by a RMIC sequence of length Ld = 1398?

According to Table 2, a PR test has been found for 33.24% of the faults having at least one PR test. Since the fraction of NR testable faults having at least one PR test is 0.959, a PR test has been found for 33.24 × 0.959 = 31.88% of the NR testable faults. According to the basic property of a robust test (detection of a fault irrespective of the presence of other delay faults), at least 31.88% of the NR testable faults will Since we have no answer to this question, assumptions can be stated. Let us assume that the success rate is S NR = 0.5 (i.e., 50%), then the weighted efficiency for our example is 31.88 + (91.75 -31.88) × 0.5 = 61.82. Let us now remind the notation and generalize the calculation: q = (Fc PR/Fc NR) is the fraction of NR testable faults having at least one PR test, S NR is the assumed value of the success rate of NR tests, and Eff (S NR ) is the weighted efficiency of a test sequence for a given S NR . The weighted efficiency of a test sequence for a given success rate S NR is given by Eq. ( 1) described below: [START_REF] Abraham | Testing of MOS VLSI Circuits[END_REF] The application of Eq. ( 1) to circuit s382 for the RMIC and RSIC test sequences and for S NR = 10%, 50% and 90% are presented in Table 3 for various test lengths. These results are also presented in the graph of Fig. 3. In this figure, it clearly appears that for a long test length, the weighted efficiency is higher for the RSIC sequence than for the RMIC sequence, for all assumed values of the success rate S NR .

Eff (S NR ) = Eff PR × q + (Eff NR -Eff PR × q) × S NR
As a conclusion, these results demonstrate that, even with a lower non-robust delay fault coverage, a RSIC test sequence may often give rise to a better test quality than that obtained with RMIC delay test sequences. This is the most important result we wanted to demonstrate in this paper. Note that this conclusion drawn from a study on circuit s382 is also valid for most of the ISCAS'89 benchmark circuits.

Application to Other Circuits and Other Fault Models

An important comment on the validation results given in this paper is that the ISCAS'85 circuits family has not been used in our experiments. This is because the number of path faults in these circuits is too huge, and TestGen [START_REF] Testgen | version Tg3.0.2 User Guide[END_REF] is unable to generate the corresponding delay fault dictionaries for fault simulation. However, results similar to those reported in the previous section for circuit s382 were obtained for various circuits of the ISCAS'89 benchmark set. For example, results obtained for circuits s1238 are given in Fig. 4.

Results obtained for other ISCAS'89 circuits are listed in Table 4 (the test length considered are 100 Ld for every circuit).

A RSIC sequence has the important property that it is made of adjacent vectors. Hence, it may be efficient for testing delay faults (the aim of this paper) or for testing stuck-open faults that also require pairs of vectors. In addition, it may be interesting to study the efficiency of RSIC sequences for other fault models.

For stuck-at faults, for example, it is expected that the test length would be longer than for a RMIC sequence since a RMIC sequence samples more quickly all the input state space. However, all this space can be covered when the RSIC sequence length increases. This observation is confirmed by the example given in Fig. 5.

Concluding Remarks and Future Work

The efficiency of RSIC test sequences for delay fault testing was analyzed in this paper. The performance measurement takes into account both robust tests and non-robust tests. The main conclusion drawn from this study is that, even with a lower non-robust delay fault coverage, a RSIC test sequence may often give rise to a better test quality than that obtained with RMIC delay test sequences. This conclusion is a fundamental result that can be further used by test engineers when they will have to provide test sequences for delay fault testing. Several extension to this study have been foreseen for the near future.

A first one concerns with the development of a complete hardware generation structure providing random SIC sequences [START_REF] Abraham | Testing of MOS VLSI Circuits[END_REF]. This generation will require that: 1) a random number is hardware transformed into a random bit; 2) the bit changed at time l is completely independent from the bit changed at time l -1; 3) the period of the RSIC sequence is large enough. Solutions to points 1 and 2 are explicitly given in [START_REF] David | Random Testing of Digital Circuits: Theory and Applications[END_REF], and point 3 is to be developed. In addition, the area overhead has to be estimated and compared with other methods.

Another extension of this study concerns with the evaluation of SIC test sequences in testing faults other than delay faults, i.e. stuck-at faults, bridging faults, etc. A preliminary study has been carried out on some ISCAS circuits. Results and analysis can be found in [START_REF] Girard | Comparison Between Random and Pseudo-Random Generation for BIST of Delay, Stuck-at and Bridging Faults[END_REF]. Some discussions on the universal nature of SIC sequences for the testing purpose will be tackled in a future paper.

Fig. 1 .

 1 Fig. 1. Robust and non-robust test conditions for path delay faults.

Fig. 2 .

 2 Fig. 2. Efficiency comparison of RMIC and RSIC sequences.

Fig. 3 .

 3 Fig. 3. Comparison of weighted efficiencies of RMIC and RSIC sequences for various values of S NR : (a) 10%, (b) 50%, (c) 90%.

Fig. 4 .

 4 Fig. 4. Comparison of RSIC and RMIC sequence efficiencies for circuit s1238.

Fig. 5 .

 5 Fig. 5. Comparison between stuck-at fault coverage of RSIC and RMIC-circuit s382.

Table 1 .

 1 Results obtained on the ISCAS'89 circuits.

	Circuits	n	Ls	Fc SA	Ld	Fc PR	Fc NR
	s298	17	38 100	688	76.19	78.79
	s382	24	36 100	1398	88	91.75
	s386	13	86 100	736	100	100
	s420	35	55 100	1322	100	100
	s510	25	67 100	1366	100	100
	s526	24	68	99.82	1396	86.34	87.8
	s641	54	46 100	2778	66.16	71.02
	s713	54	90	93.46	4418	22.54	37.13
	s1238	32 196	94.80	4992	62.17	62.93
	s1494	14 158	99.19	3518	98.83	98.89
	s3330	172	238	100	15862	89.87	90.94
	s5378	214	143	99.13	29730	74.02	82.2

Table 2 .

 2 PR and NR efficiencies of given RMIC and RSIC sequences.

			RMIC sequence	RSIC sequence
		Number
	s382	of vectors	Eff PR	Eff NR	Eff PR	Eff NR
		10	12.78	16.49	1.70	1.63
		100	25.57	52.04	14.35	15.40
	Ld/4	349	30.97	73.97	24.57	25.61
	Ld/2	699	33.09	86.78	43.47	46.05
	Ld	1398	33.24	91.01	52.84	57.22
	2Ld	2796	33.52	96.86	61.36	66.49
	5Ld	6990	33.52	99.31	76.56	84.33
	10Ld	13980	33.52	99.86	81.53	90.32
	20Ld	27960	33.52	99.86	86.65	96.19
	50Ld	69900	33.52	99.86	89.06	98.64
	100Ld	139800	33.52	99.86	90.19	99.73
		1000000	33.52	99.86	90.48	100.00

Table 3 .

 3 Weighted efficiency for various assumptions about the success rate S NR of NR tests. Eff NR -Eff PR S NR 10% S NR 50% S NR 90% Eff NR -Eff PR S NR 10% S NR 50% S NR 90%

				Weighted efficiencies for circuit s382
				RMIC sequence			RSIC sequence
		Number of vectors 10	3.71	12.60	13.98	15.36	-0.07	1.63	1.60	1.58
		50	16.02	21.92	27.88	33.83	0.73	8.51	8.78	9.06
		100	26.48	26.98	36.82	46.66	1.05	13.86	14.25	14.64
	Ld/4	349	43.01	33.69	49.67	65.65	1.04	23.66	24.05	24.43
	Ld/2	699	53.69	36.72	56.67	76.62	2.58	41.92	42.88	43.84
	Ld	1398	57.77	37.24	58.71	80.18	4.38	51.08	52.71	54.34
	2Ld	2796	63.34	38.03	61.57	85.11	5.12	59.32	61.23	63.13
	5Ld	6990	65.79	38.26	62.71	87.16	7.77	74.14	77.03	79.91
	10Ld	13980	66.34	38.31	62.96	87.61	8.79	79.01	82.27	85.54
	20Ld	27960	66.34	38.31	62.96	87.61	9.54	83.98	87.53	91.07
	50Ld	69900	66.34	38.31	62.96	87.61	9.58	86.30	89.86	93.42
	100Ld	139800	66.34	38.31	62.96	87.61	9.53	87.38	90.92	94.47
		1000000	66.34	38.31	62.96	87.61	9.52	87.65	91.19	94.73
	be actually detected during testing of the CUT. On the
	other hand, at most 91.75% of the NR testable faults
	will be actually detected. As a matter of fact, there are
	91.75 -31.88 = 59.87% of the NR testable faults for
	which no PR test has been found but for which a NR
	test has been found. A question now arises: what is the
	fraction of faults that will be actually detected among
	the faults that are only NR testable? This depends on
	the number and the nature of the distributed defects that
	affect the circuit and that can invalidate the non-robust
	tests. Let us call this fraction the success rate of NR
	tests, and denote by S NR .

Table 4 .

 4 Comparison between pseudo-robust fault coverage of RSIC and RMIC testing.

		ATPG	Number of patterns	RSIC	RMIC
	Circuit	FC %	⇒	Eff %	Eff %
	s298	76.19	68800	91.48	36.65
	s382	88.00	139800	91.19	33.38
	s386	100	73600	98.79	40.82
	s420	100	132200	66.86	37.24
	s510	100	136600	91.11	31.78
	s526	86.34	139600	90.40	30.09
	s641	66.16	277800	76.01	18.12
	s713	22.54	441800	83.72	20.76
	s1238	62.17	499200	89.11	29.35
	s1494	98.83	351800	94.51	30.83
	s3330	89.87	158620	62.28	23.47
	s5378	74.02	297300	64.77	24.51

He got his PhD degree in electrical engineering in 1983 on symbolic layout for IC design. Since 1984, he was been interested in the testing domain. He obtained the "doctorat d'état" degree in 1987 for his works on switch level automatic test pattern generation. He is presently interested in delay fault testing, design for testability and power consumption optimization. He has authored and co-authored numerous papers on these fields, and has supervised several PhD dissertations. He has also participated to several European project (Microelectronic regulation, Esprit, Medea).