J A C Martins

S Barbarin

M Raous

A Pinto Da Costa

Dynamic stability of finite dimensional linearly elastic systems with unilateral contact and Coulomb friction

Necessary and sufficient conditions are established for the occurrence of dynamic instabilities in finite dimensional linearly elastic systems in unilateral frictional contact with a rigid flat surface. These conditions apply in particular to the systems that result from the finite element discretization of linearly elastic bodies. From the numerical point of view, these conditions lead to studying eigenproblems relative to a non-symmetric (tangent) stiffness matrix that incorporates the effect of the current state of the contact candidate particles. Illustrative small-sized examples are presented together with an application to the case of an experimentally tested block of polyurethane, where friction induced instability phenomena were observed.

Introduction

This paper adresses the subject of the dynamic stability of equilibrium states and some quasistatic evolutions of finite dimensional linearly elastic systems with plane motion having frictional contact with a fixed flat surface.

This topic is important for at least two reasons. First, because the instability of frictional contact systems may give rise to friction-induced oscillations, which in turn may be responsible for undesired noise and/or severe wear of the contacting surfaces. Secondly, because frictional systems provide a prototype for the non-associative elastic-plastic behavior encountered for instance in soils, rocks, concrete, ceramics, etc.

The analysis in this paper is restricted to finite dimensional systems, because at present, enormous mathematical difficulties still remain to be solved in what concerns the proper formulation 1 of dynamic unilateral contact problems involving continuum bodies with or without friction.

Even in what concerns finite dimensional dynamic problems with unilateral contact and friction, rigorous mathematical formulations and complete existence studies are available only in some particular situations [START_REF] Jean | A system of rigid bodies with dry friction[END_REF][START_REF] Manuel | Differential Inclusions in NonSmooth Mechanical Problems : Shocks and Dry Friction[END_REF][START_REF] Stewart | Existence of solutions to rigid body dynamics and the Painlevé paradoxes[END_REF]. The unilaterality of the contact and the non-symmetry (non-associativity, non-self-adjointness) induced by the Coulomb friction law are responsible for these mathematical difficulties. Of course, some of the results obtained in the present paper, in particular their variational inequality statements, can be formally(!) re-stated for the corresponding continuum problems.

Since the classical unilateral contact and Coulomb friction laws introduce multivalued operators in their mathematical formulation, frictional contact problems are highly non-linear. By restricting the scope of the present paper to linearly elastic systems, we aim at taking advantage, as much as possible, of the consequent linear structure of the problem outside the surfaces of discontinuity induced in the phase space by the above mentioned severe non-linearities (multivaluedness). Under appropriate conditions or restrictions, some of the results obtained here can of course be extended to other finite dimensional systems having a smooth non-linear elastic behavior [START_REF] Martins | Stability of finite dimensional systems with unilateral contact and friction : nonlinear elastic behavior and obstacle curvature[END_REF].

Due to its practical importance, the subject of the stability of frictional contact systems has been dealt with in numerous studies ; a recent survey on these works has been published by Ibrahim [START_REF] Ibrahim | Friction-induced vibration, chatter, squeal and chaos, Part I : Mechanics of contact and friction, Part II : Dynamics and modelling[END_REF]. However no general theoretical framework for these analyses, which may be systematically applicable to large systems (rigid multibody systems, finite element systems) has yet been developed. In fact, most of the early analyses of friction-induced instabilities and oscillations were limited to single degree-of-freedom (linearly elastic) systems, which left practically no other alternative for interpreting the observed phenomena than assuming, in some way or another, that the friction coefficient decreases with the sliding speed. By considering more complex (multi-degree-of-freedom) models of the actual systems, it became clear that the geometrical and dynamic properties of a system, the coupling between various degrees of freedom (tangent, normal, rotational, etc...) and the intrinsic non-symmetry of Coulomb friction law might give rise to various instability and oscillation phenomena [START_REF] Spurr | A theory of brake squeal[END_REF][START_REF] Jarvis | Vibrations induced by dry friction[END_REF] (for additional references see [START_REF] Ibrahim | Friction-induced vibration, chatter, squeal and chaos, Part I : Mechanics of contact and friction, Part II : Dynamics and modelling[END_REF]). These points were stressed in [START_REF] Oden | Models and computational methods for dynamic friction phenomena[END_REF] and [START_REF] Martins | A study of static and kinetic friction[END_REF], where, for metallic interfaces, the classical unilateral contact law was replaced by a physically motivated normal compliance law in the study of some friction-induced instability and oscillation phenomena. In what concerns the classical unilateral contact and Coulomb friction laws, the first note worthy result was obtained by Chateau and Nguyen [START_REF] Chateau | Buckling of elastic structures in unilateral contact with or without friction[END_REF], who established a sufficient condition for the absence of divergence type instabilities affecting the equilibrium of a continuous elastic body. An important contribution to the stability analysis of finite dimensional non-linear elastic systems with unilateral contact but without friction was made by Klarbring [START_REF] Klarbring | On discrete and discretized non-linear elastic structures in unilateral contact (stability, uniqueness and variational principles)[END_REF]. Mròz and Plaut have studied in [START_REF] Mróz | On the stability and post-critical behavior of elastic structures with dry friction[END_REF] the stability of elastic structures with dry friction, but the results apply only to systems such that the normal contact forces are a priori known, which eliminates the non-associative character of the most general frictional contact problems.

Of course the studies on the stability of frictional systems have some parallels with those on elastic-plastic bodies and, generally speaking, with those on rate-independent dissipative systems. Recent surveys on these topics can be found in [START_REF] Nguyen | Bifurcation and stability in dissipative media (plasticity, friction, fracture)[END_REF] and in the lectures in [START_REF]Bifurcation and Stability of Dissipative Systems[END_REF]. In particular, we wish to mention the early contribution of Mandel [START_REF] Mandel | Conditions de stabilité et postulat de Drucker[END_REF], because, in the case of non-associative elastic-plastic solids, it stresses the distinction between necessary and sufficient conditions for dynamic instability (by divergence) and presents simple frictional contact examples. It is also worth mentioning the contribution of Petryk [START_REF] Petryk | The energy criteria of instability in time-independent inelastic solids[END_REF][START_REF] Petryk | General theory of bifurcation and stability in time-independent plasticity, Bifurcation and Stability of Dissipative Systems[END_REF], who, in the case of associative elastic-plastic solids and in a finite dimensional context, proves a sufficient condition for divergence instability : the direct paths emanating from the equilibrium state used in his proof inspired the straight (direct) paths used with similar purposes in the present context of (non-associative) finite dimensional frictional contact systems. This paper is organized as follows. In Section 2, the notations used in the paper are introduced and the formulation of some finite dimensional frictional contact problems is recalled ; particular attention is paid to the characterization of the sets to which the various kinematic and static variables (displacements, reactions, and their rates) belong. In Section 3, a study is presented on the existence of dynamic solutions in the neighborhood of a given equilibrium state which might tend to diverge from that state : dynamic solutions beginning with a velocity discontinuity, dynamic solutions starting with acceleration and reaction discontinuities, and smooth dynamic solutions with initial conditions arbitrarily close to the equilibrium state are considered in turn. In Section 4, a study is presented on the existence of growing (oscillatory or non-oscillatory) dynamic solutions in the neighborhood of straight portions of quasistatic evolutions of the system. Some illustrative small-sized examples are presented at the end of Sections 3 and 4. Numerical results obtained with the finite element method are presented in Section 5, in the case of a block of polyurethane which was tested experimentally by Zeghloul and Villechaise [START_REF] Barbarin | Ondes de contraintes dans un contact glissant, modélisation et analyse experimentale[END_REF] : sudden changes in the stress field in the course of the tangential loading process were detected by those authors using photoelasticity techniques ; the possibility of correlating these experimental results with friction-induced instability phenomena is discussed.

Notation, formulation and preliminary results

2.1. Dynamic frictional contact problem. Admissible displacements, velocities, accelerations and reactions.

We consider a finite dimensional linearly elastic system with plane motion that may establish frictional contact with a fixed flat surface. A typical situation in mind is the one resulting from a finite element discretization of an elastic body, the possible contact being studied only at a set of nodes (particles) on the boundary of the discretized body (see Fig. 1).

Let t ≥ 0 be the time variable and let u i (t) (i ∈ S ⊂ IN ) be the generalized displacements of the system at time t ; the corresponding column vector of the N generalized displacements is denoted by u u u(t) ∈ IR N . The set of labels of the generalized displacements is the union of three disjoint subsets S D , S F and S C : for i ∈ S D , u i (t) is a prescribed Displacement component ; for i ∈ S F , u i (t) is an unknown (Free) displacement function which is not subjected to any kinematic constraint ; and for i ∈ S C , u i (t) is an unknown displacement component of a particle of the system that may establish Contact with the fixed surface. The latter displacement components may be either Tangential or Normal to the fixed surface, so that S C = S T ∪ S N and, for each contact candidate particle p (p ∈ P C ⊂ IN ) those displacement components are also denoted by u T p (t) and u N p (t), respectively. We denote by Φ(t) ∈ IR N the given vector of applied forces at time t and we denote by M and K the N x N symmetric, positive definite mass and positive semi-definite stiffness matrices, respectively. In order to homogeneize the prescribed displacements, let Ψ(t) ∈ IR N such that

Ψ i (t) =    u i (t), for all i ∈ S D 0, for all i ∈ S F ∪ S C
and let

u(t) def = u u u(t) -Ψ(t) , f (t) def = Φ(t) -KΨ(t) -M Ψ(t) ,
where ( ) denotes the time derivative d()/dt. The dynamic equations of motion of the present system are then of the form

M ü(t) + Ku(t) = f (t) + r(t) , (1)
, where r(t) ∈ IR N denotes the unknown vector of the reaction forces at time t. In view of the presence or absence of kinematic constraints in S D and S F , respectively, we have

u C P p Np u F S i , u i D S i i u Tp
u(t) ∈ V u , r(t) ∈ V r , (2)
where V u and V r denote the spaces of admissible displacement and reaction vectors,

V u def = u ∈ IR N : u i = 0, for all i ∈ S D , V r def = r ∈ IR N : r i = 0, for all i ∈ S F . (3)
In addition, the classical unilateral contact conditions,

u N p (t) ≤ 0 , r N p (t) ≤ 0 , u N p (t) r N p (t) = 0 , for all p ∈ P C , (4)
and the friction law of Coulomb,

r T p (t) ∈ µ r N p (t) σ( uT p (t)) , for all p ∈ P C , (5)
must hold for all the contact candidate particles; µ ≥ 0 is the coefficient of friction and σ denotes the multivalued application such that, for each x ∈ IR,

σ(x) def =    x/|x| , if x = 0, [-1, 1], if x = 0.
Initial displacement and velocity conditions have also to be satisfied :

u(0) = u 0 , u(0) = v 0 ; (6)
some conditions on the data u 0 and v 0 will be specified later.

REMARK 2.1. The unilateral contact conditions (4) and the friction law of Coulomb [START_REF] Clough | Dynamics of structures[END_REF] have forms that will be frequently encountered in the rest of this paper when dealing with the kinematic and static variables in (4) and ( 5) and with their time rates of change. Conditions of the type in (4) or ( 5) lead to constraints on the global vectors of the relevant kinematic and static variables which, in the present paper, will be established in the form of appropriate inclusion or variational statements : some (primal) unknown kinematic variable must belong to some admissible set (which depends at most on other already known kinematic or static variables) while the corresponding (dual) unknown static variable must belong to some admissible set which depends on the unknown value of the (primal) kinematic variable ; in these circumstances, the set of the admissible values of the unknown static variable will also be characterized by a variational inequality statement.

In what concerns complementarity conditions of the type in (4),

x ≤ 0, y ≤ 0, xy = 0, [START_REF] Cocu | Formulation and approximation of quasistatic frictional contact[END_REF] these are obviously equivalent to inclusion statements of the form

x ∈ K x def = {x : x ≤ 0} , y ∈ K y (x) def = {y : y ≤ 0 and yx = 0} , ( 8
)
and it is well known that ( 7) or ( 8) are equivalent to

x ∈ K x and y(x ′ -x) ≥ 0 , ∀x ′ ∈ K x . (9)
On the other hand, the inclusion statement in ( 5) is of the general form

y ∈ zσ(x) (10)
where z is independently known to satisfy z ≤ 0. In these circumstances ( 10) is known to be equivalent to

y(x ′ -x) ≥ z(|x ′ | -|x|) , ∀x ∈ IR . ( 11
)
Of course other dual statements could be given for ( 8)- [START_REF] Jarvis | Vibrations induced by dry friction[END_REF], and also some mixed-type statements ; note also that, by appropriate changes of variables, the statement (10) could be transformed into complementarity statements of the type in [START_REF] Cocu | Formulation and approximation of quasistatic frictional contact[END_REF]. The interested reader may find the development of these alternative formulations for some related problems in e.g. [START_REF] Klarbring | Contact problems with friction by linear complementarity, Unilateral problems in structural mechanics[END_REF][START_REF] Klarbring | Contact friction, discrete mechanical structures and mathematical programming, Contact problems : theory, methods and applications[END_REF][START_REF] Pfeiffer | Multibody dynamics with unilateral contacts[END_REF][START_REF] Trinkle | On dynamic multi-rigid-body problems with Coulomb friction[END_REF][START_REF] Pang | Complementarity formulation and existence of solutions of dynamic multi-rigid-body contact problem with Coulomb friction[END_REF].

2

The interpretation of the governing equations and conditions (1)-( 6) needs some care since in general the displacements u will not be twice continuously differentiable functions of time. For instance, at each transition of a contact particle from sliding to sticking, the corresponding tangential acceleration and reaction may be discontinuous. In other circumstances the velocities may be discontinuous : when a collision occurs between a contact candidate particle and the obstacle surface, or when a frictional catastrophe occurs (a velocity jump not associated to a collision [START_REF] Moreau | Unilateral contact and dry friction in finite freedom dynamics[END_REF]).

The mathematical framework suggested by Moreau and Jean [START_REF] Moreau | Unilateral contact and dry friction in finite freedom dynamics[END_REF][START_REF] Jean | Dynamics of elastic or rigid bodies with frictional contact : numerical methods, Rencontres Scientifiques du Cinquantenaire du L[END_REF][START_REF] Jean | Unilaterality and dry friction in the dynamics of rigid body collections[END_REF] to study these dynamic contact problems with friction follows from the assumption that u is an absolutely continuous function with a time derivative u with locally bounded variation. This implies that a left velocity u-(t) and a right velocity u+ (t) exist at every t ≥ 0. In these circumstances, for all times t ≥ 0, the displacement u(t) must belong to the set of admissible displacements K u :

u(t) ∈ K u def = {u ∈ V u : u N p ≤ 0 , if p ∈ P C } . (12)
In addition, for each t ≥ 0, the right and left velocities satisfy

u+ (t) ∈ K u (u(t)) , u-(t) ∈ -K u (u(t)) , (13)
where, for each u ∈ K u , the decomposition of the set P C of contact candidate particles is introduced

P C = P f (u) ∪ P c (u) , P f (u) def = {p ∈ P C : u N p < 0} [particles currently not in contact (free)] P c (u) def = {p ∈ P C : u N p = 0} [particles currently in contact]
and

K u(u) def = {v ∈ V u : v N p ≤ 0 , for each p ∈ P c } , (14)
is the (displacement dependent) set of admissible right velocities, while -K u(u) is the corresponding set of admissible left velocities.

If right (or left) limits exist for the reaction forces at some time t, r + (t) (or r -(t)), it follows from ( 4) and ( 5) that

r + (t) ∈ K r u(t), u+ (t) r -(t) ∈ K r (u(t), u-(t)) , (15)
where, for

u ∈ K u and v ∈ K u(u) (or v ∈ -K u(u)), K r (u, v) def = {r ∈ V r :r N p = r T p = 0, for p ∈ P f ; r N p ≤ 0 , r N p v N p = 0 and r T p ∈ µ r N p σ(v T p ) , for p ∈ P c } (16)
is the (displacement and velocity dependent) set of admissible right (or left) reaction forces, and

K r (u) def = {r ∈ V r :r N p = r T p = 0, for p ∈ P f ;
r N p ≤ 0, and

|r T p | ≤ -µ r N p , for p ∈ P c } (17)
is the (displacement dependent) set of admissible reaction forces.

Note that K r (u, v) ⊂ K r (u) and that, for simplicity of notation, the dependence of the sets P c and P f on the displacements u is omitted in ( 14), ( 16) and [START_REF] Klarbring | Examples of non-uniqueness and non-existence of solutions to quasistatic contact problems with friction[END_REF] ; similar simplifications will be made in the rest of the paper, whenever confusion is not likely to arise.

At an instant of a velocity discontinuity, an impulsive interpretation has to be given to the dynamic equation (1) :

M u+ (t) -u-(t) = I(t) , (18)
where I(t) is the impulse or shock percussion of the reaction forces; it represents the "integral" of the reaction forces over the "very short interval of time" during which the velocity jump occurs. Indeed this interpretation is consistent with a mathematical reinterpretation of the differential equation ( 1) as a measure differential equation (see [START_REF] Moreau | Unilateral contact and dry friction in finite freedom dynamics[END_REF][START_REF] Manuel | Differential Inclusions in NonSmooth Mechanical Problems : Shocks and Dry Friction[END_REF] for details). In what concerns the restrictions posed by the frictional contact conditions on the shock percussion I, Moreau and Jean [START_REF] Moreau | Unilateral contact and dry friction in finite freedom dynamics[END_REF][START_REF] Jean | Dynamics of elastic or rigid bodies with frictional contact : numerical methods, Rencontres Scientifiques du Cinquantenaire du L[END_REF][START_REF] Jean | Unilaterality and dry friction in the dynamics of rigid body collections[END_REF] propose that I(t) must satisfy a condition of the type in [START_REF] Klarbring | Contact problems with friction by linear complementarity, Unilateral problems in structural mechanics[END_REF] but involving some weighted average velocity during the velocity jump :

I(t) ∈ K r (u(t), ua (t)) , (19)
where

ua (t) = 1 -δ 2 u-(t) + 1 + δ 2 u+ (t) (20)
is the average velocity, and

δ ∈ [0, 1] (21)
is the dissipation index. Justifications and interpretations for these definitions can be found in [START_REF] Moreau | Unilateral contact and dry friction in finite freedom dynamics[END_REF][START_REF] Jean | Dynamics of elastic or rigid bodies with frictional contact : numerical methods, Rencontres Scientifiques du Cinquantenaire du L[END_REF][START_REF] Jean | Unilaterality and dry friction in the dynamics of rigid body collections[END_REF]. Here we only mention the following.

REMARK 2.2. The above setting for the dynamic contact problem with friction requires the initial data to be given in appropriate sets and the initial conditions (6) to be satisfied in an appropriate sense :

u(0) = u 0 ∈ K u , u-(0) = v 0 ∈ -K u(u 0 ) . 2 (22)
REMARK 2.3. Because of the definition [START_REF] Klarbring | On discrete and discretized non-linear elastic structures in unilateral contact (stability, uniqueness and variational principles)[END_REF] for the set K r (u(t), ua (t)), only the definition of the average velocity components for the particles currently in contact is relevant : the particles p ∈ P c (u(t)). In addition, a non-null shock percussion at some particle p ∈ P C may occur only if that particle is currently in contact and has a null average normal velocity, i.e.

for all p ∈ P C :

I p (t) = 0 ⇒ u N p (t) = ua N p (t) = 0 . 2 LEMMA 2.4.
The change in the kinetic energy of the system (T ) in the course of a velocity jump is equal to

T ( u+ (t)) -T ( u-(t)) =- 1 2 δM ( u+ (t) -u-(t)).( u+ (t) -u-(t)) + p∈P C µI N p (t)| ua T p (t)|, (23)
which is a non-positive quantity.

PROOF. For any right or left velocity v, the kinetic energy is given by the positive definite

quadratic form T (v) = 1 2
Mv.v . The result follows by doing the inner product of (18) with ua (t) and by using the definitions [START_REF] Klarbring | Contact friction, discrete mechanical structures and mathematical programming, Contact problems : theory, methods and applications[END_REF] and [START_REF] Klarbring | On discrete and discretized non-linear elastic structures in unilateral contact (stability, uniqueness and variational principles)[END_REF] of ua (t) and K r (u(t), ua (t)). See [START_REF] Moreau | Unilateral contact and dry friction in finite freedom dynamics[END_REF][START_REF] Jean | Dynamics of elastic or rigid bodies with frictional contact : numerical methods, Rencontres Scientifiques du Cinquantenaire du L[END_REF][START_REF] Jean | Unilaterality and dry friction in the dynamics of rigid body collections[END_REF] for details. 2 REMARK 2.5. In a more recent work, Moreau [START_REF] Moreau | Some numerical methods in multibody dynamics : application to granular materials[END_REF] proposes a more general form for the frictional impact conditions in [START_REF] Klarbring | Mathematical programming in contact problems[END_REF], ( 20) and ( 21) ; other approaches have also been proposed in the literature (see, e.g., Pfeiffer and Glocker [START_REF] Pfeiffer | Multibody dynamics with unilateral contacts[END_REF] and Frémond [START_REF] Frémond | Rigid bodies collisions[END_REF]). However, what is relevant for the present paper (see Section 3.1) is that, whatever the adopted model is, the fundamental dissipativity property in Lemma 2.4 holds : the change in kinetic energy in the course of a velocity jump must be non-positive.

Let us assume now that the right acceleration and the right reaction are well defined at some time t : ü+ (t) and r + (t). In order to characterize the sets to which these quantities belong we consider the following decomposition of the set P c of the particles currently in contact. For each admissible displacement u ∈ K u and each admissible right velocity v we let

P c (u) = P cf (u, v) ∪ P 0 (u, v) ∪ P v (u, v) ,
where

P cf (u, v) def = {p ∈ P c (u) : v N p < 0}
[particles in contact with negative normal velocity (free in the near future)] P 0 (u, v) and we define the (displacement and velocity dependent) set of admissible right accelerations,

K ü(u, v) def = {a ∈ V u : a N p ≤ 0 , for each p ∈ P 0 ∪ P v } (24)
and [having also a ∈ K ü(u, v)] we define the (displacement, velocity and acceleration dependent) set of admissible right reactions, K r (u, v, a) def = {r ∈ V r :r N p ≤ 0 and r N p a N p = 0, for p ∈ P 0 ∪ P v ; r T p ∈ µr N p σ(a T p ), for p ∈ P 0 ;

r T p = µr N p σ(v T p ), for p ∈ P v ; r N p = r T p = 0, for p ∈ P f ∪ P cf } ⊂ K r (u, v) ⊂ K r (u) . (25)
The conditions on the right acceleration and reaction vectors are then written in the form

ü+ (t) ∈ K ü(u(t), u(t)) , (26) r
+ (t) ∈ K r (u(t), u(t), ü+ (t)) . ( 27
)
Note that the conditions on a N p and r N p in [START_REF] Maier | Incremental plastic analysis in the presence of large displacements and physical instabilizing effects[END_REF] and [START_REF] Mandel | Conditions de stabilité et postulat de Drucker[END_REF] at the particles that are currently in contact with vanishing normal right velocity (u N p = v N p = 0, p ∈ P 0 ∪ P v ) are complementarity conditions of the type [START_REF] Cocu | Formulation and approximation of quasistatic frictional contact[END_REF] (see [START_REF] Maier | Incremental plastic analysis in the presence of large displacements and physical instabilizing effects[END_REF] and the first line in [START_REF] Mandel | Conditions de stabilité et postulat de Drucker[END_REF]) ; on the other hand the conditions on a T p and r T p at the particles that are currently in contact with vanishing normal and tangential right velocities (u N p = v N p = v T p = 0, p ∈ P 0 , see the second line in [START_REF] Mandel | Conditions de stabilité et postulat de Drucker[END_REF]) have a form of the type [START_REF] Ibrahim | Friction-induced vibration, chatter, squeal and chaos, Part I : Mechanics of contact and friction, Part II : Dynamics and modelling[END_REF] ; finally note that the remaining conditions on the right reactions (see the two last lines in [START_REF] Mandel | Conditions de stabilité et postulat de Drucker[END_REF]) are simple equality constraints. As a consequence of these observations and of the observations in Remark 2.1, the set of the admissible right reactions can be characterized in the following manner.

LEMMA 2.6. Let U ∈ K u , V ∈ K u(U ) and A ∈ K ü(U , V ). Then R ∈ K r (U , V , A) if and only if R.(a -A) ≥ p∈Pv µR N p σ(V T p )(a T p -A T p ) + p∈P 0 µR N p (|a T p | -|A T p |), ∀a ∈ K ü(U , V ) . 2 (28)
A detailed proof of this result can be found in [START_REF] Martins | Stability of finite dimensional systems with unilateral contact and friction : flat obstacle and linear elastic behavior[END_REF]. Related discussions on the conditions satisfied by the right accelerations and reactions in finite dimensional frictional contact problems (assuming that those quantities are well defined) can be found in [START_REF] Lötstedt | Coulomb friction in two-dimensional rigid body system[END_REF][START_REF] Pang | Complementarity formulation and existence of solutions of dynamic multi-rigid-body contact problem with Coulomb friction[END_REF][START_REF] Trinkle | On dynamic multi-rigid-body problems with Coulomb friction[END_REF].

Note that in the most common situations at which no velocity discontinuity exists, the normal velocity uNp (t) is null at each particle p currently in contact [p ∈ P c (u(t))] :

u+ (t) = u-(t) = u(t)∈K u(u(t)) ∩ (-K u(u(t))) = {v ∈ V u : v N p = 0 , for each p ∈ P c (u(t))} . (29)
Consequently, P cf (u(t), u(t)) = ∅ and P c (u(t)) = P 0 (u(t), u(t)) ∪ P v (u(t), u(t)).

2.2. Static and quasistatic frictional contact problems. Admissible first order displacement and reaction rates.

In the following sections we shall be concerned with the stability analysis of some particular types of solutions to the dynamic equations and conditions (1)- [START_REF] Cocu | Existence d'une solution du problème quasi-statique de contact unilatéral avec frottement non local[END_REF]. Our first concern will be the static equilibrium positions, i.e. the dynamic solutions with null velocity and acceleration. For given forces f 0 independent of time, the static equilibrium solutions u 0 satisfy thus

u 0 ∈ K u , r 0 def = Ku 0 -f 0 ∈ K r (u 0 ) . ( 30
)
We shall also be concerned with particular classes of quasistatic evolutions of the system. For given time dependent applied forces f and given initial displacements u 0 such that [cf. [START_REF] Jean | A system of rigid bodies with dry friction[END_REF] and [START_REF] Klarbring | Examples of non-uniqueness and non-existence of solutions to quasistatic contact problems with friction[END_REF]]

u 0 ∈ K u , Ku 0 -f (0) ∈ K r (u 0 ) , (31)
a quasistatic evolution is characterized by absolutely continuous functions u 0 and r 0 such that the balance equations

Ku 0 (t) = f (t) + r 0 (t) (32)
and the unilateral contact conditions (4) hold for all times t ≥ 0, the friction law (5) holds for a.e. t ≥ 0 and the initial condition

u 0 (0) = u 0 (33)
is satisfied at t = 0. It is now well-known that unless the friction coefficient is sufficiently small [START_REF] Licht | Remarks on a numerical method for unilateral contact including friction[END_REF], simple examples of non-existence or non-uniqueness of solution in the above mathematical framework can be found [START_REF] Klarbring | Examples of non-uniqueness and non-existence of solutions to quasistatic contact problems with friction[END_REF].

If u 0 and r 0 are right-differentiable at some τ ≥ 0, the right-derivatives u0+ (τ ) and ṙ0+ (τ ) must belong to appropriate sets. In order to characterize such sets we consider the following decomposition of the set of the particles that are currently in contact. For each u ∈ K u and each r ∈ K r (u), we let

P c (u) = P z (u, r) ∪ P d (u, r) ∪ P s (u, r) ,
where the above disjoint subsets are

P z (u, r) def = {p ∈ P c (u) : r N p = r T p = 0}
[particles in contact with zero reaction] Then for u ∈ K u , r ∈ K r (u) and v ∈ K u(u, r), we define

K u(u, r) def = {v ∈ V u :v N p ≤ 0 , for p ∈ P z ;
v N p = 0 , for p ∈ P d ∪ P s ;

v T p σ(r T p ) ≤ 0 , for p ∈ P s ;

v T p = 0 , for p ∈ P d } (34)
as the (displacement and reaction dependent) set of admissible right velocities, and

K ṙ (u, r, v) def = {w ∈ V r : w N p ≤ 0 , w N p v N p = 0, for p ∈ P z ;
w N p = 0, for p ∈ P f ;

w T p σ(r T p ) + µ w N p ≤ 0 and [w T p σ(r T p ) + µ w N p ] [v T p σ(r T p )] = 0, for p ∈ P s ; w T p ∈ µ w N p σ(v T p ), for p ∈ P z ; w T p = 0, for p ∈ P f } (35)
as the (displacement, reaction and velocity dependent) set of admissible right reaction rates, so that

u0+ (τ ) ∈ K u(u 0 (τ ), r 0 (τ )) , ṙ0+ (τ ) ∈ K ṙ (u 0 (τ ), r 0 (τ ), u0+ (τ )) . (36)
A visualization of the decomposition of the set of contact candidate particles P C into P f , P z , P d and P s is shown in Fig. 2, together with the corresponding admissible values of u N , r N and r T as well as their admissible (right) rates of change.

Note that the conditions on v N p and w N p in ( 34) and ( 35) at the particles that are currently in contact with vanishing reaction (u N p = r N p = r T p = 0, p ∈ P z , see the first lines in [START_REF]Bifurcation and Stability of Dissipative Systems[END_REF] and ( 35)) are complementarity conditions of the type [START_REF] Cocu | Formulation and approximation of quasistatic frictional contact[END_REF] ; note also that the conditions on x = v T p σ(r T p ) and y = w T p σ(r T p ) + µ w N p at the particles that are currently in contact with non-vanishing reaction on the boundary of the friction cone (u N p = 0, |r T p | = -µr N p > 0, p ∈ P s , see the third lines in [START_REF]Bifurcation and Stability of Dissipative Systems[END_REF] and [START_REF] Nguyen | Bifurcation and stability in dissipative media (plasticity, friction, fracture)[END_REF]) are also complementarity conditions of the type [START_REF] Cocu | Formulation and approximation of quasistatic frictional contact[END_REF] ; the condition on w N p , w T p and v T p in P z (see the fourth line in [START_REF] Nguyen | Bifurcation and stability in dissipative media (plasticity, friction, fracture)[END_REF]) is of the type [START_REF] Ibrahim | Friction-induced vibration, chatter, squeal and chaos, Part I : Mechanics of contact and friction, Part II : Dynamics and modelling[END_REF] ; and all the remaining conditions in [START_REF]Bifurcation and Stability of Dissipative Systems[END_REF] and [START_REF] Nguyen | Bifurcation and stability in dissipative media (plasticity, friction, fracture)[END_REF] are simple equality constraints. As a consequence of these observations and of the observations in Remark 2.1, the set (35) of the admissible right reaction rates admits the variational characterization given below. Essentially the same result is proved in [START_REF] Chateau | Buckling of elastic structures in unilateral contact with or without friction[END_REF] and [START_REF] Klarbring | Derivation and analysis of rate boundary-value problems of frictional contact[END_REF].

LEMMA 2.7. Let U ∈ K u , R ∈ K r (U ) and V ∈ K u(U , R). Then W ∈ K ṙ (U , R, V ) if and only if W .(v -V )≥ p∈Pc µW N p (|v T p | -|V T p |), ∀v ∈ K u(U , R) . 2 (37)
Note that the right-hand side of (37) can be written in more detail in the form :

p∈Pc µW N p (|v T p | -|V T p |) = p∈Pz µW N p (|v T p | -|V T p |) - p∈Ps µW N p σ(R T p )(v T p -V T p ) . ( 38
)
2.3. Quasistatic evolutions with straight portions of constant velocity. Admissible second order displacement and reaction rates.

It is important to observe that in general a quasistatic evolution is not a solution of the dynamic problem ( 1)-( 6) : if the quasistatic velocity u0 (t) is not constant, then (u 0 (t), r 0 (t)) satisfying [START_REF] Mouwakeh | Etude quantitative des phénomènes de glissements dans un contact sec à deux corps par comparaison avec la propagation d'une fissure interfaciale[END_REF] will not satisfy the dynamic equation [START_REF] Barbarin | Ondes de contraintes dans un contact glissant, modélisation et analyse experimentale[END_REF]. For this reason we shall be interested in analysing what happens in the neighborhood of some portions of quasistatic evolutions along which the quasistatic velocity is constant.

Let τ ≥ 0 and let the applied force rates be constant for t ≥ τ :

f (t) = f (τ ) + ḟ (τ )(t -τ ) , for t ≥ τ. (39)
Let us also assume that there exist quasistatic solution vectors of displacements, reactions, right velocities and right reaction rates at time τ :

u 0 (τ ) ∈ K u , r 0 (τ ) ∈ K r (u 0 (τ )), u0+ (τ ) ∈ K u(u 0 (τ ), r 0 (τ )) and ṙ0+ (τ ) ∈ K ṙ (u 0 (τ ), r 0 (τ ), u0+ (τ ))
, respectively. In other words these vectors satisfy at time τ the frictional contact conditions ( 4), [START_REF] Clough | Dynamics of structures[END_REF] or the corresponding rate form associated with the definitions ( 34), [START_REF] Nguyen | Bifurcation and stability in dissipative media (plasticity, friction, fracture)[END_REF], and also the balance equation [START_REF] Mouwakeh | Etude quantitative des phénomènes de glissements dans un contact sec à deux corps par comparaison avec la propagation d'une fissure interfaciale[END_REF] or its rate form

K u0+ (τ ) = ḟ (τ ) + ṙ0+ (τ ) . (40)
In these circumstances the following result holds.

LEMMA 2.8. There exists ∆ 0 τ > 0 such that

u 0 (t) = u 0 (τ ) + u0+ (τ )(t -τ ) , r 0 (t) = r 0 (τ ) + ṙ0+ (τ )(t -τ ) (41)
solve [START_REF] Mouwakeh | Etude quantitative des phénomènes de glissements dans un contact sec à deux corps par comparaison avec la propagation d'une fissure interfaciale[END_REF] and consequently (1) for all t ∈ [τ, τ + ∆ 0 τ [ .

PROOF. Satisfaction of the force balance equations [and consequently of the dynamic equations (1)

] is checked by direct substitution of ( 41) in [START_REF] Mouwakeh | Etude quantitative des phénomènes de glissements dans un contact sec à deux corps par comparaison avec la propagation d'une fissure interfaciale[END_REF] with u0+ (τ ) and ṙ0+ (τ ) being the above solution of the rate problem [START_REF] Oden | Models and computational methods for dynamic friction phenomena[END_REF], [START_REF] Pfeiffer | Multibody dynamics with unilateral contacts[END_REF]. Satisfaction of the frictional contact conditions results from the fact that the path defined by (u 0 (t),r 0 (t)) in the space of the displacements and reactions is straight and has the direction of the admissible displacement and reaction rates ( u0+ (τ ), ṙ0+ (τ )) at (u 0 (τ ),r 0 (τ )), while the boundaries of the regions in the displacement-reaction space where each decomposition of P C is valid are hyperplanes. The smallness of ∆ 0 τ > 0 is required to keep the solution pair (u 0 (t),r 0 (t)) for t ∈ τ, τ + ∆ 0 τ , within the admissible region selected by the rates ( u0+ (τ ), ṙ0+ (τ )). More precisely, ∆ 0 τ is sufficiently small that ∀t ∈ τ, τ + ∆ 0 τ :

u 0 N p (t) < 0, for all p ∈ P f (u 0 (τ )) ; r 0 N p (t) < 0 and |r 0 T p (t)| < -µ r 0 N p (t), for all p ∈ P d (u 0 (τ ), r 0 (τ )) ; r 0 N p (t) < 0, for all p ∈ P s (u 0 (τ ), r 0 (τ )) ; |r 0 T p (t)| < -µ r 0 N p (t), for all p ∈ P s (u 0 (τ ), r 0 (τ )) such that u0+ T p (τ ) = u0+ N p (τ ) = 0 , ṙ0+ T p (τ )σ(r 0 T p (τ )) + µ ṙ0+ N p (τ ) < 0 .
Note that no condition is imposed on the particles p ∈ P z (u 0 (τ ), r 0 (τ )) because for any straight non-trivial path of such a particle, it will remain indefinitely in one and only one of P f (u 0 (τ )), P d (u 0 (τ ), r 0 (τ )) or P s (u 0 (τ ), r 0 (τ )). 2

In Section 4 the possible existence of growing dynamic solutions in the neighborhood of these straight portions of quasistatic evolutions will be studied. Since the first order displacement and reaction rates are constant along these straight quasistatic paths, the increments of smooth dynamic solutions with respect to the quasistatic one will be assumed to belong to the sets of admissible right second order displacement and reaction rates. and their first order (right) rates for a contact candidate particle p ∈ P C ; nomenclature for the subsets of P z and P s .

In order to characterize these sets, we start by considering, for each u ∈ K u , r ∈ K r (u), v ∈ K u(u, r) and w ∈ K ṙ (u, r, v), the following decomposition of P z and P s :

P z =P zf ∪ P zz ∪ P zd ∪ P zs , P s =P sd ∪ P ss ,
where :

P zf (u, r, v) = {p ∈ P z (u, r) : v N p < 0}
[particles in P z with negative normal velocity]

P zz (u, r, v, w) = {p ∈ P z (u, r) : v N p = 0 , w T p = w N p = 0} [
P zd (u, r, v, w) = {p ∈ P z (u, r) : v T p = v N p = 0 , w N p < 0 , |w T p | < -µ w N p }
[particles in P z with non-vanishing reaction rate directed towards the interior of the friction cone]

P zs (u, r, v, w) = {p ∈ P z (u, r) : v N p = 0 , w N p < 0 , w T p = µ w N p σ(v T p )}
[particles in P z with non-vanishing reaction rate on the boundary of the friction cone]

P sd (u, r, v, w) = {p ∈ P s (u, r) : v T p = v N p = 0 , w T p σ(r T p ) + µ w N p < 0}
[particles in P s with non-vanishing reaction rate directed towards the interior of the friction cone]

P ss (u, r, v, w) = {p ∈ P s (u, r) : v N p = 0 , σ(v T p ) = -σ(r T p ) , w T p σ(r T p ) + µ w N p = 0}
[particles in P s with reaction rate on the boundary of the friction cone].

A visual interpretation of the nomenclature of these sets as the result of the first order rates of the normal displacements and of the normal and tangential reactions is presented in Fig. 2.

The decomposition of each of the sets P zz , P zs and P ss into their (disjoint) intersections with P 0 and with P v will be also important in the following : it will be important to know if the particles in those sets have a vanishing or a non-vanishing sliding velocity.

Then we denote by

K ü(u, r, v, w) ={a ∈ V u :a N p ≤ 0, for p ∈ P zz ;
a N p = 0, for p ∈ P d ∪ P s ∪ P zs ∪ P zd ;

a T p σ(r T p ) ≤ 0, for p ∈ P ss ∩ P 0 ; a T p σ(w T p ) ≤ 0, for p ∈ P zs ∩ P 0 ;

a T p = 0, for p ∈ P d ∪ P sd ∪ P zd }

the set of admissible right second order displacement rates (accelerations) ; for u, r, v and w in the same sets and a ∈ K ü(u, r, v, w), the set of the admissible right second order reaction rates is the set

K r (u, r, v, w, a) = {b ∈ V r : b N p ≤ 0, b N p a N p = 0, for p ∈ P zz ; b N p = 0, for p ∈ P f ∪ P zf ; b T p σ(r T p ) + µ b N p ≤ 0, [b T p σ(r T p ) + µ b N p ] [a T p σ(r T p )] = 0, for p ∈ P ss ∩ P 0 ; b T p σ(w T p ) + µ b N p ≤ 0, [b T p σ(w T p ) + µ b N p ] [a T p σ(w T p )] = 0, for p ∈ P zs ∩ P 0 ; b T p ∈ µ b N p σ(a T p ), for p ∈ P zz ∩ P 0 ; b T p = µ b N p σ(v T p ), for p ∈ P v ; b T p = 0, for p ∈ P f ∪ P zf } . ( 43
)
The conditions in ( 42) and ( 43) deserve a more detailed discussion. Concerning the (right) second order rates of the normal displacements and reactions, these satisfy complementarity conditions of the type (7) whenever the normal displacement and reaction, as well as their first order rates are null (p ∈ P zz , see the first lines in ( 42) and ( 43)) ; on the other hand (see the second line in [START_REF] Rajakumar | The Lanczos algorithm applied to unsymmetric generalized eigenvalue problem[END_REF]), the normal acceleration has to be null whenever the normal reaction or its first order rate is different from zero (p ∈ P d ∪ P s ∪ P zs ∪ P zd ) ; finally (see the second line in [START_REF] Ralston | A first course in numerical analysis[END_REF]) the second order rate of the normal reaction has to be null whenever the normal displacement or its first order rate indicate a near future absence of contact (p

∈ P f ∪ P zf ).
In what concerns the tangential variables, it follows from the friction law that, whenever the current sliding speed is null and the current or the (first order) near future reactions do not vanish and are on the boundary of the friction cone (p ∈ P ss ∩ P 0 or p ∈ P zs ∩ P 0 ), then complementarity conditions of the type [START_REF] Cocu | Formulation and approximation of quasistatic frictional contact[END_REF] are satisfied which guarantee that the near future sliding speed can only oppose the near future tangential reaction, that (up to second order approximations) the near future reactions remain in the friction cone, and sliding occurs only if those reactions remain on the boundary of that cone (see the third and fourth lines in [START_REF] Rajakumar | The Lanczos algorithm applied to unsymmetric generalized eigenvalue problem[END_REF] and ( 43)) ; on the other hand, if the current reactions or their first order rates indicate a near future stick (p ∈ P d ∪ P zd ∪ P sd , see the fifth line in ( 42)), then the tangential acceleration has to be null ; furthermore, if the sliding speed is null and the normal displacement and reaction, as well as their first order rates, are also null (p ∈ P zz ∩ P 0 , see the fifth line in ( 43)), then the lowest order reaction and tangential displacement rates that may be different from zero (the second order ones) must satisfy a condition of the type of the friction law ; finally, if the current sliding speed is different from zero at the particles in contact with null normal velocity (p ∈ P v = (P z ∪ P s ) ∩ P v = (P zz ∪ P zs ∪ P ss ) ∩ P v , see the sixth line in ( 43)) then, up to second order approximations, the reactions must remain at the vertex of the friction cone or on the appropriate side of its boundary, while if the normal displacement or its first order rate indicate a near future absence of contact (p ∈ P f ∪ P zf , see the last line in ( 43)), then the second order rate of the tangential reaction vanishes.

Similarly to Lemmas 2.6 and 2.7, K r (u, r, v, w, a) can be characterized in the following manner.

LEMMA 2.9. Let U ∈ K u , R ∈ K r (U ), V ∈ K u(U , R), W ∈ K ṙ (U , R, V ) and A ∈ K ü(U , R, V , W ). Then B ∈ K r (U , R, V , W , A) if and only if B.(a -A) ≥ p∈Pv µB N p σ(V T p )(a T p -A T p ) + p∈P 0 µB N p (|a T p | -|A T p |) ∀a ∈ K ü(U , R, V , W ) . 2 (44)
Note that the right hand side of (44) can be written in more detail in the form :

p∈Pv µB N p σ(V T p )(a T p -A T p ) + p∈P 0 µB N p (|a T p | -|A T p |) = p∈Pv µB N p σ(V T p )(a T p -A T p ) + p∈Pzz∩P 0 µB N p (|a T p | -|A T p |) - p∈Pss∩P 0 µB N p σ(R T p )(a T p -A T p ) - p∈Pzs∩P 0 µB N p σ(W T p )(a T p -A T p ) . (45)
In some circumstances the above sets are actually subspaces of IR N : when

P zz ∪ [(P zs ∪ P ss ) ∩ P 0 ] = ∅ , (46)
and, consequently, P z = P zf ∪ P zd ∪ P zs , P zs = P zs ∩ P v , P ss = P ss ∩ P v ,

P v = P zs ∪ P ss , (47)
it can be easily seen that no inequality restriction remains on the admissible second order displacement and reaction rates. This is relevant for the discussion that motivates the definition of these sets (the dynamic solutions in the neighborhood of quasistatic paths), because it means that oscillating solutions belong to those subspaces and can be easily considered in the analysis.

Since it is useful to use complex notation to represent such oscillating solutions, we define the following subspaces of admissible complex second order displacement and reaction rates :

V l C ü = V l C ü (u, r, v, w) ={a ∈ l C N :a i = 0, for i ∈ S D ;
a N p = 0, for p ∈ P d ∪ P s ;

a T p = 0, for p ∈ P d ∪ P zd ∪ P sd } , (48)
V l C r = V l C r (u, r, v, w) ={b ∈ l C N :b i = 0, for i ∈ S F ; b N p = 0 , for p ∈ P f ∪ P zf ; b T p = µ b N p σ(v T p ) , for p ∈ P v ; b T p = 0 , for p ∈ P f ∪ P zf } . (49)

Summary of most important notations

For future reference, we summarize in Table 1 the most important notations relative to the subsets of the set P C of the contact candidate particles. The notation of the sets of admissible (right) values of the (real) kinematic and static variables is summarized in Table 2.

Dynamic divergence from an equilibrium state

In this section we wish to study the existence of dynamic solutions in the neighborhood of a given equilibrium state which might tend to diverge from that state. First we shall consider the possibility that such dynamic solutions might initiate with a velocity discontinuity relatively to the equilibrium state. Then we consider the possibility of acceleration and reaction discontinuities at the initiation of such dynamic solutions. Finally we shall study dynamic solutions along straight paths in the neighborhood of the static equilibrium state with initial conditions arbitrarily close to that state.

Absence of initial velocity discontinuities PROPOSITION No dynamic solution with initial conditions

u(τ ) = u 0 ∈ K u , u-(τ ) = 0 , (50)
may be initiated with a velocity jump, i.e., u+ (τ ) = 0 for all dynamic solutions that satisfy the initial conditions [START_REF] Zeghloul | Phénomènes de glissements partiels découlant de l'usage de la loi de coulomb dans un contact non lubrifié[END_REF].

PROOF. The result follows from Lemma 2.4, by observing that T ( u-(τ )) = 0 and that the mass matrix is positive definite. 2

Initial acceleration and reaction discontinuities

In what concerns acceleration and reaction discontinuities we can prove the necessary conditions given next, where, for each (A, a) ∈ K ü(u 0 , 0) × K ü(u 0 , 0), we use the notations :

m(A, a) def = MA.a , m * (A, a) def = MA.a - p∈Pc(u 0 ) µ[MA] N p |a T p | . ( 51
)
PROPOSITION 3.2. A dynamic solution with initial conditions (50) may be initiated with an acceleration and reaction discontinuity only if there exists an admissible nontrivial right acceleration

A [ A = ü+ (τ ) ∈ K ü(u 0 , 0) , A = 0] such that MA + r 0 ∈ K r (u 0 , 0, A) , (52)
i.e., such that

[m * (A, a) -m * (A, A)] + [r 0 .(A -a) - p∈Pc(u 0 ) µr 0 N p (|a T p | -|A T p |)] ≥ 0 , ∀a ∈ K ü(u 0 , 0) . ( 53
)
PROOF. The right accelerations and reactions satisfy ü+ (τ ) ∈ K ü(u 0 , 0) and r + (τ ) ∈ K r (u 0 , 0, ü+ (τ )), so that

r + (τ ).(a -ü+ (τ )) ≥ p∈Pc(u 0 ) µr + N p (τ )(|a T p | -|ü + T p (τ )|) , ∀a ∈ K ü(u 0 , 0) . (54)
The result follows then by observing that

r + (τ ) = M ü+ (τ ) + r 0 . 2 (55)
COROLLARY 3.3. If (53) holds then :

(i) MA ∈ K ṙ (u 0 , r 0 , 0) , (56)
and consequently 25) and ( 17)]. Since the static reactions r 0 = Ku 0f 0 also belong to the same set, it follows that the reaction jump r + (τ )r 0 = M ü+ (τ ) belongs to K r (u 0 , 0)r 0 , which is necessarily contained in the set K ṙ (u 0 , r 0 , 0) of the admissible right reaction rates at (u 0 , r 0 ). Then the inequality (57) follows immediatly from [START_REF] Pang | Complementarity formulation and existence of solutions of dynamic multi-rigid-body contact problem with Coulomb friction[END_REF] with V = 0.

m * (A, a) ≥ 0 , ∀a ∈ K u(u 0 , r 0 ) . ( 57
) (ii) m * (A, A) ≤ 0 . ( 58
) PROOF. (i) The right reactions r + (τ ) belong to the set K r (u 0 , 0, ü+ (τ )) which is contained in K r (u 0 , 0, 0) = K r (u 0 , 0) [cf. (
(ii) Since, as already mentioned, r 0 ∈ K r (u 0 , 0, 0), it follows from (28) that

r 0 .a ≥ p∈Pc(u 0 ) µr 0 N p |a T p | , ∀a ∈ K ü(u 0 , 0) . ( 59
)
Then taking a = ü+ (τ ) in (59) and taking successively a = 0 and a = 2ü + (τ ) in (54) we get

-r 0 .ü + (τ ) ≤ - p∈Pc(u 0 ) µr 0 N p |ü + T p (τ )| r + (τ ).ü + (τ ) = p∈Pc(u 0 ) µr + N p (τ )|ü + T p (τ )| .
The result (58) follows by adding the above expressions and taking (55) into account. 2 REMARK 3.4. (i) Since [cf. Proposition 3.2] A ∈ K ü(u 0 , 0) and r 0 + MA ∈ K r (u 0 , 0, A), it follows from the definitions ( 24), ( 25) and ( 34) that A necessarily belongs to the set K u(u 0 , r 0 + ∆R), where ∆R[= MA ∈ K r (u 0 , 0)r 0 ] denotes the reaction jump. However A does not have to belong to K u(u 0 , r 0 ), because the right reaction (r 0 + ∆R) may be quite distinct from the (left) static reaction (r 0 ). However if ∆R[= MA] is sufficiently small that r 0 N p + ∆R N p < 0 and

|r 0 T p + ∆R T p | < -µ(r 0 N p + ∆R N p
) , for all p ∈ P d (u 0 , r 0 ) ; r 0 N p + ∆R N p < 0 , for all p ∈ P s (u 0 , r 0 ) ;

|r 0 T p + ∆R T p | < -µ(r 0 N p + ∆R N p ) , for all p ∈ P s (u 0 , r 0 ) such that ∆R T p σ(r 0 T p ) + µ∆R N p < 0 ; then A ∈ K u(u 0 , r 0 + ∆R) ⊂ K u(u 0 ,
m * (A, a) -m * (A, A) ≥ 0 , ∀a ∈ K u(u 0 , r 0 ) . (61)
Consequently, for such small jumps,

A ∈ K u(u 0 , r 0 ) and MA ∈ K ṙ (u 0 , r 0 , A) . (62)
(ii) Clearly no dynamic solution with initial conditions (50) may be initiated with an acceleration and reaction discontinuity if

m * (a, a) > 0 , ∀a ∈ K ü(u 0 , 0) , a = 0 . ( 63
)
(iii) The condition in (ii) is satisfied if the mass matrix M is diagonal. Actually it suffices that (Ma) N p = γ p a N p , with γ p > 0, for all p ∈ P c (u 0 ) and all a ∈ K ü(u 0 , 0). 2

It should be noted that the condition given in Proposition 3.2 is necessary, but not sufficient, for the occurrence of initial acceleration and reaction discontinuities, because it has not been checked so far that such discontinuities can actually be followed by (smooth) dynamic solutions diverging from the equilibrium configuration. In some circumstances it is possible to guarantee this. A result of this type is given next, another one will be given later in Section 3.4.

PROPOSITION 3.5. Let A satisfy (53). If, in addition, for all particles p ∈ P c (u 0 ), either :

A N p < 0 and [r 0 + MA] N p = [r 0 + MA] T p = 0 , or :

A N p = A T p = 0, [r 0 + MA] N p < 0 and |[r 0 + MA] T p | < -µ[r 0 + MA] N p ,
then there exists a dynamic solution with an initial acceleration and reaction discontinuity followed by a smooth dynamic evolution.

PROOF. It suffices to observe that the direction A of the non-vanishing right acceleration ü+ (τ ) is directed towards the interior of a region in the displacement-reaction space where the system is governed by a system of linear ordinary differential equations of the form

M # ü# + K # u # = f 0# ,
where () # denotes a submatrix or a subvector associated with the Free degrees of freedom [i ∈ S F ], plus the normal and tangential degrees of freedom of the currently free contact candidate particles [p ∈ P f (u 0 )] and of those that are currently in contact with null reactions at τ + but have strictly negative normal acceleration at τ + . 2

Divergence along straight paths of smooth dynamic solutions with perturbed initial conditions

Having discussed the occurrence of dynamic solutions initiating at an equilibrium position with a velocity or an acceleration discontinuity, we discuss now the existence of smooth dynamic solutions (u(t), r(t)) starting from (perturbed) initial conditions arbitrarily close to the equilibrium state (u 0 , r 0 ). For t in some right neighborhood of some instant τ ≥ 0 (t ∈ [τ, τ + ∆τ [), we consider such perturbed dynamic solutions in the form

u(t) = u 0 + α(t)V , r(t) = r 0 + β(t)W , (64)
where

V ∈ K u(u 0 , r 0 ) , W ∈ K ṙ (u 0 , r 0 , V ) , (65)
define constant directions in the sets of right admissible displacement and reaction increments ; the function of time α is twice continuously differentiable, α and α are non-negative and nondecreasing in [τ, τ + ∆τ [ ; the function β is continuous, non-negative and non-decreasing in the same interval ; the initial values α(τ ) ≥ 0 and α(τ ) ≥ 0 are arbitrarily small. In these circumstances we have the results given next, where, for each (V , v) ∈ K u(u 0 , r 0 ) × K u(u 0 , r 0 ), the notations (51) are still used and we further denote :

a(V , v) def = KV .v , a * (V , v) def = KV .v - p∈Pc(u 0 ) µ[KV ] N p |v T p | . ( 66
) PROPOSITION 3.6. If ∃λ ≥ 0 and V ∈ K u(u 0 , r 0 ) , V = 0 , such that (λ 2 M + K)V ∈ K ṙ (u 0 , r 0 , V ) , (67)
i.e., such that

λ 2 [m * (V , v) -m * (V , V )] + [a * (V , v) -a * (V , V )] ≥ 0 , ∀v ∈ K u(u 0 , r 0 ) , (68)
then there exists a dynamic solution of the form (64) in [τ, τ + ∆τ [ , with

W = (λ 2 M + K)V (69)
β(t) = α(t) =      α(τ ) cosh [λ(t -τ )] + α(τ ) λ sinh [λ(t -τ )] , if λ > 0 , α(τ ) + α(τ ) (t -τ ) , if λ = 0 , (70)
with ∆τ > 0 sufficiently small, for all arbitrary sufficiently small α(τ ) ≥ 0, α(τ ) ≥ 0. Hence, the equilibrium state corresponding to u 0 and r 0 is dynamically unstable (a divergence instability).

PROOF. The equivalence between (67) and ( 68) results from the characterization (37) of the set K ṙ (u 0 , r 0 , V ). By direct substitution of (64) in ( 1) with α, β and W given by ( 70) and ( 69), it follows that u and r satisfy the dynamic equations of motion. Satisfaction of the frictional contact conditions follows by the same reasons pointed out in the proof of Lemma 2.8. The smallness of α(τ ) ≥ 0, α(τ ) ≥ 0 and ∆τ > 0 is required to keep the solution pair (u(t), r(t)) for t ∈ [τ, τ +∆τ [ within the admissible region of the displacement-reaction space that corresponds to the decomposition of P C at the perturbed state (u(τ ), r(τ )) = (u 0 , r 0 )+α(τ )(V , (λ 2 M +K)V ), this region being the same as or neighbour to the one that corresponds to the decomposition of P C at the unperturbed state (u 0 , r 0 ). See the proof of Lemma 2.8 for the discussion of similar restrictions. 2 PROPOSITION 3.7. A dynamic divergence instability due to perturbed dynamic solutions of the form (64) occurs only if ∃V ∈ K u(u 0 , r 0 ) , V = 0 , such that either : m * (V , V ) = 0 and m * (V , V ).a * (V , V ) ≤ 0 ;

(71) or :

m * (V , V ) = a * (V , V ) = 0 . (72)
PROOF. Let u and r of the form (64) with non-negative α, α and β, with V = 0 and W in the sets (65) correspond to a dynamic divergence instability of the equilibrium state (u 0 , r 0 ). Substituting these functions in the dynamic equation ( 1) we get

[ α(t)M + α(t)K]V = β(t)W . ( 73
)
Using now the characterization (37) of β(t)W ∈ K ṙ (u 0 , r 0 , V ) with v ∈ K u(u 0 , r 0 ) successively equal to 0 and 2V , we get

α(t)m * (V , V ) + α(t)a * (V , V ) = 0 . ( 74
)
The conditions (71), (72) are then necessary conditions for this differential equation to have nontrivial monotonically increasing solutions α with arbitrarily small non-negative initial conditions.

2 COROLLARY 3.8. If m * (v, v) > 0 , ∀v ∈ K u(u 0 , r 0 ) , v = 0 , (75)
then a dynamic divergence instability due to perturbed dynamic solutions of the form (64) : (i) occurs only if

∃V ∈ K u(u 0 , r 0 ) , V = 0 , such that a * (V , V ) ≤ 0 ; ( 76
) (ii) cannot occur if a * (v, v) > 0 , ∀v ∈ K u(u 0 , r 0 ) , v = 0 . 2 (77)
PROPOSITION 3.9. The conditions (67) and ( 68) are necessary and sufficient for a dynamic divergence instability due to perturbed dynamic solutions of the form (64).

PROOF. Only the necessity remains to be proved.

Let V ∈ K u(u 0 , r 0 ), V = 0, such that m * (V , V ) = 0 and m * (V , V ).a * (V , V ) ≤ 0 [cf. (71)].
Then the solution α(t) of ( 74) is necessarily of the form (70) with

λ 2 =      - a * (V , V ) m * (V , V ) , if m * (V , V ).a * (V , V ) < 0 , 0 , if a * (V , V ) = 0 . Hence (67) holds in both cases. If m * (V , V ) = a * (V , V ) = 0 [cf. ( 72
)], then, doing the inner product of (73) with V , we get

α(t)m(V , V ) + α(t)a(V , V ) = β(t)W .V , (78)
where m(V , V ) > 0, a(V , V ) ≥ 0 and, without loss of generality, β(t)W .V may be assumed to satisfy either

(i) β(t)W .V = 0 , ∀t ∈ [τ, τ + ∆τ [ , or (ii) β(t)W .V = 0 , ∀t ∈ [τ, τ + ∆τ [ .
In case (i), α(t) can be monotonically increasing function only if a(V , V ) = 0, and then KV = 0, α(t) is of the form (70) with λ = 0 and (67) holds with λ = 0. In case (ii), it follows from (78) that β(t) is given by

β(t) = α(t)m(V , V ) + α(t)a(V , V ) W .V , with W .V > 0 ,
and then (73) becomes

α(t) MV - m(V , V ) W .V W + α(t) KV - a(V , V ) W .V W = 0 .
It follows that either

MV - m(V , V ) W .V W = 0, and, consequently, KV - a(V , V ) W .V W = 0 , or MV - m(V , V ) W .V W = 0, and α(t) α(t)
is a non-negative constant independent of time.

In both cases (67) holds, which completes the proof. 2

In the following we wish to interpret the above necessary and sufficient conditions for divergence dynamic instabilities in terms of the properties of some tangent or effective mass and stiffness matrices.

For u ∈ K u and r ∈ K r (u) we denote by S * = S * (u, r) the set of degrees of freedom that may be (right) active in that state, i.e., those corresponding to Free velocity components [v i , i ∈ S F ], to normal and tangential velocities of the contact candidate particles that are currently free [v T p and v N p , p ∈ P f (u)], to normal and tangential velocities of the contact candidate particles that are currently in contact with zero reaction [v T p and v N p , p ∈ P z (u, r)] and to tangential (slip) velocity components of the contact candidate particles that are currently in contact with a reaction on the friction cone [v T p , p ∈ P s (u, r)]. Those (right) velocities must belong to the set

K * u(u, r) = {v * ∈ IR N * : v N p ≤ 0 , for p ∈ P z ; v T p σ(r T p ) ≤ 0 , for p ∈ P s } , (79)
where N * = N * (u, r) is the number of degrees of freedom in S * = S * (u, r). We also define the N * × N * matrices M * = M * (u, r) and K * = K * (u, r) with the structure

M * =        M F,F M F,f M F,z M F,sT M f,F M f,f M f,z M f,sT M z,F M z,f M z,z M z,sT M * sT,F M * sT,f M * sT,z M * sT,sT        , K * =        K F,F K F,f K F,z K F,sT K f,F K f,f K f,z K f,sT K z,F K z,f K z,z K z,sT K * sT,F K * sT,f K * sT,z K * sT,sT        (80)
i.e., for all j ∈ S * :

M * ij = M ij and K * ij = K ij ,
for all degrees of freedom i ∈ S * that are not tangential degrees of freedom of particles p ∈ P s ;

(81) M * T p,j = M T p,j + µσ(r T p )M N p,j and K * T p,j = K T p,j + µσ(r T p )K N p,j , for all particles p ∈ P s .

(82)

With these definitions the functionals ( 51) and (66) evaluated at (v, v) may be written :

m * (v, v) = M * S v * .v * - p∈Pz µ[Mv] N p |v T p | , (83)
a * (v, v) = K * S v * .v * - p∈Pz µ[Kv] N p |v T p | , (84)
where

v =              v D v d v sN v *              =              0 0 0 v *              ∈ K u(u, r) , v * =              v F v f v z v sT              ∈ K * u(u, r) (85)
and

K * S = 1 2 (K * + K * T ) , M * S = 1 2 (M * + M * T ) . (86)
We may easily establish the following results.

PROPOSITION 3.10. If

∃λ ≥ 0 and V * ∈ K * u(u 0 , r 0 ) , V * = 0 , such that (λ 2 M * + K * )V * = 0 , (87)
then (67) holds with the same value of λ and

V =              V D V d V sN V *              =              0 0 0 V *              . 2 (88) PROPOSITION 3.11. If M * S (u 0 , r 0 ) is positive semi-definite, and (89)
K * S (u 0 , r 0 ) is positive definite, (90)
then a divergence instability due to perturbed dynamic solutions of the form (64) cannot occur.

PROOF. It suffices to observe that, if λ ≥ 0 and V ∈ K u(u 0 , r 0 ), V = 0, are such that (68) holds, then

0 = λ 2 m * (V , V ) + a * (V , V )=λ 2 M * S V * .V * + K * S V * .V * - p∈Pz µ[λ 2 MV + KV ] N p |V T p | ≥λ 2 M * S V * .V * + K * S V * .V * , (91)
because [λ 2 MV + KV ] N p = W N p ≤ 0, for all p ∈ P z . 2

A stability criterion was derived earlier by Chateau and Nguyen [START_REF] Chateau | Buckling of elastic structures in unilateral contact with or without friction[END_REF] which coincides with the single condition (77) in the present finite dimensional context when the assumption (75) is made. It can be seen (cf. [START_REF] Chateau | Buckling of elastic structures in unilateral contact with or without friction[END_REF][START_REF] Maier | Incremental plastic analysis in the presence of large displacements and physical instabilizing effects[END_REF][START_REF] Bazant | Stability of Structures[END_REF]) that a * (v, v) represents the initial power rate of some external perturbation forces which are added to the given system so as to keep its equilibrium along an evolution initiated at an equilibrium configuration in the direction of an admissible velocity v, while the given applied forces f are kept constant. Power rate statements of the type (77) are sometimes taken as definitions of stability [START_REF] Maier | Incremental plastic analysis in the presence of large displacements and physical instabilizing effects[END_REF][START_REF] Bazant | Stability of Structures[END_REF] and their failure for some admissible velocity direction V [cf. (76)] is then interpreted as an instability of that equilibrium state.

It is important to observe that the scalar (power rate) conditions (71), (72) are only necessary conditions for divergence instabilities of the type (64). The (necessary and) sufficient conditions for such instabilities are the conditions (67), (68) that require the resolution of the inclusion eigenproblem (67) or, equivalently, the variational inequality eigenproblem (68). The reason for the distinction between (68) and ( 71), (72) is related with the non-symmetry of m * (., .) and a * (., .) with respect to the first and second variables. The full consequences of that non-symmetry cannot be perceived with scalar conditions involving only m * (V , V ) and a * (V , V ).

The above considerations have consequences on the distinction between the onset of (87) and the failure of (89), (90). However, for these matrix conditions, there exist additional reasons that make them distinct. We present and discuss all those reasons in the following.

(i) The construction of the specific matrices M * and K * is done by tacitely assuming that all the candidate particles that are currently in contact with zero reaction [p ∈ P z (u 0 , r 0 )] will have a null reaction rate for the tested admissible velocity direction V , and all the contact candidate particles that are currently in contact with a non-vanishing reaction on the friction cone [p ∈ P s (u 0 , r 0 )] will remain so for the tested admissible velocity direction V ; this is done because it eliminates any possible positive (stabilizing) contributions of the particles in P z (u 0 , r 0 ) to λ 2 m * (V , V ) + a * (V , V ), and because instabilizing effects may come only from the particles in P s (u 0 , r 0 ) that, for the considered direction V , do slide and do remain in P s (u 0 , r 0 ) ; this has also the advantage that the reaction rate obtained by computing W = (λ 2 M + K)V , with V given by (88), does belong to the set K ṙ (u 0 , r 0 , V ), no additional checks on this being needed; of course, only the resolution of the full non-linear eigenproblem (67) (i.e., the consideration of all admissible velocities and their associated effective mass and stiffness matrices) leads to a fully satisfactory sufficient (and necessary) condition.

(ii) Assuming (for simplicity of the discussion) that M is diagonal (hence M * S = M * is diagonal and positive definite and the analysis of the conditions in Propositions 3.10 and 3.11 reduces to the analysis of the eigenvalues of K * and K * S ), the distinction between the failure of (90) and the onset of (87) results also from the inequality [START_REF] Bazant | Stability of Structures[END_REF] minimum eigenvalue of K * S ≤ minimum real part of any eigenvalue of K * .

(iii) The vectors V * considered in (87) are real eigenvectors which in addition must belong to the set K * u(u 0 , r 0 ) ; consequently, minimum real part of any eigenvalue of K * ≤ any real eigenvalue of K * corresponding to a real eigenvector V * ∈ K * u(u 0 , r 0 ) .

As a result of (i), (ii) and (iii) above, failure of (90) cannot occur later (and frequently occurs much earlier) than the onset of (87). This will be illustrated in the numerical examples given later. In those examples, the main reason for the above distinction results from (ii) because, in typical situations, the peculiar set P z (u 0 , r 0 ) is empty, because in those examples the instabilizing effects of the initiation of sliding are clear for the particles in P s (u 0 , r 0 ), and because the orientation of the vectors V * obtained from the resolution of the linear eigenproblem (87) can be chosen such that the restrictions in (79) can be satisfied.

Divergence instabilities, initial acceleration and reaction discontinuities and their smooth continuation.

It is important to observe the following.

REMARK 3.12. The sufficient condition (63) for the absence of reaction and acceleration discontinuities also implies the condition (75) that reduces the necessary condition for divergence instability to (76).

2

In addition, it is possible to establish a result that combines possible reaction and acceleration discontinuities with their smooth continuation by trajectories that diverge from the equilibrium state.

PROPOSITION 3.13. If

∃V ∈ K u(u 0 , r 0 ) , V = 0 , such that MV ∈ K ṙ (u 0 , r 0 , V ) and KV = γMV , for some γ ≥ 0 ,

then there exists an infinite number of solutions of the form (64), with an initial acceleration and reaction discontinuity, followed by a smooth trajectory that diverges from the equilibrium state. For t ∈ [τ, τ + ∆τ [, with ∆τ > 0 sufficiently small, these solutions are characterized by

α(t) = (1/2)(t -τ ) 2 α+ (τ ) + ϕ(t) , (93)
β(t) = [ α+ (τ ) + φ(t)] + γα(t) , ( 94
) W = MV , ( 95
)
where α+ (τ ) > 0 is sufficiently small but arbitrary, and ϕ may be any twice continuously differentiable function satisfying ϕ(τ ) = φ(τ ) = φ(τ ) = 0, and ϕ(t) ≥ 0, φ(t) ≥ 0, for all t ∈ [τ, τ + ∆τ [.

PROOF.

The right acceleration at t = τ is ü+ (τ ) = α+ (τ )V and the corresponding initial reaction discontinuity is r + (τ )r 0 = α+ (τ )MV . The smallness of α+ (τ ) > 0 is required to guarantee a small reaction jump, as discussed in Remark 3.4(i). It is easy to check that the smooth continuation given by ( 64) and ( 93)-(95) does indeed satisfy the governing equations and inequalities (1)-( 5) in [τ, τ + ∆τ [ with ∆τ > 0 sufficiently small. 2

Illustrative examples of small dimension.

EXAMPLE 3.14. Initial acceleration and reaction discontinuities. The two d.o.f. system of Fig. 3 is a modified version of the well known example presented by Klarbring [START_REF] Klarbring | Examples of non-uniqueness and non-existence of solutions to quasistatic contact problems with friction[END_REF]. It consists of a particle of mass m restrained by vertical and inclined springs with stiffness constants K V and K I , respectively. The main modification to the system discussed in [START_REF] Klarbring | Examples of non-uniqueness and non-existence of solutions to quasistatic contact problems with friction[END_REF] is that an additional mass M is connected to the previous one by a massless and inextensible string, which (see below) yields a non-diagonal mass matrix. The reference configuration of the system is its equilibrium configuration under gravity action alone, and the particle displacements u are measured from that reference state. The mass and stiffness matrices of the system are then

M =   m + M sin 2 θ M -M sin θ M cos θ M -M sin θ M cos θ M m + M cos 2 θ M   , K =   K I sin 2 θ K -K I sin θ K cos θ K -K I sin θ K cos θ K K V + K I cos 2 θ K   .
A frictional flat obstacle is considered such that the kinematically admissible half-plane is defined by u N ≤ 0. The equilibrium state studied hereafter results from the application of additional constant forces f 0 to the particle m. We consider equilibrium states such that 3: Modified form of Klarbring's example [START_REF] Klarbring | Contact problems with friction by linear complementarity, Unilateral problems in structural mechanics[END_REF] leading to a non-diagonal mass matrix.

u 0 T = u 0 N = 0 and r 0 = (α, -1) with 0 ≤ α ≤ µ, m f 0 f 0 u u N θ K T T N M K I K M g θ V Fig.
i.e. with reactions in the interior of the friction cone (Case 1, 0 ≤ α < µ) or with non-vanishing reaction on the friction cone (Case 2, α = µ). Note that the analysis for -µ ≤ α ≤ 0 would be similar, and that no initial acceleration and reaction discontinuities are possible when r 0 = 0 and u 0 N ≤ 0.

CASE 1

[static reaction in the interior of the friction cone, 0 ≤ α < µ]. First we observe that it is possible to show that no reaction jump is possible towards a contact state with non-vanishing reaction and possible slip to the right (r + N < 0 and r + T = µr + N ). The possible cases of initial acceleration and reaction discontinuities are discussed next (see Fig. 4).

M < 0 NT TT α jump if M r + > 0 µ + 0 N 0 < NT α < 0 r N < α 0 µ 0 T 0 r TT K + + M r growth if exponential < 0 M α jump if K NT TT T r N =-

CASE 1A

[jump towards a contact state with non-vanishing reaction and possible slip to the left :

ü+ T < 0 , ü+ N = 0 , r + N < 0 , r + T = -µr + N ]. From the jump equations   M T T M T N M N T M N N     ü+ T ü+ N   =   r + T -r 0 T r + N -r 0 N   we get ü+ T = µ -α M T T + µM N T and r + N = - M T T + αM N T M T T + µM N T ,
so that a necessary condition for such a jump is that M T T + αM N T < 0. Notice that this implies that α > 0. Notice also that the jump occurs with a reduction of the absolute values of both the tangential and the normal reaction. For the example we are dealing with, condition

M T T + αM N T < 0 gives m M + sin θ M (sin θ M -α cos θ M ) < 0 ,
which is satisfied for angles θ M in an interval strictly contained in ]0, tan -1 µ[. The size of this interval increases as m/M decreases. Considering the equations of motion for a smooth near future dynamic evolution of the system we get, after eliminating the normal degree of freedom,

(M T T + µM N T )ü T (t) + (K T T + µK N T )u T (t) = µ -α .
Since M T T + µM N T < 0, the particle m begins to slide towards the left with a sinusoidal, an exponential, or a quadratic evolution in time, when K T T + µK N T < 0, > 0, or = 0, respectively. Note that for the system of Fig. 3,

K T T + µK N T = K I sin θ K (sin θ K -µ cos θ K ) is negative for θ K ∈]0, tan -1 µ[.

CASE 1B

[jump towards a contact state with no reaction : ü+ N ≤ 0 , r + = (0, 0)]. The acceleration discontinuity is characterized by

ü+ T = - M T N + αM N N det M and ü+ N = M T T + αM N T det M ,
with, necessarily, M T + αM N T ≤ 0 (det M > 0). From the equations of motion for a smooth evolution in the near future we conclude that : (1) if M T T + αM N T < 0 the normal acceleration is strictly negative (ü + N < 0) and the particle m looses contact (Proposition 3.5 could be applied to this case); (2) if M T T + αM N T = 0 the particle m cannot begin to slide towards the right but it can loose contact or begin sliding towards the left with a sinusoidal evolution in time if

K T T + αK N T < 0.

CASE 2

[non-vanishing reaction on the friction cone, α = µ]. Similarly to Case 1, it can be shown that no jump is possible towards a contact state with non-vanishing reaction and possible slip to the right (r + N < 0 and r + T = µr + N ). The possible cases of initial acceleration and reaction discontinuities are the following (see Fig. 5).

CASE 2A

[jump towards a contact state with non-vanishing reaction and possible slip towards the left : states having reactions with absolute values larger than the equilibrium values. Moreover, ü+ T < 0 is arbitrary but sufficiently small that the normal reaction r + N remains negative. After elimination of the normal degree of freedom, the equations of motion for a smooth near future evolution yield

ü+ T < 0 , ü+ N = 0 , r + N < 0 , r + T = -µr + N ].
(K T T + µK N T )u T (t) = 0 .
If K T T +µK N T = 0, no such smooth near future evolution is possible, while if K T T +µK N T = 0, there are infinitely many dynamic evolutions along which the system may leave its equilibrium (Proposition 3.13 could be applied to this case).

CASE

[jumps towards a no reaction contact state :

ü+ N ≤ 0 , r + = (0, 0)]. The discontinuity is ü+ T = - M T N + µM N N det M and ü+ N = M T T + µM N T det M ,
with, necessarily, M T T + µM N T ≤ 0. From this and the equations of motion for a near future smooth evolution we conclude that : (1) if M T T + µM N T < 0 the particle m looses contact (Proposition 3.5 could be applied to this case); (2) if M T T + µM N T = 0 the particle m cannot begin sliding towards the right, but it can loose contact if K T T + µK N T < 0, or it can begin sliding towards the left in an arbitrary way if K T T + µK N T = 0. For the structure in Fig. 3 and for fixed m/M and µ a possible configuration that fulfills the latter conditions is qualitatively represented in Fig. 6. Finally we observe that the equilibrium state of the system of Fig. 3 may also be unstable by divergence due to smooth dynamic solutions of the form (64) with perturbed initial conditions (see Proposition 3.6). This can happen only when the system is in a state of impending slip with non-vanishing normal reaction (r 0 N < 0 and |r 0 T | = -µr 0 N ). In these circumstances, this kind of instability occurs in the above case 2 (α = µ) if and only if the system parameters satisfy either :

M T T + µM N T = 0 and (K T T + µK N T ).(M T T + µM N T ) ≤ 0 , or : M T T + µM N T = K T T + µK N T = 0 .
u D , u F , u T and u N . Following the notation of Section 2.1 :

Ψ =              2P/K 0 0 0              , Φ =              0 0 µP 0             
, the stiffness and mass matrices, and the vector of applied forces (after homogeneization of the displacement u D ) are, respectively,

K = K        1 -1 0 0 -1 2 0 -1 0 0 1 0 0 -1 0 1        , M = M        1 0 0 0 0 1 0 0 0 0 1 0 0 0 0 1        , f = Φ -KΨ =              -2P +2P µP 0              .
We study the stability of the equilibrium state

u 0 =              0 P/K 0 0              , r 0 =              P 0 -µP -P              .
The corresponding sets of admissible right velocities and right reaction rates are here :

K u(u 0 , r 0 ) = {(v D , v F , v T , v N ) : v D = 0, v N = 0, v T ≥ 0} , K ṙ (u 0 , r 0 , V ) = {(w D , w F , w T , w N ) : w F = 0, -w T + µw N ≤ 0, (-w T + µw N )|V T | = 0} .
Application of Propositions 3.6 and 3.9 leads here to the following system of equations and conditions on the unknowns λ, V F , V T , W D , W T and W N (note that

V D = V N = W F = 0) -KV F = W D , (λ 2 M + 2K)V F = 0 , (λ 2 M + K)V T = W T , -KV F = W N , λ ≥ 0 , V F = 0 or V T = 0 , V T ≥ 0 , -W T + µW N ≤ 0 , (-W T + µW N )|V T | = 0 .
It can be easily checked that no solution to this system of equations and conditions exists, so that divergence instability due to perturbed solutions of the form (64) is excluded.

In this simple example we illustrate the distinction between the sufficient (and necessary) condition (67) for divergence instability and the (necessary) power rate condition (76) [note that, since M is diagonal, condition (75) necessarily holds]. In this example, the comparison between those conditions reduces to the analysis of the smallest real eigenvalue of K * and K * S , which are here [recall (81)-( 82) and (86)]

K * = K   2 0 µ 1   , K * S = K   2 µ/2 µ/2 1   .
Note that in this example no a priori assumption is involved in letting V T = 0 in the construction of K * [recall observation (i) in Section 3.3], because the stiffnes matrix K is positive definite in K u(u 0 , r 0 ) and because there is only one particle in contact, so that the only way divergence instability might occur would be with V T = 0. Note also that, again because there is only one particle in contact, the only sign restriction on the admissible eigenvectors of K * and K * S (V T ≥ 0) can be immediately satisfied [recall observation (iii) in Section 3.3]. The minimum eigenvalue of K * is 1, independently of µ, and the minimum eigenvalue of

K * S is (3 -1 + µ 2 )K/2 which becomes non-positive for µ ≥ 2 √ 2. The corresponding eigenvector V * S is V * S =    -β 1    , β = µ 1 + 1 + µ 2 .
We stress the fact that the direction

V S =              0 -β 1 0              ,
for which the necessary condition (76) is satisfied is not a direction along which a divergence instability might initiate for µ ≥ 2 √ 2 : the eigenvalues of K * are always positive and, as seen above, the sufficient condition (67) is never satisfied. 2

Growth of perturbed dynamic solutions in the neighborhood of quasistatic straight paths

We consider now the quasistatic straight portions of quasistatic evolutions considered in Section 2.3. Throughout this section we assume thus the linear variation (39) of the applied forces and that the quasistatic evolution (u 0 (t), r 0 (t)) has the linear variation given by [START_REF] Progri | Analyse de glissements dans un contact sec[END_REF]

in some interval [τ, τ + ∆ 0 τ [.
In some right neighborhood [τ, τ + ∆τ [ of τ , with 0 < ∆τ ≤ ∆ 0 τ , we seek perturbed dynamic solutions in the form

u(t) = u 0 (t) + α(t)A , r(t) = r 0 (t) + β(t)B , (96)
where

A ∈ K ü(u 0 (τ ), r 0 (τ ), u0+ (τ ), ṙ0+ (τ )) , B ∈ K r (u 0 (τ ), r 0 (τ ), u0+ (τ ), ṙ0+ (τ ), A) , (97)
α is a twice continuously differentiable function such that α and α are non-negative and nondecreasing in that interval, and the function β is continuous, non-negative and non-decreasing in the same interval ; the initial values α(τ ) ≥ 0, α(τ ) ≥ 0 are arbitrarily small. For simplicity of notation, we shall denote hereafter the first set in (97) by K ü(• • •) and the second one by

K r (• • • , A).
In these circumstances we have the results given next, where, for each (A, a) 51), (66) and ( 44)]

∈ K ü(• • •) × K ü(• • •), we use the notations [recall (
m * * (A, a) def = MA.a - p∈Pv µσ( u0+ T p (τ ))[MA] N p a T p - p∈P 0 µ[MA] N p |a T p | (98)
a * * (A, a) def = KA.a - p∈Pv µσ( u0+ T p (τ ))[KA] N p a T p - p∈P 0 µ[KA] N p |a T p | . ( 99
) PROPOSITION 4.1. If ∃λ ≥ 0 and A ∈ K ü(• • •) , A = 0 , such that (λ 2 M + K)A ∈ K r (• • • , A) , (100)
i.e., such that

λ 2 [m * * (A, a) -m * * (A, A)] + [a * * (A, a) -a * * (A, A)] ≥ 0 , ∀a ∈ K ü(• • •) , (101)
then there exists a dynamic solution of the form (96) in [τ, τ + ∆τ [ , with

B = (λ 2 M + K)A , (102)
β(t) = α(t) =      α(τ ) cosh [λ(t -τ )] + α(τ ) λ sinh [λ(t -τ )] , if λ > 0 , α(τ ) + α(τ ) (t -τ ) , if λ = 0 , (103)
with ∆τ ≤ ∆ 0 τ , positive and sufficiently small, for all arbitrary, sufficiently small α(τ ) ≥ 0, α(τ ) ≥ 0.

PROOF. The proof follows the same lines of Proposition 3.6 and is omitted.

2 REMARK 4.2.
Similarly to what was done in Section 3.3 (Propositions 3.10 and 3.11), it is also possible to establish necessary and sufficient conditions for the existence of these growing solutions in the neighborhood of the quasistatic path, in terms of the properties of some matrices M * * and K * * . These are here N * * × N * * matrices, where N * * is the total number of degrees of freedom that may be active in the interval [τ, τ + ∆ 0 τ [, i.e., those corresponding to free degrees of freedom (i ∈ S F ), normal and tangential degrees of freedom of the particles p ∈ P f ∪ P zf ∪ P zz , tangential degrees of freedom of the particles p ∈ P zs ∪ P ss . 2

We assume now that P zz ∪ [(P zs ∪ P ss ) ∩ P 0 ] = ∅ [recall [START_REF] Raous | Stress waves in a sliding contact, Part 2: modelling[END_REF] and the related definitions ( 48), (49) of the subspaces of admissible complex second order displacement and reaction rates].

In this situation we seek growing dynamic solutions in the neighborhood of the quasistatic path (u 0 (t), r 0 (t)) [START_REF] Progri | Analyse de glissements dans un contact sec[END_REF] in the form

u(t) = u 0 (t) + Re[α(t)A] , r(t) = r 0 (t) + Re[β(t)B] , (104)
with

A ∈ V l C ü , B ∈ V l C r . (105)
In the first of these subspaces of l C N we define the sesquilinear forms given by m(A, a) ü . In the above expressions, the inner product in the first term of the right hand side is the inner product of complex vectors [x.y = N i=1 x i ȳi , for each (x, y) ∈ l C N × l C N ] and, for each y ∈ l C, ȳ denotes the complex conjugate of y. In the same circumstances we denote by S = S(u 0 (τ ), r 0 (τ ), u0+ (τ ), ṙ0+ (τ )) the set of the active degrees of freedom along the time interval [τ, τ + ∆ 0 τ [, i.e., those corresponding to :

• free degrees of freedom [i ∈ S F ] ;

• normal and tangential degrees of freedom of the contact candidate particles that are free at time τ , or that are in contact with zero reaction at time τ but have a negative normal velocity [p ∈ P f ∪ P zf ] ;

• tangential degrees of freedom of the contact candidate particles that, at time τ , are in contact with a non-vanishing tangential velocity and have either a non-vanishing reaction on the friction cone, or a vanishing reaction and a non-vanishing reaction rate on the friction cone [p ∈ P zs ∪P ss ].

Following a procedure similar to the one used in Section 3.3, we consider the N × N matrices constructed in the following manner ( N denotes the number of active degrees of freedom in S).

For all j ∈ S : M i,j = M i,j and K i,j = K i,j , for all degrees of freedom i ∈ S that are not tangential degrees of freedom of particles p in P zs ∪ P ss ; M T p,j = M T p,jµσ( u0+ T p (τ ))M N p,j , and K T p,j = K T p,jµσ( u0+ T p (τ ))K N p,j , for all particles p in P zs ∪ P ss .

In these circumstances we have the result. C with Re(λ) > 0 such that one of the following equivalent conditions holds :

(i) ∃A ∈ V l C ü , A = 0 , such that (λ 2 M + K)A ∈ V l C r ; ( 106
) (ii) ∃A ∈ V l C ü , A = 0 , such that λ 2 m(A, a) + a(A, a) = 0 , ∀a ∈ V l C ü ; ( 107
) (iii) ∃ A ∈ l C N , A = 0 , such that [λ 2 M + K] A = 0 ; ( 108
)
then there exists a dynamic solution of the form (104) in [τ, τ + ∆τ [, with

B = (λ 2 M + K)A (109) β(t) = α(t) = α(τ ) cosh [λ(t -τ )] + α(τ ) λ sinh [λ(t -τ )] , (110)
with ∆τ ≤ ∆ 0 τ , positive and sufficiently small, for all arbitrary, sufficiently small |α(τ

)| > 0 and | α(τ )| > 0. 2
We observe that the above statement is not interpreted as an instability result because it refers only to a portion of a quasistatic evolution that is also a dynamic solution in some possibly finite time interval [τ, τ + ∆ 0 τ [. The above statement does correspond to an instability if (u 0 (t), r 0 (t)) given by ( 41) solves the quasistatic problem ( 32), ( 2)-( 6) 1 , and consequently the dynamic problem ( 1)-( 6), for all t ∈ [τ, +∞[ . 111) and ( 46) hold are those corresponding to steady sliding solutions (u 0 (t), r 0 (t)) for t ≥ τ , with non-vanishing reaction forces and non-vanishing sliding speeds at the particles actually in contact. In such situations, the applied forces Φ(t) (recall Section 2.1) are constant, the prescribed displacement rates Ψ(t) are also constant, r 0 (t) is a constant vector such that P z (u 0 (t), r 0 (t)) = P d (u 0 (t), r 0 (t)) = ∅, and, for all t ≥ τ , u u u 0 (t) = Ψ(t) + u0 (t) is a constant vector that belongs to the kernel of the stiffness matrix K : the quasistatic rate balance equations (40) (recall also Section 2.1) effects due to the axial force on the beam and due to gravity action on the mass M are neglected a priori. The mass and stiffness matrices are :

K u0 = ḟ + ṙ0 = -K Ψ reduce to K u u u 0 = 0 .
M =                 M D M 0 0 0 0 H 2 M + R 2 M L 2 - H M L 0 0 - H M L 1 0 0 0 0 1                 , K =                 1 H B L -1 0
H B L 4R 2 B + H 2 B L 2 - H B L 6 R B L 2 -1 - H B L 1 0 0 6 R B L 2 0 12 R B L 2                
.

We consider now the dynamic stability of the steady sliding solution

u u u 0 (t) =        u 0 D (τ ) + u0 D .(t -τ ) u 0 θ (τ ) u 0 T (τ ) + u0 D .(t -τ ) 0        , r 0 (t) =        µr 0 N (τ ) 0 -µr 0 N (τ ) r 0 N (τ )        , where u 0 θ (τ ) = µ M g H B 4EAL 1 + 3 2 µ H B L R B L 2 , r 0 N (τ ) = - M g EA 1 + 3 2 µ H B L , u 0 D (τ ) -u 0 T (τ ) = - µ M g 1 + 1 4
H B R B 2 EA 1 + 3 2 µ H B L .
Note that the prescribed displacement u D (t) has not been homogeneized in the expressions above. One of the initial displacements u 0 D (τ ) or u 0 T (τ ) is arbitrary. After elimination of the D and N degrees of freedom, we get the eigenproblem of Eq. ( 108) in Proposition 4.3, where the 2x2 matrices M and K are here

M =    1 - H M L - H M L H 2 M + R 2 M L 2    , K =     1 - H B L -6µ R B L 2 σ( u0 D ) - H B L 4R 2 B + H 2 B L 2     .
The characteristic equation for this system is a bi-quadratic equation in λ. For instance, in the specific simpler case of H M = 0, the sufficient conditions for the occurrence of dynamic instabilities are :

-Flutter : 1 + 4 R B R M 2 + H B R M 2 2 -8 R B R M 2 . 2 -3µ H B L σ( u0 D ) < 0, and -Divergence : 2 -3µ H B L σ( u0 D ) ≤ 0.
From these inequalities we conclude that flutter occurs only when σ( u0 D ) = -1, i.e. only when the beam is pulled towards the left, and that divergence occurs only when σ( u0 D ) = +1, i.e. only when the beam is pushed towards the right. Note however that in the latter case, it is possible to show that when the condition for divergence is satisfied, a steady sliding solution with (downward) gravity and prescribed displacement u D does not exist : the situation is similar to the one reported by Nguyen [34] for Klarbring's example with prescribed displacement. Just for illustration purposes we present in Fig. 9 a phase plane plot of the evolution of the tangential displacement and velocity of the contact point of the pin as a result of the flutter instability of the steady sliding solution ; note that the variables plotted are not (u T (t), uT (t)) but the translated variables, (u T (t)-u0 D .(t-τ ), uT (t)-u0 D ), i.e., the displacements and velocities relative to an observer that follows the motion of the block of mass M D . In Fig. 10 the phase plane plot of the variables (u N (t), uN (t)) is shown. The development of the flutter instability and the convergence of the solution to a limit cycle is clear in both of those plots. The impacts are assumed to be perfectly plastic : the dissipation index δ [recall [START_REF] Klarbring | Contact friction, discrete mechanical structures and mathematical programming, Contact problems : theory, methods and applications[END_REF], [START_REF] Lebon | Multibody contact problem including friction in structure assembly[END_REF]] is taken to be 1. The other data used is listed in the caption of Fig. 9. ) (and previously with R. Progri and M. Mouwakeh [START_REF] Progri | Analyse de glissements dans un contact sec[END_REF][START_REF] Mouwakeh | Etude quantitative des phénomènes de glissements dans un contact sec à deux corps par comparaison avec la propagation d'une fissure interfaciale[END_REF]) a deformable block of polyurethane is first pressed against a rigid plane Araldite surface (see Fig. 11) with a prescribed displacement (U 0 N ) which insures the initial resultant normal force F 0 N to be -55N. The block is then loaded in the tangential direction by prescribing a slow tangential motion to the right (≃ 0.05 mm/s) to the Araldite base plate. The normal resultant will change during the tangential motion while the normal prescribed displacement U 0 N is kept constant. In the course of this tangential loading process, fast isolated stress waves can be observed in the block by photoelasticity. These waves move from right to left at a speed that is nearly a thousand times faster than the prescribed speed of the Araldite base plate. The fast propagation of these waves induces jumps in the measured tangential displacements on the left and right ends of the block and in the total tangential and normal forces as given on Fig. 12 (see Zeghloul [START_REF] Zeghloul | Étude des phénomènes d'adhérences et de glissements dans un contact entre solides: approche expérimentale et modélisation[END_REF] and Zeghloul-Villechaise [START_REF] Zeghloul | Stress waves in a sliding contact, Part 1: experimental study[END_REF]). Waves of two kinds may occur :

-partial waves which do not reach the left end of the block and correspond to a partial release of the stresses : this effect is emphasized either with other dimensions of the block (H/L=0.7 for example) or with other velocities of the plate (v=4.1 10 -2 mm/s for example), -total waves that are characterized by a jump of the tangential displacement of the left end of the block, and which lead to a larger release of global tangential forces (see Fig. 12) and of global stresses and displacements.

Fig. 12: Experimental measurement of the total tangential force (from Zeghloul [START_REF] Progri | Analyse de glissements dans un contact sec[END_REF]).

A response involving jumps is thus obtained for a quasistatic loading.

The possible correlation between these experimental observations and friction-induced instability phenomena is discussed next, using the finite element method.

Finite element model

The block is assumed to be in a plane stress state and to have linear elastic behaviour. In all the meshes, linear P1 finite elements are used because in this case the convex set K u provides an internal approximation for the corresponding set of admissible displacements of the continuum body, but higher order elements could also be used. Various meshes are constructed in the same manner so that the symmetry in all directions is preserved. The number of contact nodes varies from 11 to 81. A 21 contact node mesh is shown in Fig. 13.

A B

Fig. 13: Mesh with 21 contact nodes.

Quasistatic solution

To simplify the following formulae, the quasistatic solution, previously denoted by u 0 (t) at time t, is now denoted by u k at the discrete time instants t k .

Recalling that V u = {u ∈ IR N : u i = 0 , for i ∈ S D }, the discretized quasistatic problem (2)-(5) (31)-( 33) is equivalent to the discrete incremental variational formulation presented by

The memory due to the velocity formulation of the friction is taken into account in (117) by u k T p , the previous tangential contact displacement.

To ensure the maximum accuracy, the contact forces are directly computed from the residues of equilibrium. Details about the algorithm can be found in [START_REF] Raous | Numerical methods for frictional contact problems and application[END_REF] and [START_REF] Lebon | Multibody contact problem including friction in structure assembly[END_REF]. A diagonal process is used for the fixed point method on g, i.e., coarse resolutions are carried out with the first values of g. The number of iterates on g is usually less than 8 and the total number of S.O.R.P. iterates is between one and two times the number of degrees of freedom.

The stability analysis presented next is done at each time t k of the quasistatic incremental solution (u k , r k ).

Numerical characterization of instabilities

For each straight portion of the quasistatic solution, i.e. between two successive changes of contact status, the tangent mass and stiffness matrices are constructed. Within such portions and in the typical situations which occur during the numerical resolution of the quasistatic problem, the following sets of contact candidate nodes are empty : the peculiar set P z (u k , r k ) of the nodes currently in contact with null reaction ; the peculiar set P ss0 (u k , r k , uk , ṙk ) of the nodes currently in contact with non-null reaction on the friction cone, with reaction rates tangent to that cone and null sliding speed ; and the set P sd (u k , r k , uk , ṙk ) of the nodes currently in contact with non-null reaction on the friction cone, with null sliding speed and reaction rates directed towards the interior of the friction cone. Note that in practice, because of the space and time discretizations involved in the numerical resolution of the quasistatic problem and because of the finite digit accuracy we never need to worry about the above sets, i.e., at each discrete time point t k we "always" have

P z (u k , r k ) = ∅ P s (u k , r k ) = P ss (u k , r k , uk , ṙk ) ∩ P v (u k , uk ) i.e. P C = P f (u k ) ∪ P d (u k , r k ) ∪ [P ss (u k , r k , uk , ṙk ) ∩ P v (u k , uk )]
and there is no distinction either between the matrices M * , M * * and M or between the matrices K * , K * * and K. Hereafter we shall denote them simply by M * and K * .

These matrices are obviously constructed from the finite element stiffness and mass matrices K and M by removing the equations corresponding to the non-active degrees of freedom [recall the definition of S * = S * (u k , r k ) given in Section 3] and by modifying the coefficients of the rows of the matrices K and M associated with the tangential displacement of the sliding nodes [recall (82)]. Matrices K * and M * are therefore non symmetrical and dependent on µ and they are modified whenever the contact conditions change. When desired, the symmetrized matrices M * S and K * S [see (86)] are also constructed. We observe that matrices M * and M * S were never found to be singular (M * S is always positive definite) in any of the computations involving finite element discretizations performed in the course of this work. Therefore during the quasistatic loading process, i.e., at each time t k , we search for the following : this is necessary for a divergence instability ;

(ii) non-positive real eigenvalues of the matrix K * , or, equivalently, of M * -1 K * satisfying the same constraint as above [recall (87) and (79) with M * non-singular and P z = ∅] ; this is sufficient for a divergence instability ;

(iii) complex eigenvalues of M * -1 K * [recall (108) and [START_REF] Stewart | Existence of solutions to rigid body dynamics and the Painlevé paradoxes[END_REF] with P z = ∅, P ss ∩ P v = P s , M = M * , K = K * and M * non-singular] ; this is sufficient for the existence of growing flutter-type oscillations in the neighborhood of the quasistatic solution.

The search for (i) is done by three different methods ; they all give the same results :

1. minimizing the bilinear form a * (v, v) under constraints (118) and v = 1 and checking the sign of the minimum, 2. applying Cholesky's method to K * S and checking when it fails (K * S is no longer positive definite), 3. computing the smallest eigenvalue of K * S and checking its sign. Method 2 is a direct one and therefore the fastest. The last method is presented because it could be used to search for other eigenvalues, and also because it gives the eigenvector of K * S corresponding to its smallest eigenvalue.

The smallest eigenvalue of K * S is computed as follows : using the Power Method, we first compute the eigenvalue Λ L of K * S with largest modulus. If Λ L were negative, that would be the desired result. If Λ L > 0, we shift all the eigenvalues towards the negative axis by constructing the modified matrix K * S -Λ L I ; a computation of the eigenvalue Λ M with largest modulus of this modified matrix gives the smallest eigenvalue Λ m of K * S as Λ m = Λ M + Λ L and we check its sign.

The search for (ii) and (iii) is done by solving the generalized eigenproblem K * V * = ΛM * V * (Λ = -λ 2 ) with the Lanczos algorithm [START_REF] Rajakumar | The Lanczos algorithm applied to unsymmetric generalized eigenvalue problem[END_REF] along with the double QR algorithm [START_REF] Ralston | A first course in numerical analysis[END_REF]. Most of these computations are performed with the consistent mass matrix M * ; in order to speed up some computations with the most refined meshes, a diagonal mass matrix M is used which also leads to a diagonal matrix M * .

The results will be presented in Section 5.6. Note that in cases (i) and (ii) condition (118) has always been satisfied by the computed eigenvectors.

Results of the quasistatic loading

For µ = 1.1, H/L=0.5 , U 0 N = -0.5mm (which corresponds to F 0 N = -55N) and a mesh with 21 contact nodes, the contact stresses are plotted in Fig. 14 at various stages along the process of increasing the imposed tangential displacement of the araldite base. It can be seen in that figure that, in the course of that tangential loading process, some nodes on the right of the contact surface successively loose contact and also that the sliding zone expands from right to left. In the final phase of the tangential loading process all the nodes still in contact attain sliding and, thereafter, all deformations and stresses in the block remain constant : a steady sliding equilibrium state is reached.

Note that for different friction coefficients µ and different aspect ratios H/L quite different patterns of evolution can be obtained for the contact states. The computed quasistatic evolution of the total tangential contact force (T) is plotted in Fig. 15 relatively to the time variable t. It is a smooth evolution on which it is possible to distinguish a straight initial portion (while all the contact nodes are still stuck) and a final flat portion corresponding to the steady sliding state. The material characteristics used in the numerical calculations are : Poisson's ratio ν = 0.48, and density ρ=1.2 10 -3 g/mm 3 . Note that the value of Young's modulus used in the numerical calculations was identified by comparing the experimental and theoretical results for the early stages of the quasistatic evolution of the total tangential contact force, when all the nodes are still stuck (see Fig. 12 and Fig. 15) : E=5MPa. evolution with loss of energy into some external sink. That mode shape clearly suggests the waves experimentally observed in the polyurethane block. However, as seen in Table 1, the sufficient condition for a divergence instability (a negative eigenvalue of M * -1 K * ) never gets to be satisfied.

In what concerns flutter, it occurs for very high modes as soon as any node begins sliding. However it is likely that significant consequences will arise only when some lower order mode is affected by flutter. Table 1 shows that flutter occurs for a low order mode (the 3rd mode) at the same time as the necessary condition for divergence is first satisfied (a little earlier in other examples). Table 1 also shows that the steady-sliding state attained at the end of the theoretical quasistatic solution is unstable by flutter (see the last row in that table). In Fig. 17 the time evolution of both experimental and theoretical tangential forces is plotted. On the theoretical curve, the dot indicates the time at which the necessary condition for divergence instability and the sufficient condition for growing flutter-type oscillations are first satisfied. Clearly those conditions are satisfied too early as compared to the first experimental observations of jumps : both those conditions are satisfied for U T =-1.35mm while the first small jump is observed for U T =-2.2mm and a second jump is observed for U T =-2.71mm (see also Table 1).

In the next section we present a preliminary study on the effect that the presence of viscous damping may have on the stability analysis and results.

Effect of viscous damping

As seen above, flutter occurs too early in comparison with the experiment. On the other hand, a simplified form to introduce viscous damping effects in the problem is to consider a damping matrix of the Rayleigh type αM + βK in the dynamic equations of motion and the dynamic stability analysis. In this manner it is possible to delay the flutter occurrence (and even to eliminate it with sufficiently strong damping).

The frequencies of the 18 first natural modes of the structure range between 454 Hz and 1530 Hz in its initial phase of elastic behavior, i.e., when all the nodes are stuck. The damping ratio ξ depends on the frequencies and is known to be decreasing in that range for the polyurethane.

Fig. 1 :

 1 Fig. 1: Finite dimensional system in frictional contact with a fixed flat surface ; prescribed Displacements (i ∈ S D ), Free displacements (i ∈ S F ), Normal and Tangential displacements of the Contact candidate particles (p ∈ P C ) .

def=

  {p ∈ P c (u) : v N p = 0, v T p = 0} [particles in contact with normal and tangential velocities equal to 0] P v (u, v) def = {p ∈ P c (u) : v N p = 0, v T p = 0} [particles in contact with vanishing normal velocity and non-vanishing tangential velocity].

P

  d (u, r) def = {p ∈ P c (u) : r N p < 0 and |r T p | < -µ r N p } [particles in contact with reaction strictly inside the friction cone and consequent vanishing (right) displacement rate] P s (u, r) def = {p ∈ P c (u) : r N p < 0 and |r T p | = -µ r N p } [particles in contact with non-vanishing reaction on the boundary of the friction cone and consequent possible slip in the near future].

Fig. 2 :

 2 Fig.2: Admissible normal displacement (u N ), normal and tangential reactions (r N and r T ), and their first order (right) rates for a contact candidate particle p ∈ P C ; nomenclature for the subsets of P z and P s .

  r 0 ) and, consequently [cf. (57),(58)], m * (A, A) = 0 (60) and also [cf. (57),(60)],

Fig. 4 :

 4 Fig. 4: Possible initial reaction jumps followed by smooth near future evolutions in Case 1 [static reaction in the interior of the friction cone].

Fig. 5 :

 5 Fig. 5: Possible initial reaction jumps followed by smooth near future evolutions in Case 2 [non-vanishing static reaction on the friction cone].

def=

  MA.a -p∈Pzs∪Pss µσ( u0+ T p (τ ))[MA] N p āT p = MA.a + p∈Pzs µσ( ṙ0+ T p (τ ))[MA] N p āT p + p∈Pss µσ(r 0 T p (τ ))[MA] N p āT p a(A, a) def = KA.a -p∈Pzs∪Pss µσ( u0+ T p (τ ))[KA] N p āT p = KA.a + p∈Pzs µσ( ṙ0+ T p (τ ))[KA] N p āT p + p∈Pss µσ(r 0 T p (τ ))[KA] N p āT p for each pair (A, a) ∈ V l C ü × V l C

PROPOSITION 4 . 3 .

 43 [START_REF] Raous | Stress waves in a sliding contact, Part 2: modelling[END_REF] holds and if ∃λ ∈ l

2 EXAMPLE 4 . 6 .Fig. 8 :

 2468 Fig.8: A pin-on-flat system subjected to gravity forces, prescribed displacements (u D ) and frictional contact reactions on the boundary of the kinematically admissible half-plane u N ≤ 0.

Fig. 9 :

 9 Fig. 9: Phase plane plot of the evolution of the tangential displacement and velocity of the contact point of the pin for u0 D =-0.015, µ=2.2, M g EA =0.01, R B L =0.07, H B L =0.1, R M L =0.15 and H M L =0.05.

Fig. 10 :

 10 Fig. 10: Phase plane plot of the evolution of the normal displacement and velocity of the pin.

2 5.Fig. 11 :

 211 Fig. 11: Experimental assembly scheme

  (i) non-positive eigenvalues of the matrix K * S [recall (76) and (84), with positive definite M * S and with P z = ∅] corresponding to real eigenvectors V * satisfying the constraint [cf. (79) with P z = ∅] σ(V * T p ) = -σ(r k T p ) for all p ∈ P s (u k , r k ); (118)

Fig. 14 :

 14 Fig. 14: Evolution of the contact stresses and of the contact conditions along the contact zone.

Fig. 15 :

 15 Fig. 15: Numerical total tangential force.

Fig. 17 :

 17 Fig. 17: Evolution in time of experimental and theoretical total tangential forces.

  PROPOSITION 4.4. If (111) holds and the assumptions of Proposition 4.1 or 4.3 hold, then the solution (u 0 (t), r 0 (t)), t ≥ τ , is dynamically unstable. If Im(λ)=0, we have a divergence instability and if Im(λ) = 0 we have a flutter instability. 2 REMARK 4.5. Typical situations in which (

	 
		(111)

Acknowledgements : We are grateful to Prof. Bernard Villechaise and D r . T. Zeghloul, from LMS Poitiers-Angoulême, for the experimental results given on figure 12 and 17 and for the discussions we had about the experimental study.

This study has been partially supported by a project CNRS/JNICT and a Réseau Formation Recherche (MENESR and JNICT).

 6:A possible configuration for the system of Fig. 3 to have an infinite number of smooth evolutions after an initial reaction and acceleration jump. The value θ M corresponding to M T T +µM N T = 0 is a solution (if it exists) to the equation m/M +sin θ M (sin θ M -µ cos θ M ) = 0 ; depending on m/M and µ, this equation may have at most two solutions in the interval [0, π/2]. 2 EXAMPLE 3.15. Necessary and sufficient conditions for (smooth) divergence instability. We consider the plane system represented in Fig. 7, where each of the three particles has mass M and each of the three linear springs has stiffness K. The displacement u D has the constant prescribed value 2P/K and the constant force Φ T = µP is applied on the direction of u T . The two vertical springs remain orthogonal to the fixed contact surface when the system deforms, and gravity acceleration is not considered. The generalized displacements are, by this order, Cocu et al ( [START_REF] Cocu | Existence d'une solution du problème quasi-statique de contact unilatéral avec frottement non local[END_REF][START_REF] Cocu | Formulation and approximation of quasistatic frictional contact[END_REF]) :

where the bilinear form a(., .) was already defined by (66), j(u, v) = -

r N p v N p , and the following backward difference approximation is used for the velocity at time t k+1 :

Another equivalent incremental formulation is (see [START_REF] Cocu | Existence d'une solution du problème quasi-statique de contact unilatéral avec frottement non local[END_REF]) :

As shown in [START_REF] Cocu | Existence d'une solution du problème quasi-statique de contact unilatéral avec frottement non local[END_REF], two algorithms can be deduced : one is set on displacement variables and the other on displacement increments. The first algorithm is used here. At each step t k , the problem to be solved has the same form as a static frictional contact problem with an extra term which corresponds to the tangential contact displacement u k T at the previous step. It therefore characterizes the loading history. This is natural because we have a path dependent solution.

Various methods have been implemented for solving problem (115). Three of them are based on projection techniques for solving the minimization problem under constraints obtained by using a fixed point method on the friction force limit g = µ|r N | (sequence of Tresca problems) : an Over Relaxation method with Projection, a Projected Conjugate Gradient method with preconditioning or a projected Gauss Seidel algorithm with an Aitken acceleration (see Raous and co-workers [START_REF] Raous | Numerical methods for frictional contact problems and application[END_REF][START_REF] Raous | Preconditioned conjugate gradient method for a unilateral problem with friction[END_REF][START_REF] Chabrand | Comparison of various numerical methods for solving unilateral contact problems with friction[END_REF]). Another one is a direct mathematical programming method used to solve the problem written in linear complementarity form (see [START_REF] Klarbring | Mathematical programming in contact problems[END_REF]) : this is the Lemke method (see [START_REF] Raous | Numerical methods for frictional contact problems and application[END_REF][START_REF] Chabrand | Comparison of various numerical methods for solving unilateral contact problems with friction[END_REF]).

With the examples presented below, the quasistatic solution is computed using the incremental formulation (115) and the S.O.R.P. (Successive Over Relaxation with Projection) method coupled with the fixed point on g. This is written as follows for each time t k+1 :

where r N ((g k+1 ) ν ) is the normal contact force associated to the solution (u k+1 ) ν the problem :

The results of the stability analyses

The results of the stability analyses performed at each step of the quasistatic evolution are presented in Table 1. For each value of the prescribed tangential displacement U T , the total number of stick, slip or no-contact nodes is given in that table together with the minimum real eigenvalues of K * S and M * -1 K * and the lowest mode number corresponding to a complex eigenvector of M * -1 K * . The necessary condition for divergence instability is first satisfied for U T =-1.35mm when a second node begins sliding. The eigenvector associated with this negative eigenvalue of K * S is plotted in Fig. 16. It represents a direction along which an unstable evolution would be energetically admissible : an Therefore we compute different values of α and β (see Clough [START_REF] Clough | Dynamics of structures[END_REF]) so as to obtain the various decreasing evolutions of the damping ratio given on Fig. 18. When referring to each of those curves the value of ξ for 454 Hz will be used. The damping ratio ξ is defined as -Re(λ) |Im(λ)| ( λ is the eigenvalue of the eigenproblem below, with all the contact nodes stuck). Stability has been studied with these 6 different damping levels. The eigenvalues are now computed from :

The results are given on Fig. 19. The values of the imposed tangential displacement U T at which flutter occurs are given in that figure for the six damping levels of Fig. 18. The sudden change in U T occurring between the values of ξ=0.1475 and ξ=0.15 on Fig. 19 is due to the change in the excited mode. For the first value, mode 3 is concerned and for the second one, mode 5 is concerned. That means that when the damping is strong enough, the mode 3 is not excited anymore by flutter.

Therefore, it can be observed on Fig. 19 that the occurrence of instability is greatly delayed by the introduction of viscous damping (from U T =-1.35mm (the no-viscosity case) to U T =-6.2mm (for the maximum damping considered here)).

A more precise estimation of the occurrence of flutter when viscosity is taken into account would require a precise measurement of the loss angle of the material in the appropriate frequency range.

Conclusions

The steady sliding solution of the block is found to be unstable by flutter.

During the tangential quasistatic loading process, we observe that :

• the necessary condition for divergence is satisfied too early in comparison with the occurrence of the first experimental jump ; the corresponding eigenmode suggests the initiation of the waves experimentally observed, but

• the sufficient condition for divergence is never satisfied ;

• growing flutter-type oscillations for the lower order modes initiate here at the same time as the necessary condition for divergence is satisfied, which is too early in comparison with the experiment ;

• consideration of damping in the stability analysis delays the occurrence of flutter, probably leading to a better agreement with the experimental results.