Study of some rheological models with a finite number of degrees of freedom

Eur. J. Mech. A/Solids 19 (2000) 277–307

2000 Éditions scientifiques et médicales Elsevier SAS. All rights reservedS0997-7538(00)00163-7/FLA

Study of some rheological models with a finite number of degrees of freedom

Jérôme Bastiena, Michelle Schatzmanb, Claude-Henri Lamarquea

a URA 1652 CNRS, département génie civil et bâtiment, laboratoire géomatériaux, École Nationale des travaux publics de l’Etat,rue Maurice Audin, 69518 Vaulx-en-Velin cedex, France

b UMR 5585 CNRS, analyse numérique, université Claude Bernard Lyon I, 69622Villeurbanne cedex, France

(Received 6 April 1999; revised and accepted 2 November 1999)

Abstract – A large number of rheological models can be covered by the existence and uniqueness theory for maximal monotone operators. Numericalsimulations display hysteresis cycles when the forcing is periodic. A given shape of hysteresis cycle in an appropriate class of polygonal cycles canalways be realized by adjusting the physical parameters of the rheological model. 2000 Éditions scientifiques et médicales Elsevier SAS

friction / St-Venant element / elastoplasticity / maximal monotone graph

1. Introduction

In the literature, many studies about the behaviour of nonlinear models can be found within the last twentyyears. B. Brogliato (1996) and V. Palmov (1998) give in their respective books numerous nonlinear mechanicalmodels. Among these models, we are interested in those involving friction laws.

Authors especially study mechanical models with a finite number of degrees of freedom involving frictionterms and submitted to dynamical solicitations. Some of these papers provide mathematical results of existenceand uniqueness (Jean and Pratt, 1985; Laghdir and Monteiro Marques, 1995; Matrosov and Finogenko, 1995,1996a, b; Monteiro Marques, 1994; Trinkle et al., 1997). Some others investigate physical behaviours withoutdealing with theoretical results: some works present friction laws issued from experiments (Anderson and Ferri,1990 and Ferri and Bindemann, 1995); stick-slip phenomena is the main topic investigated in (Awrejcewicz andDelfs, 1990a, b; Baumberger et al., 1995; Ionescu and Paumier, 1993; Pratt and Williams, 1981) and (Stelter,1992). The study of friction may be based upon analytical calculation of solutions (Capecchi and Vestroni,1995). Some works describe experiments, identification and modelling of friction (Dowell and Schwartz,1983a, b) and (Tomlinson and Chen, 1996). Numerical experiments have been made by many authors usingclassical numerical schemes (Stewart, 1996; Stewart and Trinkle, 1996, 1997). Chaotic behaviour has beenexhibited by nonlinear models including friction (Chua et al., 1986a, b; Madan, 1993; Popp and Stelter, 1990)and (Shaw, 1986). Except for the paper of Monteiro Marques (1994), these works do not present convergenceresults of the approximate solution.

In older literature, we can read works coping with continuous models involving elastic linear behaviourand friction terms at boundary conditions (for example (Duvaut and Lions, 1972)). In this case, one hasexistence and uniqueness results for the mechanical problems which are expressed via variational inequalities.These inequalities are discretized by numerical schemes with good mathematical properties; one has results ofconvergence of the approximate solution described by the numerical schemes to the exact solution (Glowinskiet al., 1976a, b). Nevertheless these works are related to stationary problems. So they are useless for the studyof our rheological problems.

278 J. Bastien et al.

In order to study micro-plasticity, Fougeres and Sidoroff (1989) investigates the classical elastoplasticMasing model, with a finite or infinite number of degrees of freedom; but, this model is a stationary one.

In this paper, we intend to generalize the study by M.D.P. Monteiro Marques (1994): he examines the modelcomposed of one linear spring and one St-Venant element connected in parallel to a material point; he providesexistence and uniqueness results. This proof is based on the properties of two numerical schemes. This modelcan be used to process the case of ann-degrees-of-freedom model of elements connected in parallel. But hisresults cannot be used in order to deal with connections in series of spring and St-Venant elements. So wedescribe a general mathematical frame which permits us to process many rheological models including springsand St-Venant elements connected either in series or in parallel. Moreover this mathematical frame provideswell adapted numerical schemes and their properties of convergence to the exact solution of the model. In amathematical point of view, linear viscosity and linear elasticity are similar: so, our developments are valid forboth visco-plastic or elasto-visco-plastic models.

This paper is organized as follows: in Section 2, we will describe some models involving springs, dashpotsand St-Venant elements connected in series or in parallel; we give a mathematical theory of the dynamics ofthis rheological models. In Section 3, we can see that the motion of these dynamical systems is governed bydifferential inclusions of the form

X(t)+MA(X(t)) 3 F (t,X(t)), (1.1)

whereM is an invertible matrix,X is a function from[0, T ] in RN , A is a maximal monotone graph onRNandF a function from[0, T ]×RN in RN . The existence and uniqueness of solutions are consequences of veryclassical results. In Section 4, we give a numerical scheme and we report the results of numerical simulations forelastoplastic models. We find in particular that under periodic loading, the solutions tend to a limit cycle whichin the displacement-force plane is a hysteresis cycle. We will prove that there is a one-to-one correspondencebetween the geometrical shapes of these cycles and the mechanical characteristics of the elastoplastic models.

2. The physical models

In this section, we will study combinations of springs, dry friction elements, dashpots and material points.In Section 2.1 we study two elementary models: we recall first the model composed of one spring and one St-Venant element connected in parallel, which has been studied by Monteiro Marques (1994). Then we presentthe association of one spring, one St-Venant element and one dashpot, connected in series. In Section 2.2 wewill consider elementary associations (spring, St-Venant element and material points) and their combinations.

2.1. Two elementary models

2.1.1. One spring and one St-Venant element connected in parallel

We have to keep in mind that the following mechanical system has been previously studied by MonteiroMarques. A material point of massm submitted to an external forceF is connected in parallel to a spring withstiffnessk and a St-Venant element with thresholdα; let x be the abscissa of this material point. This system isgoverned by the equation

mx + kx + ασ(x) 3 F, (2.1)

with initial data

x(0)= x0, x(0)= x0; (2.2)

Study of some rheological models with a finite number of degrees of freedom 279

Figure 1. The graphσ (a) and the graphβ (b).

Figure 2. A viscoelastoplastic model.

the graphσ is defined by (figure 1a):

σ (x)=−1 if x < 0,

1 if x > 0,

[−1,1] if x = 0.

(2.3)

2.1.2. A viscoelastoplastic model and the Prandtl rheological model

We present now a viscoelastoplastic model, described infigure 2: a material point of massm is submitted toan external forceF and is connected in series to a spring with stiffnessk, a St-Venant element with thresholdαand a dashpot with viscosityc. Letu be the displacement of the extremityA of the spring, relative to a referenceposition, letv be the difference betweenAB and a reference length and letw be the difference betweenBCand a reference length. Denotex the abscissa of the material point of massm, relative to a reference position,so that

x = u+ v+w. (2.4)

Denotef the force exerted byA on B. We make the assumption thatk 6= 0, c 6= 0 andα 6= 0. We write theconstitutive law of the spring under the form

f =−ku, (2.5)

the constitutive law of the St-Venant element under the form

f ∈−ασ(v), (2.6)


and the constitutive law of the dashpot under the form

f =−cw. (2.7)

The fundamental theorem of dynamics and (2.5) give:

mx = F − ku. (2.8)

We rewrite (2.4), (2.5), (2.6), (2.7) and (2.8) to better understand the mathematical structure of the problem.Differentiating (2.4) with respect to time, we infer from (2.5), (2.6) and (2.7) the following relation

ku ∈ ασ(x − u− ku/c). (2.9)

We can remark that (2.9) implies

ku ∈ [−α,α]. (2.10)

Let us define the graphβ by (figure 1b):

β(x)=

∅ if x ∈ ]−∞,−1[ ∪ ]1,+∞[,{0} if x ∈ ]−1,1[,R− if x =−1,

R+ if x = 1.

(2.11)

With this notation, the graphβ is the inverse ofσ and (2.9) is equivalent to

u+ β(ku/α) 3 x − ku/c. (2.12)

The initial datax(0), u(0) andx(0) are given. Therefore, the system describing the mechanical set-up offigure 2is (2.8) together with (2.12) and the initial data

x(0)= x0, x(0)= x0, u(0)= u0. (2.13)

Conversely, from the knowledge of (u, x) we can recover the functions (v,w) thanks to (2.4) and (2.7). Bysetting

η= α/k and y = x,we can see that the system (2.8), (2.12) and (2.13) is equivalent to the system

x = y,y = (F − ku)/m,u+ β(u/η) 3 y − ku/c,

(2.14)

with the initial data

x(0)= x0, y(0)= y0, u(0)= u0 ∈ [−η, η]. (2.15)


Figure 3. The Prandtl rheological model.

Figure 4. Two springs and two St-Venant elements.

For the system (2.14) and (2.15), the limit casec → +∞ corresponds to the Prandtl rheological model,governed by equation (2.15) and

x = y,y = (F − ku)/m,u+ β(u/η) 3 y.

(2.16)

This can be explained on physical grounds: if the viscosity of the dashpot is very large, it works as a solidand the model described infigure 2 is equivalent to the association of one spring and one St-Venant element(figure 3). The addition of dashpots does not change the mathematical theory; thus, we study henceforth onlyelastoplastic models, composed of springs and St-Venant elements.

2.2. Other models

Now let us consider models made out of a number of Hooke and St-Venant elements and material points. Ingeneral, the elements that we consider satisfy the following constitutive laws

fi =−kiui for thei-th Hooke element, (2.17)

with ki > 0 and

gi ∈−αiσ (vi) for thei-th St-Venant element, (2.18)

with αi > 0. In order to obtain interesting and new models, we cannot associate arbitrarily these components;moreover, we show the equivalence between certain classes of models, and this reduces the number of cases toconsider. Thus, in Section 2.2.1 we present two ill-posed models: the first one is an undetermined model; forthe second one, we present two physical models which are mathematically equivalent. Then, in Section 2.2.2we present all the models that we study in this paper.

2.2.1. Two ill-posed models

We consider first the connection with two springs (k1, k2) and two St-Venant elements (α1, α2) described infigure 4. From (2.17) and (2.18), we have

f =−k1u1=−k2u2 and f ∈−α1σ (v1), f ∈ −α2σ (v2). (2.19)


Figure 5. The Prandtl model with linear hardening (a) and model with three springs and one St-Venant element (b).

If α1< α2, we have according to (2.19)

|f |6 α1<α2,

and

v2= 0;hence, the St-Venant element 2 is always locked and does not alter the mechanical system. Ifα1 = α2, thissystem is undetermined from a mechanical point of view: thus, if

|f |< α1= α2,

both St-Venant elements are locked; if

|f | = α1= α2,

we cannot determine which St-Venant element will run. Therefore this model is neither mechanically soundnor interesting.

Let us show now an equivalence between the two models given byfigure 5aandfigure 5b.

We can observe that the model offigure 5ais governed by the system (withη= α/k)x = y,y = (F − k0x − ku)/m,u+ β(u/η) 3 y,

(2.20)

with initial data

x(0)= x0, y(0)= y0, u(0)= u0 ∈ [−η, η]. (2.21)

The model offigure 5bis governed by the same system by setting

k0= K2K3

K2+K3, k = K1K3

2

(K1+K2+K3)(K2+K3)and α =A K3

K2+K3.

We remark that, if we chosek0= 0 in (2.20), we find again the system (2.16). The model described infigure 5ais the Prandtl model with linear hardening.

2.2.2. The studied models

This study shows that we may restrict ourselves to two kinds of elementary pairs: thePi-elements (figure 6a)are parallel associations of one spring of stiffnesski > 0 and one St-Venant element of thresholdαi > 0; the


Figure 6.Pi -element (a) andSi -element (b).

constitutive law ofPi is mathematically described by:

fi =−kiui and gi ∈−αiσ (ui). (2.22)

TheSi-elements (figure 6b) are serial associations of one spring of stiffnesski > 0 and one St-Venant elementof thresholdαi > 0; the constitutive law ofSi is mathematically described by:

fi =−kiui ∈−αiσ (vi). (2.23)

We assume that there is no secular term of deformation when the external forcefi is vanishing: according tothe model described in Section 2.1.2 we have chosen a reference position of the material point and a referenceposition of the St-Venant element such thatfi is vanishing whenui is equal to zero.

We present now some models withPi -elements,Si-elements and several material points. According to theforegoing study of the examples in Section 2.2.1 we consider particular combinations ofPi -elements orSi-elements and material points. We will give successively:

• the model withn Pi -elements connected in parallel with one material point (see Section 2.2.3),• the model withn Si-elements connected in parallel with one material point (see Section 2.2.4),• the model withn Si-elements andn material points connected in series (see Section 2.2.5),• the model withn Pi -elements andn material points connected in series (see Section 2.2.6),• the model withn Pi -elements, one spring and one material point connected in series (see Section 2.2.7),• and the mixed models with oneS1-element, oneP2-element and material points (see Section 2.2.8).

2.2.3. Parallel association ofn Pi-elements with one material point

We consider the association ofn Pi-elements, fori = 1, . . . , n, connected in parallel (figure 7). The materialpoint has abscissax and massm. This material point is submitted to an external forceF .

We have the constitutive laws (2.22) fori = 1, . . . , n. As previously, we have the differential system

mx +(

n∑i=1

ki

)x +

(n∑i=1

αi

)σ (x) 3 F, (2.24)

which is the inclusion (2.1) in Section 2.1. Thus, we have the equivalence betweenn Pi -elements connected inparallel and onePi -element.

2.2.4. Parallel association ofn Si-elements and one material point

This system is described infigure 8a; the material point has abscissax and massm.


Figure 7. Parallel association ofn Pi -elements and one material point.

Figure 8. The generalized Prandtl rheological model (a) and the generalized Prandtl rheological model with linear hardening (b).

This mechanical system is described by the differential inclusion

{mx =−∑n

i=1 kiui + F,∀i ∈ {1, . . . , n}, ui + β(kiui/αi) 3 x. (2.25)


Figure 9. Association in series ofn Si -elements andn material points.

Setting

y = x and ∀i ∈ {1, . . . , n}, ηi = αi/ki, (2.26)

we write equation (2.25) under the form

x = y,y = F/m− (1/m)∑n

i=1 kiui,

∀i ∈ {1, . . . , n}, ui + β(ui/ηi) 3 y;(2.27)

the initial data are

x(0)= x0, y(0)= y0, ∀i ∈ {1, . . . , n}, ui(0)= ui,0 ∈ [−ηi, ηi]. (2.28)

This system is a generalized Prandtl rheological model.

Setting one of theαis equal to infinity amounts to saying that one of theSi-elements is a pure spring; thisstrengthens the model. We obtain the equation

x = y,y = F/m− (k0/m)x − (1/m)∑n

i=1 kiui,

∀i ∈ {1, . . . , n}, ui + β(ui/ηi) 3 y.(2.29)

The relevant mechanical system is described infigure 8b(with n Si-elements and one spring); it is a generalizedPrandtl rheological model with linear hardening. This model is also called the discrete Masing model.

2.2.5. Association in series ofn Si-elements andn material points

This system is described infigure 9; for i = 1, . . . , n, xi denotes the abscissa of material point of massmi .

By setting

∀i ∈ {1, . . . , n}, yi = xi and ηi = αi/ki,


Figure 10.Association in series ofn Pi -elements andn material points.

we obtain the system

x1= y1,

y1= (−k1u1+ k2u2)/m1,

u1+ β(u1/η1) 3 y1,

∀i ∈ {2, . . . , n− 1},xi = yi,yi = (−kiui + ki+1ui+1)/mi,

ui + β(ui/ηi) 3 yi − yi−1,xn = yn,yn = (−knun + F)/mn,un + β(un/ηn) 3 yn − yn−1,

(2.30)


∀i ∈ {1, . . . , n} xi(0)= xi,0, yi(0)= yi,0, ui(0)= ui,0 ∈ [−ηi, ηi]. (2.31)

2.2.6. Association in series ofn Pi -elements andn material points

This system is described infigure 10; for i = 1, . . . , n, xi denotes the abscissa of the material point of massmi .

We have the system {X= Y,Y +D−1A(Y ) 3D−1H −D−1KX,

(2.32)


X(0)=X0 and Y (0)= Y0. (2.33)

We have set

X = (x1, . . . , xn), Y = (y1, . . . , yn), D = diag(m1, . . . ,mn),

H =

0...

0F

and K =

k1+ k2 −k2 0 · · · 0 0−k2 k2+ k3 −k3 · · · 0 0

. . .. . .

. . .

0 · · · −ki ki + ki+1 −ki+1 · · ·. . .

. . .. . .

0 · · · 0 −kn kn

.


Figure 11.Association in series ofn Pi -elements, one spring and one material point.

The maximal monotone operatorA is defined onRn by, for all Y = (y1, . . . , yn) ∈ Rp

A(Y )= {(gi − gi+1)16i6n: g1 ∈ α1σ (y1), ∀i ∈ {2, . . . , n}, gi ∈ αiσ (yi − yi−1), gn+1= 0}.

2.2.7. Association in series ofn Pi -elements, one spring and one material point

The model described infigure 11is a discrete Masing model constituted by the combination in series ofn Pielements and one spring. In the static setting, it is equivalent to the parallel combination ofn Pi -elements andone spring, described infigure 8b(Fougeres and Sidoroff, 1989). We shall not consider here the problem of thedynamical equivalence of these two kinds of models.

We assumek0 > 0. We denote byx the abscissa of the material point of massm. By choosing appropriateauxiliary functionsv1, . . . , vn, we can prove that this mechanical system is governed by

x = y,y = F/m+ γ (∑n

i=1 vi/ki − x),

V +KB(V ) 3−k0yW,

(2.34)

with initial data

x(0)= x0, y(0)= y0 and ∀i ∈ {1, . . . , n}, vi(0)= vi,0 ∈ [−αi, αi]. (2.35)

The numberγ is defined by

γ =(m

n∑i=0

1

ki

)−1

.

The vectorsV andW are given by

V = (v1, . . . , vn) and W = (1, . . . ,1).The coefficients of the matrixK are

∀i, j ∈ {1, . . . , n}, Ki,j ={k0+ ki if i = j,k0 if i 6= j ;

this matrix is symmetric, positive and definite. The maximal monotone operatorB is defined onRn by

∀X= (x1, . . . , xn) ∈Rn, B(X)= β(x1/α1)× β(x2/α2)× · · · × β(xn/αn).


Figure 12.OneS1-element, oneP2-element and two material points connected in series.

2.2.8. Two other associations

The following two kinds of systems are not redundant with any of the previous models. Their equations ofmotions simultaneously contain the graphsσ andβ. This system is described infigure 12; x1 andx2 denote theabscissa of material points of massm1 andm2.

Let us sety1= x1, y2= x2 andη1= α1/k1, we obtain

x1= y1,

y1− (1/m1) g2= (− z1+ k2(x2− x1))/m1,

x2= y2,

y2+ (1/m2) g2= (F − k2(x2− x1))/m2,

g2 ∈ α2σ (y2− y1),

u1+ β(u1/η1) 3 y1,

(2.36)

with initial data

u1(0)= u1,0 ∈ [−η1, η1] and for i = 1,2 xi(0)= xi,0, yi(0)= yi,0. (2.37)

When we choose theS1-element and theP2-element connected in parallel, then we obtain the system (withη1= α1/k1)

x = y,y + (α2/m)σ (y) 3 (F − z− x)/m,u+ β(u/η1) 3 y,

(2.38)


x(0)= x0, y(0)= y0 and z(0)= z0 ∈ [−η1, η1], (2.39)

wherex denotes the abscissa of the material point of massm.

3. Existence and uniqueness results

We observe that all the models introduced in Section 2, except the one infigure 4, can be subsumed underone form which is conveniently described in the language of maximal monotone operators. In Section 3.1 werecall some properties of maximal monotone operators. Then we give existence and uniqueness results for thesedifferential inclusions in Section 3.2.


3.1. Reminder about maximal monotone graphsσ andβ

The reader is referred to (Brézis, 1973). Let〈 , 〉 be the scalar product onRp. If φ is a convex proper andlower semi-continuous function fromRp to ]−∞,+∞], we can define its sub-differential∂φ by{

y ∈ ∂φ(x)⇐⇒∀h ∈Rp, φ(x + h)− φ(x)> 〈y,h〉,D(∂φ)= {x: ∂φ(x) 6= ∅}; (3.1)

moreover,∂φ is a maximal monotone graph inRp ×Rp. The maximal monotone graphsσ andβ of Section 2are sub-differentials of proper semi-continuous convex functions|x| andψ[−1,1] defined by

∀x ∈R, ψ[−1,1](x)={

0 if x ∈ [−1,1],+∞ if x 6∈ [−1,1], (3.2)

and the choice of the canonical scalar product inR. Therefore:

∀x ∈R, σ (x)= ∂|x| and β(x)= ∂ψ[−1,1](x). (3.3)

We observe that ifRp is equipped with its canonical scalar product, and with another scalar product

〈x, y〉M = xTM−1y, (3.4)

whereM is the symmetric positive definite, then we can relate the sub-differential∂φ of φ relative to thecanonical scalar product and the sub-differential∂Mφ relative to〈 , 〉M by

∂Mφ(x)=M∂φ(x). (3.5)

3.2. Mathematical study of differential systems

We give now the general mathematical formulation of our problem. We assume thatT is strictly positive andthatG is a function from[0, T ] ×Rp toRp which is Lipschitz continuous with respect to its second argument,i.e. there existsω > 0 such that

∀t ∈ [0, T ], ∀X1,X2 ∈Rp,∥∥G(t,X1)−G(t,X2)

∥∥6 ω‖X1−X2‖. (3.6)

Moreover, we assume that

∀Y ∈Rp, G(·, Y ) ∈L∞(0, T ,Rp). (3.7)

The matrixM is symmetric positive definite, andφ is convex proper and lower semicontinuous onRp.

PROPOSITION 3.1: Under the previous assumptions, for allξ ∈D(∂φ), there exists a unique function X inW 1,1(0, T ,Rp) such that {

X(t)+M∂φ(X(t)) 3G(t,X(t)) a.e. on]0, T [,X(0)= ξ. (3.8)

Proof. –We observed above that∂Mφ(x)=M∂φ(x); therefore Proposition 3.13, page 107 of (Brézis, 1973)applies.


Table I. The dimension of the system, the convex function and the invertible matrix used for the above described mechanicalmodels.

System p Functionφ Matrix M

(2.1)–(2.2) 2 φ(x,y)= α|y| diag(1,1/m)

(2.14)–(2.15),(2.15)–(2.16),(2.20)–(2.21)

3 φ(x,y,u)=ψ[−η,η](u) I

(2.27)–(2.28),(2.28)–(2.29)

n+ 2 φ(x,y,u1, . . . , un)=n∑i=1

ψ[−ηi ,ηi ](ui ) I

(2.34)–(2.35) n+ 2 φ(x,y, z1, . . . , zn)=n∑i=1

ψ[−αi ,αi ](zi)

(1 0 00 1 00 0 K

),

K not diagonal

(2.30)–(2.31) 3n φ(x1, y1, u1, . . . , xn, yn,un)=n∑i=1

ψ[−ηi ,ηi ](ui) I

(2.32)–(2.33) 2n φ(x1, . . . , xn, y1, . . . , yn)= α1|y1| +n∑i=2

αi |yi − yi−1|(I 00 D−1

),

D diagonal

(2.36)–(2.37) 5 φ(x1, y1, x2, y2, u1)= α2|y2− y1| +ψ[−η1,η1](u1) diag(1,1/m1,1,1/m2,1)

(2.38)–(2.39) 3 φ(x,y,u)= α2|y| +ψ[−η1,η1](u) diag(1,1/m,1)

Thus, all the systems of Section 2 can be written under the form (3.8) and have a unique solution. For allsystems,table I provides the corresponding integerp, functionφ and matrixM . It is easy to prove thatφ areconvex proper and lower semi-continuous function onRp. With this table, we can observe that there are threeclasses of mechanical systems:

• in the first class, the functionφ is a linear combination of|xi|;• in the second class, the functionφ is a sum of functionsψ[−ηi,ηi ];• in the third class,φ involves the two functions|x| andψ[−ηi,ηi ].

4. Numerical simulations

In Section 4.1 we present some numerical simulations for the Prandtl model, the Prandtl model withlinear hardening and the generalized Prandtl model with linear hardening in order to show that they candescribe materials with elastoplastic constitutive laws. In Section 4.2 we give several simulations of theviscoelastoplastic model. All the curves obtained are presented in Appendix A. We have also numericallysimulated other models, governed by systems (2.1) and (2.2), (2.30) and (2.31), (2.38) and (2.39); the numericalmethods are similar to those used for the Prandtl models and we obtain hysteresis cycles, but we were not ableto infer from them relevant information.

For all these numerical simulations, we solve the differential inclusion (3.8) by using the implicit Eulerscheme: letn ∈N∗; let h= T /n and forq ∈ {0, . . . , n}, let tq = qh. We solve{

∀q ∈ {0, . . . , n− 1}, Xq+1−Xqh+M∂φ(Xq+1) 3G(tq,Xq),

X0= ξ.(4.1)


The first equality of (4.1) is equivalent to

∀q ∈ {0, . . . , n− 1}, Xq+1= (I + hM∂φ)−1(hG(tq,Xq)+Xq). (4.2)

We denote byXh the linear interpolation of theXn’s. The functionXh converges to the solutionX of system(3.8) inC0([0, T ],Rp). This result is proved in (Crandall and Evans, 1975) which contains a much more generalresult.

4.1. Study of the rheological Prandtl models

In this section, we study the rheological Prandtl models, which are shown in Section 2, i.e. the Prandtl model(see system (2.15) and (2.16)), the Prandtl model with linear hardening (see system (2.20) and (2.21)) and thegeneralized Prandtl model with linear hardening (see system (2.28) and (2.29)). We choose

m= 1. (4.3)

We will study a harmonic forcing (Section 4.1.1) and other periodic forcings (Section 4.1.2).

4.1.1. Harmonic forcing

For this section, we choose

F(t)= f cos(ωt). (4.4)

We first present essential differences, between the Prandtl model and the Prandtl model with linear hardening:we have computed the functionsx, y and z for the Prandtl model on the interval[0, T ] with the followingparameters

T = 2000, F (t)= 20cos(0.1t), η= k = 1 and x0= y0= u0= 0, (4.5)

and for the Prandtl model with linear hardening on the interval[0, T ] with the following parameters

T = 300, F (t)= 20cos(0.1t), k0= 1, η= k = 1 and x0= y0= u0= 0. (4.6)

Theses functions are infigures 16and17.

We can see that the amplitude of the functionsx andx are larger for the Prandtl model, which correspondsto the system (2.20) withk0 = 0 than for the Prandtl model with hardening, which corresponds to the samesystem withk0 6= 0. We sketch a physical explanation: if the St-Venant element slips, then we have from (2.20)

u≡ εη,and

x(t)+ k0

mx(t)= f

mcos(ωt)− εkη

m, (4.7)

whereε ∈ {−1,1}. If we assumek0= 0, and if the slip phase starts att0, then (4.7) gives

x(t)= 2f

mω2sin(ω(t + t0))sin

(ω(t − t0))− εkη

2m(t − t0)2+

(x(t0)− f

mωsin(ωt0)

)(t − t0)+ x(t0). (4.8)


Forω� 1 andt ≈ t0, we have

x(t)≈ 2f

mω2sin(ω(t − t0)); (4.9)

thus, the approximate amplitude ofx is roughly

2f

mω2= 4000, (4.10)

as we can see infigure 16. If k0= 1 and if the slip phase starts att0, then (4.7) implies

x(t)=D cos(ω0t + φ)+ f

m(ω20− ω2)

cos(ωt)− εkη

mω20, (4.11)

whereω0 = √k0/m 6= ω and (D,φ) depends on(x(t0), x(t0)). Here, withω = 0.1, the amplitude of thesinusoidal component of periodω of the functionx is roughly

f

m(ω20−ω2)

≈ 20, (4.12)

as we can see infigure 17.

We will study rheological properties of Prandtl models. In all cases, the material point is submitted to theexternal forceF and has an abscissax. The choice of representation is crucial; we have plotted infigure 18athe curve{x(t),F (t)}t∈[3500,8000] for the Prandtl model defined by the following parameters

F(t)= 20cos(0.1t), η= k = 1, and x0= y0= u0= 0.

We observe a limit cycle, but it gives no information on the physically relevant parameters of the system; ifhowever we plot, as infigure 18b, the curve{x(t),F (t)−mx(t)}t∈[3500,8000], we observe that the slope of theoblique parts of the cycle is equal tok.

In figure 19, we can observe transients (t ∈ [0,1500]) followed by a periodic regime (t > 1500), composedof periodic hysteresis cycles. We plot only the periodic regime. Infigure 19, we can see thatx(t) is vanishingwhenF(t) − mx(t) is equal to zero (fort = 0), according to the choice of Section 2.2.2. Nevertheless, infigure 20, we can see that the cycle does not contain the point(0,0) in the plane(x,F −mx): we have notplotted transients in this figure.

In table II we define the parameters of the simulations for the generalized Prandtl model with linearhardening, which includes as particular cases, the Prandtl model and the Prandtl model with hardening: thecorresponding system of equations is (2.28) and (2.29).

In figures 20–23, we observe hysteresis cycles in the(x,F −mx) plane. Infigure 20, we can see identicalhysteresis cycles up to translation for different initial data. Infigure 21a(n= 5 andf = 15), we can see thatthe hysteresis cycle is composed of 12 line segments. But, infigure 21b(n= 5 andf = 10), we can see only6 line segments. Infigures 22and23a (large values ofn), we can distinguish line segments when magnifyingthe figure (figure 23).

We have also made simulations with other values of parameters: we obtain hysteresis limit cycles; thesecycles are periodic but they are not convex.


Table II. The parameters of the numerical simulations of (2.28) and (2.29) with forcing (4.4).

No. of figure f ω a b n k0 x0 y0 ηi ki u0,i

20a 1 0.5 50 2000 1 0 0 0 1 1 0

20b 1 0.5 50 2000 1 0 100 −200 1 1 0.99

21a 15 0.5 400 3000 5 0 0 0 i 1 0

21b 10 0.5 400 3000 5 0 0 0 i 1 0

22a 12 0.5 500 4000 30 0 0 0 i 1/n 0

22b 12 0.5 500 4000 30 1/n 0 0 i 1/n 0

23 55 0.5 500 4000 100 0 0 0 i 1/n 0

Figure 13. The functionsH0 andH1.

Table III. The parameters of the numerical simulations of (2.28) and (2.29) with forcing (4.13).

No. of figure H f ω a b n k0 x0 y0 ηi ki u0,i

24 H0 15 0.5 5 2000 1 0 0 0 1 1 0

25 H1 80 0.5 400 500 1 0 0 0 1 1 0

26 H1 80 10 400 500 1 0 0 0 1 1 0

27a H1 80 0.5 400 500 1 1 0 0 1 1 0

27b H1 80 10 400 500 1 1 0 0 1 1 0

28 H1 800 0.5 960 2000 30 1/n 0 0 i 1/n 0

4.1.2. Periodic forcing

We choose now a forcing of pulsationω and amplitudef

∀t ∈ [0, T ], F (t)= fH(ωt), (4.13)

whereH is 1-periodic and of amplitude 1. We performed simulations using the functionsH0 andH1 graphedin figure 13.

These simulations are similar to the simulations of Section 4.1.1. The parameters are defined intable III.

In figure 24, we find again the same shape of hysteresis cycle as infigure 20a. Thus, we can conclude thatwe obtain hysteresis cycles with a continuous and periodic forcing. On the contrary, we do not necessarilyobtain limit periodic hysteresis cycles with a discontinuous forcing; infigure 25, we observe a resonancewhose amplitude decreases when the pulsation increases (figure 26). For the static models this resonancephenomenon is known as the ratchet phenomenon. When we choosek0 6= 0, we obtain limit periodic hysteresis


Figure 14.Loading curve for the generalized Prandtl rheological model with linear hardening.

cycles (figures 27and28). In figure 28, we obtain periodic hysteresis cycles similar to the hysteresis cycle offigure 22b.

4.1.3. Conclusions on elastoplastic models and analysis of hyteresis cycles

Most of the responses of the generalized Prandtl model to a sinusoidal forcing (figures 20–23) or to a periodic(but non-harmonic) forcing (seefigures 24, 27and28), tend to hysteresis limit cycles, which we are going tostudy.

We observe that these hysteresis cycles have a centre of symmetry. Therefore, we study a loading phasecorresponding to a half-cycle in the(x,F −mx) plan (figure 14). Let us denote byA1,A2, . . . ,An+1 andAn+2

the ends of segments which constitute the hysteresis half-cycle. Fori ∈ {1, . . . , n}, di denotes the differencebetween abscissa ofAi+1 andA1 (positive real). Fori ∈ {1, . . . , n+ 1}, pi denotes the slope of the segment[Ai,Ai+1].

If the forcingF is not constant on any open non-empty subinterval of[0, T ], thenx has the same property;assume indeed thatx is constant on(t1, t2); thenx andx vanish on(t1, t2) and by uniqueness of the solution of{

ui + β(ui/ηi) 3 0,

ui(t1) given,(4.14)

ui(t) is equal toui(t1) on (t1, t2). This implies that the expression

F =mx + k0x +n∑i=1

kiui

is constant over(t1, t2), which contradicts our assumption.


We are now able to describe the shape of the representation of a trajectory of the system in thex andF −mxcoordinates, in some special case. We assume that theηi = αi/ki are all distinct; otherwise two elements withidentical ratioηi would enter the plasticity phase simultaneously. We reorder the indices so that

η1< η2< · · ·< ηn−1 < ηn. (4.15)

Assume thatF vanishes on no open non-empty subinterval of[0, T ], and the following properties:

(1) uj(0)=−ηj , ∀j = 1, . . . , n,(2) x is increasing on[0, T ].

Then there exists an increasing sequencet16 t26 · · ·6 tn such that for somej ∈ {1, . . . , n}t1< t2< · · ·< tj = tj+1= · · · = T (4.16)

and we have moreover

x(tj )− x(0)= 2ηj , (4.17)

F −mx −(k0+

n∑i=j+1

ki

)x is constant on]tj , tj+1]. (4.18)

Observe first thatx is strictly increasing on[0, T ], thanks to our assumption onF . Let

uj (t)=−ηj + x(t)− x(0), (4.19)

and denote byξ the inverse function ofx and define

tj ={ξ(2ηj + x(0)) if 2ηj + x(0) < x(T ),T if 2ηj + x(0)> x(T ). (4.20)

We check immediately that the function

uj (t)={uj (t) if 0 6 t < tj ,ηj if t > tj ,

(4.21)

is the (unique) solution of {uj + β(uj/ηj ) 3 x,uj (0)=−ηj . (4.22)

Therefore, on]tj , tj+1], relation (4.18) holds. Moreover, relation (4.20) is equivalent to (4.17). We writeequations (4.18) and (4.17) under the form

∀j ∈ {1, . . . , n+ 1}, pj = k0+n∑l=jkl, (4.23)

∀i ∈ {1, . . . , n}, di = 2ηi. (4.24)

On the other hand, we remark that the assumptionsk0> 0, ki > 0 and equation (4.23) imply

pn+1< pn < · · ·< p2< p1; (4.25)


Table IV. Constantsηi , ki and constantsηci

and kci

computed from the cycle of thefigure 21a.

i 0 1 2 3 4 5

ηci

− 0.999921 1.999945 3.000250 3.999900 5.000013

ηi − 1 2 3 4 5

kci

3.2210−14 0.999892 0.999959 1.000082 0.999971 0.999997

ki 0 1 1 1 1 1

Table V. Constantsηi andki computed from thecycle of thefigure 15.

i 0 1 2 3 4 5

ηi − 1/2 1 2 7/2 11/2

ki 0 1 2 2/3 1/12 1/4

equations (4.15) and (4.24) imply

d1< d2< · · ·< dn−1 < dn. (4.26)

From (4.23) and (4.24) we have a one to one correspondence between the parameters of generalized Prandtlmodelki andηi and geometrical parameterspj anddj of the hysteresis cycle in the(x,F −mx) plane.

Thus, a partial identification of the model is possible. But we cannot infer from the values of(pi)16i6n+1 and(di)16i6n the abscissa and the ordinate of the lower left point of the cycle: they depend on the initial conditionsand on the forcing: see for examplefigure 20: the two cycles are identical up to translation.

In figure 20, 21a, 22a, 22band23, we can see(n+ 1) different phases in one half-cycle. On the contrary, infigure 21bwe have just three different phases: the amplitude of forcing is too weak and the five St-Venantelements do not enter the slip phase; two of them are always locked and the generalized Prandtl modelworks as a model with threeSi-elements, defined by(k1, η1), (k2, η2) and(k3, η3) and one spring of stiffnessK = k0 + k4 + k5. Let us consider a given hysteresis cycle: we can only calculate the parameters of Prandtlelements which enter slip phase.

In the casek0 = 0, we can prove that if the amplitude of the forcing is large enough for all the St-Venantelements to be unlocked, the ordinate of the symmetry centre of the cycle vanishes. Moreover, the slopepn+1

of the last segment must be zero.

4.1.4. Numerical identification of a model from its limit cycles

If we start from the cycles already computed, we find the physical constants of the systems underconsideration. For instance, forfigure 21a, we identify five St-Venant elements. We give intable IV theconstantsηci andkci that have been obtained and the constantsηi andki that characterize the studied Prandtlmodel (table II).

We find that

max16i65

∣∣ηci − ηi ∣∣= 2.50110−4 and max06i65

∣∣kci − ki ∣∣= 1.07410−4. (4.27)

Conversely, let a polygonal convex half-cycle plotted infigure 15. Formulae (4.23) and (4.24) determine theconstantηi anddi of a generalized Prandtl model; they are given intable V.


Figure 15.Upper part of an arbitrary limit cycle.

With these values ofηi andki andn= 5, we make numerical simulations for the corresponding generalizedPrandtl model. The other parameters are defined arbitrarily by

F(t)= f cos(0.5t), m= 1, x0= x0= 0 and ∀i ∈ {1, . . . ,5} u0,i = 0. (4.28)

As is natural, the amplitude of the forcing has to be large enough to retrieve the cycle offigure 15. A forcingamplitude larger than 6.6 increases the size of the last segment of the cycle.

Forf = 6.6, we obtain a symmetrical hysteresis cycle with 12 segments (figure 29). We find

max16i65

∣∣ηci − ηi ∣∣= 5.45010−5 and max06i65

∣∣kci − ki ∣∣= 1.04210−4. (4.29)

4.2. A viscoelastoplastic model

Let us show some simulations for the system (2.14) and (2.15). In this section, we choose

F(t)= cos(0.5t), (4.30)

and the parameters have the following values

η= k =m= 1 and x0= y0= u0= 0 and c= 100 or c = 1. (4.31)

We have plotted the curves{x(t),F (t) − mx(t)}t∈[400,2000] for c = 100 in figure 30aand for c = 1 infigure 30b. In every case, we obtain hysteresis cycles. For the largest value ofc (figure 30a), the cycle issimilar to those obtained for Prandtl model with harmonic forcing (figure 20a); the behaviour is similar to theelastoplastic case. For smaller values ofc, the cycle no longer looks like an elastoplastic cycle (figure 30b). Thelimit cycle obtained does not enable us to identify the parameters of forcing.

5. Conclusion

We examine several classes of rheological models consisting of springs, St-Venant elements and dashpots.All these models are governed by the differential inclusion (3.8). Note how equation (3.8) is interesting: itpermits us to unify a large number of rheological models with a finite number of degrees of freedom. In amathematical point of view, we proved by using classical results that a unique solutionX exists. Moreover,


we used a numerical scheme which is simple and easy to realize. It avoids distinguishing between the differentphases (friction, movement, stick and slip) which could be very fastidious for systems with many degrees offreedom. Also, convergence of the approximate solution toward the exact solution, can be proved. In a physicalpoint of view, equation (3.8) can be written as{

X(t)+M∂φ(X(t))+C.X(t) 3H(t) a.e. on]0, T [,X(0)= ξ,

whereC is a matrix. This differential inclusion points out the movement of a material point with coordinatesX(t) in Rp (X(t) are generalized coordinates of the studied rheological model): this material point issubmitted to external forceH(t), to repelling forceC.X(t) (which denotes the linear part of the model,consisting of springs) and to convex potentialφ (which denotes the plastic part of the model, the St-Venantelements).

Numerous simulations have been presented for different periodic forcing. All of them exhibit a hysteresislimit cycle that occurs after a transient phase.

In the case of the generalized Prandtl model, the hysteresis cycle obtained determines the physical parametersof the model.

Appendix

Appendix A. Curves of the numerical simulations

Figure 16.The functionsx, y andu for the Prandtl model defined byT = 2000,F(t)= 20cos(0.1t), η= k = 1 andx0= y0= u0= 0.


Figure 17. The functionsx, y andu for the Prandtl model with linear hardening, defined byT = 300,F(t) = 20cos(0.1t), k0 = 1, η = k = 1 andx0 = y0 = u0= 0.

Figure 18. The curves{x(t),F (t)}t∈[3500,8000] (a) and{x(t),F (t)−mx(t)}t∈[3500,8000] (b) for the Prandtl model defined byF(t) = 20cos(0.1t),η= k = 1, andx0= y0= u0= 0. These two figures differ only by the choice of coordinates.


Figure 19.The curve{x(t),F (t)−mx(t), t}t∈[0,1800] for the Prandtl model, defined byF(t)= 20cos(0.1t), η= k = 1 andx0 = y0 = u0= 0.

Figure 20. The curve{x(t),F (t)−mx(t)}t∈[50,2000] for the generalized Prandtl rheological model with linear hardening, defined byF(t)= cos(0.5t),n= 1,k0= 0,η1= k1= 1 and(x0= y0= u0,1= 0) (a),(x0= 100,y0 =−200,u0,1= 0.99) (b). These two figures differ only by the initial conditions.


Figure 21. The curve{x(t),F (t) − mx(t)}t∈[400,3000] for the generalized Prandtl rheological model with linear hardening, defined byF(t) =f cos(0.5t), k0 = 0, n = 5, x0 = y0 = 0, ∀i ∈ {1, . . . ,5}, ki = 1, ηi = i, u0,i = 0 andf = 15 (a) andf = 10 (b). These two figures differ only by

the amplitude of the forcing.

Figure 22. The curve{x(t),F (t) − mx(t)}t∈[500,4000] for the generalized Prandtl rheological model with linear hardening, defined byF(t) =12cos(0.5t), n = 30, x0 = y0 = 0, and∀i ∈ {1, . . . ,30}, ki = 1/30, ηi = i, u0,i = 0 andk0 = 0 (a) k0 = 1/30 (b). These two figures differ only

by the value ofk0.


Figure 23. The curve{x(t),F (t) − mx(t)}t∈[500,4000] for the generalized Prandtl rheological model with linear hardening, defined byF(t) =55cos(0.5t), n = 100, k0 = 0, x0 = y0 = 0, and∀i ∈ {1, . . . ,100}, ki = 1/100, ηi = i, u0,i = 0. Figure (b) is the magnifying and dilating square

of figure (a).

Figure 24. The curve {x(t),F (t) − mx(t)}t∈[50,2000] for the generalized Prandtl rheological model with linear hardening, defined byn = 1,F(t) = 15H0(0.5t), k0 = 0, x0 = y0 = u0,10, k1 = η1 = 1. This curve can be compared to the curve offigure 20awhich corresponds to the same

values of parameters and to a harmonic forcing.


Figure 25. The curve{x(t),F (t) − mx(t), t}t∈[450,500] for the generalized Prandtl rheological model with linear hardening, defined byF(t) =80H1(ωt), n= 1, k0= 0, η1= k1= 1, x0 = y0 = u0,1= 0, andω= 0.5. This curve shows the ratchet phenomenon.

Figure 26. The curve{x(t),F (t) − mx(t), t}t∈[450,500] for the generalized Prandtl rheological model with linear hardening, defined byF(t) =80H1(ωt), n= 1, k0= 0, η1= k1= 1, x0 = y0= u0,1= 0, andω= 10. (b) is the magnifying square of figure (a).


Figure 27. The curve{x(t),F (t)−mx(t)}t∈[450,500] for the generalized Prandtl rheological model with linear hardening, defined byF(t)= 80H1(ωt),n= 1, k0= 1, η1= k1= 1, x0= y0= u0,1= 0, andω= 0.5 (a) andω= 10 (b). These two figures differ only by the pulsation of the forcing.

Figure 28. The curve{x(t),F (t) − mx(t)}t∈[960,1000] for the generalized Prandtl rheological model with linear hardening, defined byF(t) =800H1(0.5t), n= 30,k0= 1/n, x0= y0= 0, and∀i ∈ {1, . . . ,30}, ki = 1/30,ηi = i, u0,i = 0.


Figure 29. The curve{x(t),F (t) − mx(t)}t∈[750,1000] for the generalized Prandtl rheological model with linear hardening, defined bytable VandF(t)= f cos(0.5t), n= 5,m= 1, x0= y0= 0, ∀i ∈ {1, . . . ,5} u0,i = 0 andf = 6.6.

Figure 30. The curves{x(t),F (t)−mx(t)}t∈[400,2000] for the viscoelastoplastic model defined byF(t)= cos(0.5t), η= k =m= 1,x0= y0= u0= 0andc= 100 (a) andc= 1 (b). These two figures differ only by value ofc.


References

Anderson J., Ferri A, 1990. Behaviour of a single degree of freedom system with a generalized friction law. J. Sound Vibration.Awrejcewicz J., Delfs J., 1990a. Dynamics of a self-excited stick-slip oscillator with two degrees of freedom, part I. Investigation of equilibria.

European J. Mech. A/Solids 9 (4), 269–282.Awrejcewicz J., Delfs J., 1990b. Dynamics of a self-excited stick-slip oscillator with two degrees of freedom, part II, slip-stick, slip-slip, stick-slip

transitions, periodic and chaotic orbits. European J. Mech. A/Solids 9 (5), 397–418.Baumberger T., Caroli C., Perrin B., Ronsin, O., 1995. Nonlinear analysis of the stick-slip bifurcation in the creep-controlled regime of dry friction.

Physical Review E 51 (5), 4005–4010.Brézis, H., 1973. Opérateurs maximaux monotones et semi-groupes de contractions dans les espaces de Hilbert. North-Holland Mathematics Studies,

No. 5. Notas de Matemática (50), North-Holland Publishing Co., Amsterdam, 1973.Brogliato B., 1996. Nonsmooth impact mechanics. Springer-Verlag London Ltd., London, 1996. Models, dynamics and control.Capecchi D., Vestroni F., 1995. Asymptotic response of a two degrees of freedom elastoplastic system under harmonic excitation; internal resonance

case. Nonlinear Dynamics 7, 317–333.Chua L.O., Komuro M, Matsumoto T., 1986a. The double scroll family. I. Rigorous proof of chaos. IEEE Trans. Circuits and Systems 33 (11),

1072–1097.Chua L.O., Komuro M, Matsumoto T., 1986b. The double scroll family. II. Rigorous analysis of bifurcation phenomena. IEEE Trans. Circuits and

Systems 33 (11), 1097–1118.Crandall M.G., Evans, L.C., 1975. On the relation of the operator∂/∂s + ∂/∂τ to evolution governed by accretive operators. Israel J. Math. 21 (4),

261–278.Dowell E., Schwartz H., 1983a. Forced response of a cantilever beam with a dry friction damper attached, part I: Theory. J. Sound Vibration 91 (2),

255–267.Dowell E., Schwartz H., 1983b. Forced response of a cantilever beam with a dry friction damper attached, part II: Experiment. J. Sound Vibration

91 (2), 269–291.Duvaut G., Lions, J.-L., 1972. Les inéquations en mécanique et en physique. Travaux et Recherches Mathématiques, No. 21, Dunod, Paris.Ferri A., Bindemann A., 1995. Large amplitude vibration of a beam restrained by a non-linear sleeve joint. J. Sound Vibration 184 (1), 19–34.Fougeres R., Sidoroff F., 1989. The evolutive masing model and its application to cyclic plasticity and ageing. Nuclear Engineering and Design 114,

273–284.Glowinski R., Lions J.-L., Trémolières R., 1976a. Analyse numérique des inéquations variationnelles. Tome 1. Théorie générale premières

applications. Méthodes Mathématiques de l’Informatique, 5, Dunod, Paris.Glowinski R., Lions J.-L., Trémolières R., 1976b. Analyse numérique des inéquations variationnelles. Tome 2. Applications aux phénomènes

stationnaires et d’évolution. Méthodes Mathématiques de l’Informatique, 5, Dunod, Paris.Ionescu T., Paumier J.-C., 1993. Friction dynamique avec coefficient dépendant de la vitesse de glissement. C. R. Acad. Sci. Paris Sér. I Math. 316

(1), 121–125.Jean M., Pratt E., 1985. A system of rigid bodies with dry friction. Internat. J. Engrg. Sci. 23 (5), 497–513.Laghdir M., Monteiro Marques M.D.P., 1995. Dynamics of a particle with damping, friction, and percussional effects. J. Math. Anal. Appl. 196 (3),

902–920.Madan R.N. (Ed.) 1993. Chua’s Circuit: A Paradigm for Chaos. World Scientific Publishing Co. Inc., River Edge, NJ. With 1 IBM-PC floppy disk

(3.5 inch; HD).Matrosov V.M., Finogenko I.A., 1995. Right-hand solutions of the differential equations of dynamics for mechanical systems with sliding friction.

J. Appl. Math. Mech. 59 (6), 837–844.Matrosov V.M., Finogenko I.A., 1996a. On the existence of right-side solutions to differential equations of dynamics of mechanical systems with

dry friction. Differential Equations 32 (2), 186–194.Matrosov V.M., Finogenko I.A., 1996b. To the theory of differential equations arising in the dynamics of systems with friction. I. Differential

Equations 32 (5), 610–618.Monteiro Marques M.D.P., 1994. An existence, uniqueness and regularity study of the dynamics of systems with one-dimensional friction. European

J. Mech. A/Solids 13 (2), 277–306.Palmov V., 1998. Vibrations of Elasto-Plastic Bodies. Springer-Verlag, Berlin. Translated from the 1976 Russian original by A. Belyaev and revised

by the author.Popp K., Stelter P., 1990. Stick-slip vibrations and chaos. Phil. Trans. R. Soc. Lond. A 332, 89–105.Pratt T., Williams R., 1981. Non-linear analysis of stick-slip motion. J. Sound Vibration 74 (4), 531–542.Shaw S.W., 1986. On the dynamic response of a system with dry friction. J. Sound Vibration 108 (2), 305–325.Stelter P., 1992. Nonlinear vibrations of structures induced by dry friction. Nonlinear Dynamics 3, 329–345.Stewart D.E., 1996. A numerical method for friction problems with multiple contacts. J. Austral. Math. Soc. Ser. B 37, 288–308.Stewart D.E., Trinkle J.C., 1996. An implicit time-stepping scheme for rigid body dynamics with inelastic collisions and Coulomb friction. Internat.

J. Numer. Methods Engrg. 39 (15), 2673–2691.Stewart D.E., Trinkle J.C., 1997. Dynamics, friction, and complementarity problems. In: Complementarity and Variational Problems (Baltimore,

MD, 1995), SIAM, Philadelphia, PA, pp. 425–439.


Tomlinson G., Chen Q., 1996. Parametric identification of systems with dry friction and nonlinear stiffness using a time series model. Journal ofVibration and Acoustics 118 (2), 252–263.

Trinkle J.C., Pang J.S., Sudarsky S., Lo G., 1997. On dynamic multi-rigid-body contact problems with coulomb friction. Z. Angew. Math. Mech. 77(4), 267–279.

Documents

Study of some rheological models with a finite number of degrees of freedom