# An overview of the formulation, existence and uniqueness issues for the initial value problem raised by the dynamics of discrete systems with unilateral contact and dry friction

Patrick Ballard, Alexandre Charles

► To cite this version:

Patrick Ballard, Alexandre Charles. An overview of the formulation, existence and uniqueness issues for the initial value problem raised by the dynamics of discrete systems with unilateral contact and dry friction. Comptes Rendus Mécanique, 2018, 346 (3), pp.222-236. 10.1016/j.crme.2017.12.010. hal-01737939

HAL Id: hal-01737939

https://hal.sorbonne-universite.fr/hal-01737939

Submitted on 20 Mar 2018

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
The legacy of Jean-Jacques Moreau in mechanics

# An overview of the formulation, existence and uniqueness issues for the initial value problem raised by the dynamics of discrete systems with unilateral contact and dry friction

<figure><img src="image_2.png" /></figure>

Check for
updates

Une revue des questions de formulation, existence, unicité, soulevées par la dynamique des systèmes discrets en présence de contact et frottement sec

Patrick Ballard a,*, Alexandre Charles b

ARTICLE INFO

Article history :

Received 20 May 2017

Accepted 11 October 2017

Available online 15 December 2017

Keywords :

Dynamics

Rigid solid

Discrete system

Impact

Friction

Contact

Existence

Uniqueness

Mots-clés :

Dynamique

Solide rigide

Système discret

Impact

Frottement

Contact

Existence

Unicité

A B S T R A C T

In the end of the seventies, Schatzman and Moreau undertook to revisit the venerable dynamics of rigid bodies with contact and dry friction in the light of more recent mathematics. One claimed objective was to reach, for the first time, a mathematically consistent formulation of an initial value problem associated with the dynamics. The purpose of this article is to make a review of the today state-of-art concerning not only the formulation, but also the issues of existence and uniqueness of solution.

© 2017 Académie des sciences. Published by Elsevier Masson SAS. This is an open access

article under the CC BY-NC-ND license

(http://creativecommons.org/licenses/by-nc-nd/4.0/).

RÉSUMÉ

À la fin des années 70, Schatzman et Moreau entreprirent de reformuler l'antique dynamique des solides rigides en présence de contact et frottement sec à la lumière de mathématiques plus récentes. Un des objectifs revendiqués était de parvenir, pour la première fois, à la formulation d'un problème d'évolution à partir d'une condition initiale, qui soit mathématiquement cohérent. Le but de cet article est de brosser un état de l'art actuel, concernant non seulement les questions de formulation, mais également d'existence et d'unicité de solution.

© 2017 Académie des sciences. Published by Elsevier Masson SAS. This is an open access

article under the CC BY-NC-ND license

(http://creativecommons.org/licenses/by-nc-nd/4.0/).

<footnote>a Sorbonne Universités, UPMC (Université Paris-6), CNRS, UMR 7190, Institut Jean-Le-Rond-d'Alembert, 75005 Paris, France</footnote>

<footnote>b Safran Tech, rue des Jeunes-Bois, Châteaufort, CS 80112, 78772 Magny-les-Hameaux cedex, France</footnote>

<footnote>* Corresponding author.</footnote>

<footnote>E-mail addresses: patrick.ballard@dalembert.upmc.fr (P. Ballard), alexandre.charles@safrangroup.com (A. Charles).</footnote>
# 1. Introduction

The dynamics of rigid solids with contact and dry friction conditions is a venerable subject, which was developed mainly in the second half of the nineteenth century and the beginning of the twentieth century to answer some questions raised by engineering. Then, the attention of engineers began to be driven towards elasticity and continuum mechanics, and less attention was paid to frictional contact multibody dynamics. In the seventies, a renewal of interest occurred, mainly driven by the development of numerical modelling in granular dynamics and the control issues associated with robotics. For both concerns, it turned out that the foundations of the venerable theory were not firm enough and that they should be reconsidered in the light of more recent mathematics.

A fundamental impulse was given simultaneously by Michelle Schatzman [1] and Jean-Jacques Moreau [2], who first considered an evolution problem for the configuration $\mathbf{q} : [0, T] \to \mathbb{R}^d$ in the framework of functions whose second derivative (in the distributional sense) is a Radon measure. At that time, the antique point of view of different systems of equations applying to the different phases of motion (without or with active contact) was still prevailing. Their new point of view permitted to formulate, for the first time, a mathematical evolution problem associated with multibody contact dynamics. It also paved the road for the design of efficient strategies for numerical computing and enabled the first investigations on the stability and control issues that are crucial in the analysis of the problems facing frictional contact events in robotics.

The seminal work of Michelle Schatzman and Jean-Jacques Moreau also initiated a series of contributions on the general formulation and the mathematical analysis of the initial value problem associated with multibody contact dynamics. A brief sketch of the history follows.

- The first studies about the formulation were restricted to the model problem of the dynamics of a particle evolving in an admissible region of $\mathbb{R}^d$ bounded by an obstacle. In [1], Michelle Schatzman formulated a consistent evolution problem in the frictionless case, under the additional restriction that the admissible region is convex. She was able to successfully implement a penalty method to prove the existence of a solution for the initial value problem. Her original work was restricted to impacts preserving the kinetic energy (the so-called elastic impact law) and an external force depending only on time. This result was generalized later by Paoli in her PhD thesis [3] to the case of an arbitrary impact law and an external force possibly also depending on current position and velocity. In parallel, an alternative strategy for proving the existence of a solution was designed by Monteiro Marques [4]. He introduced a time-stepping approximation and proved the convergence (of a subsequence) towards a solution. It was restricted to the completely inelastic impact law (zero restitution coefficient), but he was able to relax the convexity assumption of Schatzman on the admissible region. More importantly, he was able to generalize to the case where the contact with the obstacle obeys the Coulomb law of dry friction with a given friction coefficient $\mu$ (the frictionless case is recovered by taking $\mu = 0$). One benefit of this new strategy is that it directly suggested an algorithm for numerical computations. The time-stepping approach was further developed by Paoli [5] and [6], who, in particular, extended it to the case of an arbitrary restitution coefficient. Her work, however, is up to now concerned only with the frictionless case.

- It was recognized very early by Michelle Schatzman that issues should be expected with the uniqueness of the solution for the initial value evolution problem. In particular, she exhibited in [1] a striking example of multiple solutions for the unilateral dynamics of the one-degree-of-freedom particle submitted to an external force that is a $C^\infty$ function of time. This issue was further considered by Percivale [7,8], who noticed that the uniqueness of the solution could be recovered in the one-degree-of-freedom problem considered by Schatzman, provided that the given external force was assumed to be not only a $C^\infty$, but also an analytic function of time. His work was suggesting that uniqueness could be expected in general, provided a regularity assumption of analyticity on the data. This was proved in full generality for the frictionless problem by Ballard in [9] and [10]. Local uniqueness in the analytic framework was also exploited in this work to design a third alternative strategy (in addition to penalty and time-stepping methods) to prove the existence of the solution. This new strategy turned out to yield more general (except for the additional assumption of analyticity) existence results than those which were available at that time from the penalty and time-stepping methods, encompassing the multi-constraint case with an arbitrary impact law. This strategy of proof was also adapted in [11] and [12] to the dynamics of a point particle with contact conditions and Coulomb friction. It yields a slightly more general (except for the analyticity assumption) existence result than that of Monteiro Marques, and provides, in addition, the uniqueness of the solution from a given initial value.

- In the eighties and the nineties, most of the articles that appeared on the subject of the mathematical formulation of the initial value problem and the issues of existence and uniqueness of solutions were restricted to the model problem of a point particle evolving in $\mathbb{R}^d$, or rather in an admissible region of $\mathbb{R}^d$. There is one noteworthy exception: the seminal article [2] by Jean-Jacques Moreau. In this article, Jean-Jacques Moreau addressed the formulation of the dynamics of a collection of rigid bodies submitted to frictionless unilateral constraints such as the ones arising from the non-interpenetration conditions. The framework is from the beginning that of Lagrange in which the motion is represented as a curve $\mathbf{q}(t)$ in the configuration space, identified with a subset of $\mathbb{R}^d$. In that framework and in accordance with the ideas of Lagrange, the reaction force that appears in the formulation is a *generalized reaction force*. In particular, the detailed distribution of reaction forces in the real world (meaning forces in $\mathbb{R}^3$ from one body onto another) is generally undefined. The existence and uniqueness result of Ballard in [9] and [10] applies to this general framework for frictionless unilateral multibody dynamics, under the assumption of analyticity of the data. It yields existence and
uniqueness for the motion $q(t)$ and for the generalized reaction force, but the detailed distribution of reaction forces in the real world remains undetermined in general, and there may be several such distributions that are compatible with the generalized reaction force that is associated with the solution. The use of generalized forces, originating in Lagrange's idea, has now been made systematic all over continuum mechanics under the name 'Principle of Virtual Power'. It conveys the idea that the representation of forces within a mechanical theory must be made consistent with the kinematics: forces must be taken in the dual space of the velocity space, the duality product between forces and velocities being nothing but the power. Surprisingly, coming to frictional unilateral multibody dynamics, the use (originating from the nineteenth century) was to invoke the Coulomb law of dry friction applying to reaction forces in the real world, contradicting Lagrange's point of view and the Principle of Virtual Power that stipulates that the dynamics should be formulated in terms of appropriate generalized forces. However, in some simple cases involving a small number of contacting rigid bodies, a mathematical evolution problem can still be formulated based on this historical point of view [13], although no modern formulation of the general evolution problem has ever been obtained in that setting. Unsurprisingly, adopting the historical point of view of real world forces (when it is possible) raises some inconsistencies, known as Painlevé and Kane paradoxes [19] and [22,23]. Because of this lack of consistency within the usual view about the formulation of frictional unilateral multibody dynamics, no precise formulation of the evolution problem was obtained, and the only mathematical investigations about the existence and uniqueness of solutions for frictional unilateral dynamics have been so far restricted to the case of a finite collection of point particles. That was only recently that such a consistent general formulation in the line of Lagrange's ideas about the use of generalized forces was derived by Charles [14] in his PhD thesis.

## 2. The one-degree-of-freedom problem

Consider a point particle, of unit mass, that is constrained to move along a line, the location of which is represented by the abscissa $q \in \mathbb{R}$. We assume that an obstacle is located at the origin so that the particle is constrained to remain in the half-line defined by $q \ge 0$. To enforce this constraint during an arbitrary motion $q(t)$ of the particle, an unknown reaction $r$ force must be added in the equation of motion:

$$ \ddot{q} = f + r \quad (1) $$

where the external force $f(t)$ is supposed to be a given (integrable) function of time only (for the sake of simplicity). The usual physical assumption is that the existence of an obstacle has no influence on the motion of the particle when contact is not active, and the reaction force must therefore be supported in those instants where the contact is active:

$$ \text{Supp } r \subset \{t \mid q(t) = 0\} $$

Elementary examples then show that velocity jumps cannot be avoided in this framework and the acceleration that appears in the equation of motion (1) should be understood in the sense of (Schwartz's) distributions. Hence, the reaction force $r$ should not be expected to be a function, but rather a distribution. As it is usually assumed that the obstacle is only able to repel the particle, a nonnegativity assumption also has to be required on the unknown reaction force: $r \ge 0$. In the extended framework where $r$ is a distribution, this can only mean that the distribution $r$ returns nonnegative real values when tested by means of a nonnegative $C^\infty$ trial function with compact support. But, it is a classical (and easy) result that such a nonnegative distribution must actually be a nonnegative measure. Hence, given a bounded time interval $[0, T]$, the largest possible functional space in which the motion $q(t)$ can be sought is the space $MMA([0, T])$ (the acronym standing for 'Motions with Measure Acceleration') of those distributions on $[0, T]$ whose second derivative is a Radon measure $(\dot{q} \in \mathcal{M}([0, T]))$. Distributions in that space $MMA([0, T])$ are actually continuous functions, admitting left and right derivatives $\dot{q}^{-}(t)$, $\dot{q}^{+}(t)$ (in the classical sense) at every instant $t \in ]0, T[$. The side derivatives $\dot{q}^{-}(t)$ and $\dot{q}^{+}(t)$ are actually equal, except possibly at some instants belonging to a countable subset of $[0, T]:$ the impact instants. The two functions $\dot{q}^{-}(t)$ and $\dot{q}^{+}(t)$ are functions with bounded variation.

Finally, given an initial condition $(q_0, v_0) \in \mathbb{R}^{+} \times \mathbb{R}$ compatible with the obstacle ($q_0 = 0 \Rightarrow v_0 \ge 0$), the evolution problem associated with the unilateral dynamics of the particle reads as follows.

**Problem P_1.** Find $q \in MMA([0, T])$ and $r \in \mathcal{M}([0, T])$ such that:

$$ - q(0) = q_0, \quad \dot{q}^{+}(0) = v_0, $$

$$ - \ddot{q} = f + r, $$

$$ - \forall t \in [0, T], \quad q(t) \ge 0, $$

$$ - \text{Supp } r \subset \{t \in [0, T] \mid q(t) = 0\}, $$

$$ - r \ge 0, $$

$$ - q(t) = 0 \Rightarrow \dot{q}^{+}(t) = -e\dot{q}^{-}(t). $$

Here, the last line is an additional requirement with respect to the introductory discussion. If it were not stated, a particle free of external force impacting the obstacle could either subsequently remain stuck on the obstacle, or bouncing according
to a sign-reversed velocity (infinitely many intermediate choices are possible), and all of these events would be compatible with both the equation of motion and the contact conditions. Such an indeterminacy was already noted by Newton, who introduced the concept of a restitution coefficient to recover determinism. The indeterminacy originates in the fact that real bodies are always deformable and their bouncing is actually governed by the deformation waves that travel within the body during an impact. As the model of a point particle is too coarse to describe these deformation waves, the outcome of the impact remains undetermined and should be artificially added to the equations by means of an impact law, which has the same status as constitutive laws in continuum mechanics. The restitution coefficient $e$ is supposed to be given in the interval $[0, 1]$ (the choice of $e=1$ is usually baptized the elastic impact law and that of $e=0$, the completely inelastic impact law).

Actually, it was only recently that it was discovered that the original proposition of Newton to add the impact law fails to yield the uniqueness of the solution, even in the case where $f$ is supposed to be a $C^\infty$ function of time. The first discovery of this striking fact seems to be that of Bressan [15], although he did not have such a precise statement of the evolution problem. The construction of such a $C^\infty$ function $f$ enabling multiple solutions $q$ can be found in [1] in the case of the elastic impact law $e=1$ and in [9] in the case of the completely inelastic impact law $e=0$. In both cases, the $C^\infty$ functions $f$ that are exhibited are highly oscillating (having in particular infinitely many zeroes in the bounded interval $[0, T]$).

The following results about problem $\mathcal{P}_1$ have been proved in various articles.

(i) Suppose that $f \in L^1(0, T)$ and $e \in [0, 1]$. For all $\varepsilon > 0$, the penalized initial value problem:

$$ \left| \begin{array}{l} q_{\varepsilon}(0) = q_0, \qquad \dot{q}_{\varepsilon}(0) = v_0, \\ \ddot{q}_{\varepsilon}(t) = f(t) - \frac{2|\log e|}{\sqrt{\varepsilon[\pi^2 + (\log e)^2]}} \dot{q}_{\varepsilon}(t) \operatorname{sgn}^{-1}(q_{\varepsilon}(t)) - \frac{q_{\varepsilon}(t)}{\varepsilon} \operatorname{sgn}^{-1}(q_{\varepsilon}(t)), \quad \text{for a.a. } t \in [0, T] \end{array} \right. $$

(where $\operatorname{sgn}^{-1}$ is the function taking the value 1 on $]-\infty, 0]$ and 0 on $[0, \infty[$), has a unique solution $q_{\varepsilon} \in W^{2,1}(0, T)$. As $\varepsilon$ goes to zero, one can extract a subsequence in $q_{\varepsilon}$ that converges strongly in $W^{1,1}(0, T)$ towards some $q \in \text{MMA}([0, T])$ that solves problem $\mathcal{P}_1$. This result was proved in [3], whereas the simpler case $e=1$ (where the damping term in the penalty differential equation reduces to 0) had been previously treated in [1]. Hence, although the penalty approach fits naturally with the particular case of the elastic impact law $e=1$, it can be extended to dissipative impact laws.

(ii) Suppose that $f \in L^1(0, T)$ and $e \in [0, 1]$. Picking $n \in \mathbb{N} \setminus \{0\}$ and setting $h = T/n$, we define a sequence of approximants $q_n \in C^0([0, T])$ by the following induction.

$$ - Q_n^0 = q_0, \quad V_n^0 = v_0, $$

$$ - \forall i \in \{1, 2, \dots, n\}, $$

$$ F_n^i = \frac{1}{h} \int_{(i-1)h}^{ih} f, $$

$$ V_n^i = \begin{cases} V_n^{i-1} + hF_n^i, & \text{if } Q_n^{i-1} > 0 \\ -eV_n^{i-1} + \langle (1+e)V_n^{i-1} + hF_n^i \rangle^+, & \text{if } Q_n^{i-1} \le 0, \end{cases} $$

$$ - q_n(t) = Q_n^{i-1} + (t - (i-1)h)V_n^{i-1}, \quad \forall t \in [(i-1)h, ih], $$

where $\langle x \rangle^+ = \max\{x, 0\}$ stands for the positive part function. Then, it was proved in [6] that a subsequence of $q_n(t)$ converges strongly in $W^{1,1}(0, T)$ towards some $q \in \text{MMA}([0, T])$ that solves problem $\mathcal{P}_1$. The simpler case $e=0$ had been previously treated in [4]. Hence, although the time-stepping approach fits naturally with the particular case of the completely inelastic impact law $e=0$, it can be extended to the general case of an arbitrary restitution coefficient $e \in [0, 1]$.

(iii) Uniqueness for problem $\mathcal{P}_1$ can be recovered in $\text{MMA}([0, T])$, if the function $f(t)$ is more regular than $C^\infty$. In the analytic case, uniqueness was first proved for the particular case $e=1$ in [7], and for an arbitrary restitution coefficient $e \in [0, 1]$ in [9]. Needless to say, these results apply in the case of piecewise analyticity.

(iv) If the function $f(t)$ is analytic, or piecewise analytic, then the solution mapping:

$$ \left\{ \begin{array}{l} \mathbb{R}^2 \rightarrow C^0([0, T]) \\ (q_0, v_0) \mapsto q(t) \end{array} \right. $$

is continuous. In other words, the unique solution to problem $\mathcal{P}_1$ depends continuously on the initial data. A proof is to be found in [9].
### 3. Frictional unilateral dynamics of a point particle

We now consider the motion of a point particle, of unit mass, in a region of $\mathbb{R}^d$ defined by the unilateral constraint $\varphi(\mathbf{q}) \ge 0$, where $\varphi$ is a smooth (of class $C^1$ at least) function such that:

$$ \varphi(\mathbf{q}) = 0 \quad \Rightarrow \quad d\varphi_{\mathbf{q}} \neq \mathbf{0} $$

where $d\varphi_{\mathbf{q}}$ denotes the differential of the real-valued function $\varphi$ at point $\mathbf{q}$. The equation $\varphi(\mathbf{q}) = 0$ defines the geometry of the obstacle, and $d\varphi_{\mathbf{q}}$ defines the direction of the outward normal to the obstacle. When the particle is in contact with the obstacle, the reaction force exerted by the obstacle will be denoted by $\mathbf{r} \in \mathbb{R}^d$. It can be classically split into normal and tangential parts:

$$ \mathbf{r} = \mathbf{r}_t + r_n \frac{d\varphi_{\mathbf{q}}}{|d\varphi_{\mathbf{q}}|}, \quad \text{with } r_n = \left\langle \mathbf{r}, \frac{d\varphi_{\mathbf{q}}}{|d\varphi_{\mathbf{q}}|} \right\rangle $$

where $\mathbf{q}$ denotes the location of the particle at contact and $\langle \cdot, \cdot \rangle, |\cdot|$ are the canonical scalar product and norm in $\mathbb{R}^d$. The case where dry friction between the particle and the obstacle can occur will be considered. The simplest law describing dry friction is the empirical law of Coulomb, which reads as:

$$ |\mathbf{r}_t| \le \mu r_n \quad \text{and} \quad \begin{cases} \text{if } |\mathbf{r}_t| < \mu r_n, & \text{then } \dot{\mathbf{q}}_t^+ = \mathbf{0}, \\ \text{if } |\mathbf{r}_t| = \mu r_n, & \text{then } \exists \lambda \in \mathbb{R}^+, \quad \mathbf{r}_t = -\lambda \dot{\mathbf{q}}_t^+ \end{cases} \quad (2) $$

where $\dot{\mathbf{q}}_t^+$ denotes the tangential (with respect to the obstacle) component of the right-velocity, and $\mu \ge 0$ is a given friction coefficient. The Coulomb law above can be compactly and equivalently rewritten under the weak form:

$$ \forall \mathbf{v} \in \mathbb{R}^d, \quad \left\langle \mathbf{r}_t, \mathbf{v}_t - \dot{\mathbf{q}}_t \right\rangle + \mu r_n \left( |\mathbf{v}_t| - |\dot{\mathbf{q}}_t| \right) \ge 0 $$

where $\langle \cdot, \cdot \rangle$ stands for the canonical scalar (duality) product in $\mathbb{R}^d$, as seems to have been first pointed out independently in the 1960s by Jean-Jacques Moreau and Georges Duvaut. This weak form turns out to be the appropriate form to extend the classical Coulomb law to the case where the reaction force $\mathbf{r}$ is a Radon measure (with respect to time):

$$ \forall \mathbf{v} \in C^0([0, T]; \mathbb{R}^d), \quad \int_{[0, T]} \left\langle \mathbf{r}_t, \mathbf{v}_t(t) - \dot{\mathbf{q}}_t(t) \right\rangle + \mu r_n \left( |\mathbf{v}_t(t)| - |\dot{\mathbf{q}}_t(t)| \right) \ge 0 $$

Note that if $\mathbf{q} \in \text{MMA}([0, T]; \mathbb{R}^d)$, then $\dot{\mathbf{q}}^+ : [0, T] \to \mathbb{R}^d$ is a function with bounded variation that is therefore universally integrable (integrable with respect to any measure). Hence, the integral in the above weak form of the Coulomb law is well-defined for $\mathbf{q} \in \text{MMA}([0, T]; \mathbb{R}^d)$ and $\mathbf{r} \in \mathcal{M}([0, T]; \mathbb{R}^d)$. This weak form can be assumed globally in $[0, T]$, so that the Coulomb law will be enforced both during the impacts and the possible smooth phases of motion along the obstacle.

Given an initial condition $(\mathbf{q}_0, \mathbf{v}_0) \in \mathbb{R}^d \times \mathbb{R}^d$ compatible with the obstacle:

$$ \varphi(\mathbf{q}_0) \ge 0 \quad \text{and} \quad \varphi(\mathbf{q}_0) = 0 \Rightarrow \langle d\varphi_{\mathbf{q}_0}, \mathbf{v}_0 \rangle \ge 0 $$

we can now formulate the initial value problem that governs the dynamics of the particle in frictional contact with the obstacle. The external force $\mathbf{f}(t; \mathbf{q}, \dot{\mathbf{q}}^{-})$ is now allowed to depend on the current location and velocity (here, we can equivalently use either $\dot{\mathbf{q}}^{-}$ or $\dot{\mathbf{q}}^{+}$, with no influence on the solution).

**Problem P2.** Find $\mathbf{q} \in \text{MMA}([0, T]; \mathbb{R}^d)$ and $\mathbf{r} \in \mathcal{M}([0, T]; \mathbb{R}^d)$ such that:

$$ - \mathbf{q}(0) = \mathbf{q}_0, \quad \dot{\mathbf{q}}^+(0) = \mathbf{v}_0, $$

$$ - \ddot{\mathbf{q}} = \mathbf{f}(t; \mathbf{q}, \dot{\mathbf{q}}^{-}) + \mathbf{r}, $$

$$ - \forall t \in [0, T], \quad \varphi(\mathbf{q}(t)) \ge 0, $$

$$ - \text{Supp} \mathbf{r} \subset \{t \in [0, T] \mid \varphi(\mathbf{q}(t)) = 0\}, $$

$$ - r_n \ge 0, $$

$$ - \forall \mathbf{v} \in C^0([0, T]; \mathbb{R}^d), \quad \int_{[0, T]} \left\langle \mathbf{r}_t, \mathbf{v}_t(t) - \dot{\mathbf{q}}_t(t) \right\rangle + \mu r_n \left( |\mathbf{v}_t(t)| - |\dot{\mathbf{q}}_t(t)| \right) \ge 0, $$

$$ - \varphi(\mathbf{q}(t)) = 0 \quad \Rightarrow \quad \dot{\mathbf{q}}_n^+(t) = -e \dot{\mathbf{q}}_n^-(t). $$

The following results about problem P2 have been proved in various articles.
(i) In the frictionless case $\mu = 0$, with given restitution coefficient $e \in [0, 1]$, in the case where $f$ is a continuous function that is globally Lipschitz with respect to $(\mathbf{q}, \dot{\mathbf{q}})$, then it is possible [3] to extend the penalty method described in section 2 (i), at least when the admissible subset of $\mathbb{R}^d$ defined by $\varphi(\mathbf{q}) \ge 0$ is convex. In the same way, there exists a subsequence of the penalty approximates $\mathbf{q}_\varepsilon$ that converges strongly in $W^{1,1}(0, T)$, as $\varepsilon \to 0+$, towards some $\mathbf{q} \in \text{MMA}([0, T])$ that solves problem $\mathcal{P}_2$ (a proof is to be found in [3] in the case of a convex admissible subset of $\mathbb{R}^d$).

(ii) In the frictionless case $\mu = 0$, with given restitution coefficient $e \in [0, 1]$, in the case where $f$ is a continuous function that is globally Lipschitz with respect to $(\mathbf{q}, \dot{\mathbf{q}})$, then it is possible [6] to extend the time-stepping method described in section 2 (ii). In the same way, there exists a subsequence of the time-stepping approximates $\mathbf{q}_n$ that converges strongly in $W^{1,1}(0, T)$, as $n \to +\infty$, towards some $\mathbf{q} \in \text{MMA}([0, T])$ that solves problem $\mathcal{P}_2$ (a proof is to be found in [6]).

(iii) In the frictional case $\mu \ge 0$, and the completely inelastic impact law $e = 0$, in the case where $f$ is a continuous bounded function of $t$ and $\mathbf{q}$ only, Monteiro Marques [4] defines a sequence of approximants $\mathbf{q}_n \in C^0([0, T])$ ($n \in \mathbb{N} \setminus \{0\}$) by the following induction.

$$- \mathbf{Q}_n^0 = \mathbf{q}_0, \quad \mathbf{V}_n^0 = \mathbf{v}_0,$$

$$- \forall i \in \{1, 2, \dots, n\},$$

$$\mathbf{Q}_n^i = \mathbf{Q}_n^{i-1} + h\mathbf{V}_n^{i-1},$$ $$\mathbf{V}_n^i = \mathbf{V}_n^{i-1} + h\mathbf{f}(ih, \mathbf{Q}_n^i),$$

$$\mathbf{V}_n^i = \begin{cases} \mathbf{V}_n^i, & \text{if } \mathbf{V}_n^i \in \mathcal{V}(\mathbf{Q}_n^i) \\ \operatorname{Proj}(\mathbf{0}; [\mathbf{V}_n^i + \mathcal{C}(\mathbf{Q}_n^i)] \cap \mathcal{T}(\mathbf{Q}_n^i)), & \text{if } \mathbf{V}_n^i \notin \mathcal{V}(\mathbf{Q}_n^i) \end{cases}$$

$$- \mathbf{q}_n(t) = \mathbf{Q}_n^{i-1} + (t - (i-1)h)\mathbf{V}_n^{i-1}, \quad \forall t \in [(i-1)h, ih],$$

where $h = T/n$, $\mathcal{V}(\mathbf{q})$ is the set of admissible right-velocities at location $\mathbf{q}$:

$$\mathcal{V}(\mathbf{q}) = \begin{cases} \{\mathbf{v} \in \mathbb{R}^d \mid \langle d\varphi_{\mathbf{q}}, \mathbf{v} \rangle \ge 0\}, & \text{if } \varphi(\mathbf{q}) \le 0, \\ \mathbb{R}^d, & \text{if } \varphi(\mathbf{q}) > 0 \end{cases}$$

($d\varphi_{\mathbf{q}}$ is assumed to be nowhere $\mathbf{0}$), $\mathcal{T}(\mathbf{q})$ is the tangent hyperplane:

$$\mathcal{T}(\mathbf{q}) = \{\mathbf{v} \in \mathbb{R}^d \mid \langle d\varphi_{\mathbf{q}}, \mathbf{v} \rangle = 0\}$$

and $\mathcal{C}(\mathbf{q})$ is the friction cone:

$$\mathcal{C}(\mathbf{q}) = \{\mathbf{r} \in \mathbb{R}^d \mid \langle \mathbf{r}, d\varphi_{\mathbf{q}} \rangle \ge |\mathbf{r}| |d\varphi_{\mathbf{q}}| / \sqrt{1 + \mu^2}\}$$

Then, there exists a subsequence in the sequence $(\mathbf{q}_n)$ that converges strongly in $W^{1,1}(0, T)$, as $n \to +\infty$, towards some $\mathbf{q} \in \text{MMA}([0, T])$ that solves problem $\mathcal{P}_2$ (a proof is to be found in [4]).

(iv) In the frictional case $\mu \ge 0$, with given restitution coefficient $e \in [0, 1]$, in the case where $\varphi$ is an analytic function and $f$ is an analytic function that is globally Lipschitz with respect to $(\mathbf{q}, \dot{\mathbf{q}})$, then there exists a unique solution to problem $\mathcal{P}_2$ (a proof is to be found in [12]).

# 4. Frictionless unilateral multibody dynamics

In the preceding sections, only the dynamical evolution of a point particle has been considered. However, the applicability of the theory to situations of practical interest requires to extend it to the case of one or several rigid bodies. The framework must therefore be extended to that of Lagrange about discrete mechanical systems (that encompasses in particular any finite collection of rigid bodies, some of them being possibly connected by so-called perfect joints).

## 4.1. Virtual power and Lagrange equations

A discrete mechanical system is a mechanical system whose arbitrary configuration in the space can be described by a finite number $d$ of independent real numbers: the generalized coordinates denoted by $q^1, q^2, \dots, q^d$. Here, "independent" means that it is always possible to conceive a motion for which all generalized coordinates but an arbitrary one remain fixed. The number $d$ is called the number of degrees of freedom of the discrete mechanical system. The notation $\mathbf{q} = (q^1, q^2, \dots, q^d)$ will be used and $\mathbf{q}$ will be called the abstract configuration (or sometimes simply, the configuration) of the system. A motion of the system is simply a mapping $\mathbf{q}(t)$ defined on some time interval and taking abstract configuration values. Its derivative with respect to time is denoted by $\dot{\mathbf{q}}(t)$ and the vector $\dot{\mathbf{q}}(t)$ is called the generalized velocity at time $t$. The generalized velocity is a convenient mathematical representation of the whole velocity field over the body (or bodies). The kinetic energy $K(\mathbf{q}, \dot{\mathbf{q}})$ is quadratic with respect to the generalized velocity and defines the so-called kinetic matrix $\mathbf{M}(\mathbf{q})$:
$$K(\mathbf{q}, \dot{\mathbf{q}}) = \frac{1}{2} \dot{\mathbf{q}} \cdot \mathbf{M}(\mathbf{q}) \cdot \dot{\mathbf{q}}$$

The kinetic matrix is always symmetric and positive. It is actually positive definite, if the system does not involve any massless components, which will always be assumed in the sequel.

More precisely, the set $Q$ of all abstract configurations is endowed with the structure of a differentiable manifold. This point of view makes it possible to consider an abstract configuration $\mathbf{q}$ with no need for a prior definition of generalized coordinates (that is, a parameterization of the system) and is particularly useful when it is important to distinguish intrinsic quantities, that is, quantities not relying on a specific choice of generalized coordinates on the configuration manifold. In this context, the generalized velocities belong to tangent spaces to the configuration manifold. The kinetic matrix endows each tangent space with a scalar product, so that the configuration manifold is actually a Riemannian manifold. When emphasis is to be put on intrinsic quantities, the following alternative notations will be possibly encountered in the sequel:

$$(\dot{\mathbf{q}}_1, \dot{\mathbf{q}}_2)_q = \dot{\mathbf{q}}_1 \cdot \mathbf{M}(\mathbf{q}) \cdot \dot{\mathbf{q}}_2, \qquad |\dot{\mathbf{q}}|_q = \sqrt{\dot{\mathbf{q}} \cdot \mathbf{M}(\mathbf{q}) \cdot \dot{\mathbf{q}}}$$

In any case, the point of view of Riemannian manifold on the configuration space, although enlightening, is not strictly needed, and there will be no harm in identifying in the sequel the configuration manifold with (an open subset of) $\mathbb{R}^d$ and the tangent space $T_q Q$ at configuration $\mathbf{q}$, with $\mathbb{R}^d$ endowed with the scalar product $\mathbf{M}(\mathbf{q})$. This identification is nothing but considering a particular choice of generalized coordinates.

The generalized force $\mathbf{f}$ of Lagrange is defined by means of the virtual power that the internal and external forces develop in any virtual generalized velocity $\hat{\mathbf{v}}$. This virtual power reads as:

$$\langle \mathbf{f}, \hat{\mathbf{v}} \rangle = f_i \hat{v}^i$$

(with the usual convention of summation on repeated indices) and shows that the generalized force belongs to the dual space $T_q^*Q$ of the vector space $T_qQ$ of all generalized velocities at configuration $\mathbf{q}$. In particular, the generalized force for a mechanical system with $d$ degrees of freedom has $d$ components. The virtual power point of view is usually used to compute the generalized force from a given distribution of forces $\tilde{\mathbf{f}}(\mathbf{x})$ in the real world (here $\mathbf{x}$ denotes the space variable in the three-dimensional space). First, the real world velocity $\hat{\mathbf{v}}(\mathbf{x})$ is computed in terms of the generalized velocity $\hat{\mathbf{v}}$:

$$\tilde{\mathbf{v}}(\mathbf{x}) = \mathbf{I}(\mathbf{q}, \mathbf{x}) \cdot \hat{\mathbf{v}}$$

where $\mathbf{I}(\mathbf{q}, \mathbf{x}) : T_q Q \mapsto \mathbb{R}^3$ is a linear mapping depending in general on the current configuration $\mathbf{q}$ and $\mathbf{x}$. Then, the virtual power paradigm:

$$\forall \hat{\mathbf{v}}, \quad \langle \mathbf{f}, \hat{\mathbf{v}} \rangle = \int \tilde{\mathbf{f}}(\mathbf{x}) \cdot \tilde{\mathbf{v}}(\mathbf{x}) \, d\mathbf{x} = \left\langle \int t \mathbf{I}(\mathbf{q}, \mathbf{x}) \cdot \tilde{\mathbf{f}}(\mathbf{x}) \, d\mathbf{x}, \hat{\mathbf{v}} \right\rangle \quad \Rightarrow \quad \mathbf{f} = \int t \mathbf{I}(\mathbf{q}, \mathbf{x}) \cdot \tilde{\mathbf{f}}(\mathbf{x}) \quad (3)$$

provides the expression of the generalized force $\mathbf{f}$ in terms of the real world force distribution $\tilde{\mathbf{f}}(\mathbf{x})$ (here, in the case where $\tilde{\mathbf{f}}(\mathbf{x})$ consists in finitely many point forces, the above integral reduces to a finite sum). Let us point out once more than it is always possible to compute the generalized force from a real world force distribution but that the real world force distribution cannot be recovered in general from the generalized force.

The generalized acceleration $\gamma$ of Lagrange is defined by means of the virtual power it develops in any virtual generalized velocity $\hat{\mathbf{v}}$.

$$\langle \gamma, \hat{\mathbf{v}} \rangle = \left( \frac{D\dot{\mathbf{q}}}{Dt}, \hat{\mathbf{v}} \right)_q = \left( \frac{d}{dt} \frac{\partial K}{\partial \dot{\mathbf{q}}^i} - \frac{\partial K}{\partial q^i} \right) \hat{\mathbf{v}}^i$$

where $D/Dt$ stands for the covariant derivative along the motion. It satisfies:

$$\frac{dK}{dt} = \left( \frac{D\dot{\mathbf{q}}}{Dt}, \dot{\mathbf{q}} \right)_q$$

The fundamental principle of classical dynamics asserts the Lagrange equation of motion $\gamma = \mathbf{f}$, which is equivalent to the principle of virtual power:

$$\forall \hat{\mathbf{v}} \in T_q Q, \quad \left( \frac{D\dot{\mathbf{q}}}{Dt}, \hat{\mathbf{v}} \right)_q = \langle \mathbf{f}, \hat{\mathbf{v}} \rangle$$
## 4.2. Frictionless unilateral constraints in Lagrange's setting

In the case of several rigid bodies, it must be expected, in general, that there will be several unilateral constraints: $\varphi_{\alpha}(\mathbf{q}) \ge 0$ ($\alpha = 1, 2, \dots, n$), as the non-interpenetration conditions give rise to at least one such unilateral constraint for each pair of rigid bodies. The set of all the indices for which the corresponding constraint is active in the configuration $\mathbf{q}$ will be denoted by $J(\mathbf{q})$:

$$J(\mathbf{q}) = \{\alpha | \varphi_{\alpha}(\mathbf{q}) = 0\}$$

The unilateral constraints $\varphi_{\alpha}$ are supposed to be such that, for any admissible configuration $\mathbf{q}$, the $d\varphi_{\alpha, \mathbf{q}}$ ($\alpha \in J(\mathbf{q})$) are linearly independent. The realization of the unilateral constraints requires to complement the equation of motion with a (generalized) reaction force $\mathbf{r}$. In the particular case of frictionless unilateral constraints, it takes the form:

$$\mathbf{r} = \sum_{\alpha=1}^{n} \lambda_{\alpha} d\varphi_{\alpha, \mathbf{q}}$$

with $\lambda_{\alpha} \ge 0$ and $\text{Supp} \lambda_{\alpha} \subset \{t \in [0, T] | \varphi_{\alpha}(\mathbf{q}(t)) = 0\}$.

Although the sole restitution coefficient $e$ is enough to convey the more general frictionless impact law in the case of a single unilateral constraint, a great deal of additional complexity of the impact law is permitted in the case of several active unilateral constraints at an impact. The impact law will be postulated under its more general possible form:

$$\dot{\mathbf{q}}^{+} = \mathcal{F}(\mathbf{q}, \dot{\mathbf{q}}^{-})$$

where $\mathbf{q}$ is the configuration of the system at the impact time, $\dot{\mathbf{q}}^{-}$ the left (impacting) velocity and $\dot{\mathbf{q}}^{+}$ the right (outgoing) velocity. The function $\mathcal{F}$ is a datum of the problem. In a practical situation, it has to be identified by using either experiments or a refined theory such as the mechanics of deformable bodies. Naturally, the function $\mathcal{F}$ cannot be arbitrary, it has to fulfil some compatibility conditions in order to be compatible with the equation of motion. More precisely, it must satisfy the three following conditions.

(i) $\forall \mathbf{q}, \dot{\mathbf{q}}^{-}, \forall \alpha \in J(\mathbf{q}), \quad \langle d\varphi_{\alpha, \mathbf{q}}, \mathcal{F}(\mathbf{q}, \dot{\mathbf{q}}^{-}) \rangle \ge 0$, \\
(the post-impact velocity must not violate the unilateral constraints),

(ii) $\forall \mathbf{q}, \dot{\mathbf{q}}^{-}, \quad \mathbf{M}(\mathbf{q}) \cdot (\mathcal{F}(\mathbf{q}, \dot{\mathbf{q}}^{-}) - \dot{\mathbf{q}}^{-}) \in \sum_{\alpha \in J(\mathbf{q})} \mathbb{R}^{+} d\varphi_{\alpha, \mathbf{q}}$, \\
(no friction: the generalized reaction force impulse is directed along the normal),

(iii) $\forall \mathbf{q}, \dot{\mathbf{q}}^{-}, \quad |\mathcal{F}(\mathbf{q}, \dot{\mathbf{q}}^{-})|_{\mathbf{q}} \le |\dot{\mathbf{q}}^{-}|_{\mathbf{q}}$, \\
(the kinetic energy cannot be increased by an impact).

There exist many functions $\mathcal{F}$ satisfying these requirements. One canonical example is the Moreau impact law [2]. It is based on the decomposition of an arbitrary vector $\mathbf{v}$ of the Euclidean vector space $T_{\mathbf{q}}Q \simeq \mathbb{R}^d$, endowed with the scalar product $\mathbf{M}(\mathbf{q})$:

$$\mathbf{v} = \text{Proj}_{\mathbf{M}(\mathbf{q})}(\mathbf{v}; \mathcal{V}(\mathbf{q})) + \text{Proj}_{\mathbf{M}(\mathbf{q})}(\mathbf{v}; \mathcal{N}(\mathbf{q}))$$

on the two mutually polar cones:

$$\mathcal{V}(\mathbf{q}) = \{ \mathbf{v} \in T_{\mathbf{q}}Q \mid \forall \alpha \in J(\mathbf{q}), \langle d\varphi_{\alpha, \mathbf{q}}, \mathbf{v} \rangle \ge 0 \}, \\ \mathcal{N}(\mathbf{q}) = \sum_{\alpha \in J(\mathbf{q})} \mathbb{R}^{-} \nabla \varphi_{\alpha, \mathbf{q}} = \sum_{\alpha \in J(\mathbf{q})} \mathbb{R}^{-} \mathbf{M}^{-1}(\mathbf{q}) \cdot d\varphi_{\alpha, \mathbf{q}}$$

where $\text{Proj}_{\mathbf{M}(\mathbf{q})}$ is the orthogonal projection operator with respect to the scalar product $\mathbf{M}(\mathbf{q})$ of $T_{\mathbf{q}}Q \simeq \mathbb{R}^d$. Then, the Moreau impact law [2] reads as:

$$\mathcal{F}(\mathbf{q}, \dot{\mathbf{q}}^{-}) = \text{Proj}_{\mathbf{M}(\mathbf{q})}(\dot{\mathbf{q}}^{-}; \mathcal{V}(\mathbf{q})) - e \text{ Proj}_{\mathbf{M}(\mathbf{q})}(\dot{\mathbf{q}}^{-}; \mathcal{N}(\mathbf{q})) \quad (4)$$

in which the restitution coefficient $e$ can be chosen arbitrarily in $[0, 1]$. It is an easy matter to prove that it fulfils requirements (i), (ii), and (iii) above. Reciprocally, at a location $\mathbf{q}$ where only one unilateral constraint is active ($\text{card } J(\mathbf{q}) = 1$), requirements (i), (ii) and (iii) above, entail that the most general impact law must take the form (4) of that of Moreau, for some $e \in [0, 1]$.
<figure><img src="image_3.png" /><figcaption>Fig. 1. Non-continuous dependence of the solution on the initial datum.</figcaption></figure>

4.3. The evolution problem associated with frictionless unilateral multibody dynamics

We are given an initial condition $(\mathbf{q}_0, \mathbf{v}_0) \in TQ \simeq \mathbb{R}^d \times \mathbb{R}^d$ assumed to be compatible with the unilateral constraints:

$$ \forall \alpha, \quad \varphi_{\alpha}(\mathbf{q}_0) \ge 0, \quad \text{and} \quad \forall \alpha \in J(\mathbf{q}_0), \quad \left\langle d\varphi_{\alpha, \mathbf{q}_0}, \mathbf{v}_0 \right\rangle \ge 0 $$

and an impact law function $\mathcal{F}$ fulfilling the three requirements (i), (ii), (iii) above. The initial value problem that governs frictionless unilateral multibody now reads as follows.

**Problem P3.** Find $\mathbf{q} \in \text{MMA}([0, T]; Q)$ and $\lambda_{\alpha} \in \mathcal{M}([0, T]; \mathbb{R})$ ($\alpha \in \{1, 2, \dots, n\}$) such that:

$$ - \mathbf{q}(0) = \mathbf{q}_0, \quad \dot{\mathbf{q}}^+(0) = \mathbf{v}_0, \quad \text{(initial condition)} $$

$$ - \mathbf{M}(\mathbf{q}) \cdot \frac{D\dot{\mathbf{q}}}{Dt} = \mathbf{f}(t; \mathbf{q}, \dot{\mathbf{q}}^{-}) + \sum_{\alpha=1}^{n} \lambda_{\alpha} d\varphi_{\alpha, \mathbf{q}}, \quad \text{(equation of motion)} $$

$$ - \forall \alpha, \quad \varphi_{\alpha}(\mathbf{q}) \ge 0, \quad \lambda_{\alpha} \ge 0, \quad \text{Supp} \, \lambda_{\alpha} \subset \{t \in [0, T] \mid \varphi_{\alpha}(\mathbf{q}(t)) = 0\}, \quad \text{(contact conditions)} $$

$$ - \forall t, \quad \dot{\mathbf{q}}^+ = \mathcal{F}(\mathbf{q}, \dot{\mathbf{q}}^{-}), \quad \text{(impact constitutive law).} $$

The following results about problem P3 have been proved in various articles.

(i) If the functions $\varphi_{\alpha}$, $\mathbf{M}$ and $\mathbf{f}$ are analytic functions of their arguments, then problem P3 has a unique maximal solution [9,10]. Here, 'maximal solution' means that the solution is not necessarily defined all over $[0, T]$, but rather on a subinterval that cannot be extended, since a possible blow up at finite time can occur. This blow up at finite time is classically known to be possible in the case where there is no unilateral constraint, but is dismissed whenever $\mathbf{f}$ satisfies a global Lipschitz condition with respect to $(\mathbf{q}, \dot{\mathbf{q}})$. The same is true for problem P3: if $\mathbf{f}$ is, in addition, Lipschitz with respect to $(\mathbf{q}, \dot{\mathbf{q}})$, then the maximal solution is defined all over $[0, T]$, and problem P3 has truly one and only one solution [9,10].

(ii) The above general existence and uniqueness result of solution for problem P3 should not put a mask on a major ill-posedness issue that is generally encountered by the solution: it does not depend continuously on the data, in general. This unpleasant feature can be observed in a simple example in $\mathbb{R}^2$, with $\mathbf{f} \equiv 0$, two affine functions $\varphi_{\alpha}$ (representing two intersecting straight obstacles), and an elastic impact law. In Fig. 1, two arbitrarily closed initial locations of a point particle are considered, associated with the same initial velocity, and they are seen to produce post-impact velocities that are orthogonal, so that the subsequent motions diverge. This pathology is very unpleasant because it precludes any rational attempt to compute an approximation of the unique solution to problem P3, in general. However, this issue is specific to the situation where the unique solution to problem P3 experiments at least one multiple impact (that is, passing through a location $\mathbf{q}$ at which several unilateral constraints are simultaneously active: card $J(\mathbf{q}) \ge 2$). Hence, the difficulty is connected with multiple impacts. However, all the multiple impacts do not give rise to this pathology. In particular, the counter-example of Fig. 1 would be impossible to build in the case where the two straight obstacles would intersect with a right angle. This was generalized as follows by Ballard.

**Proposition 1.** [9] Consider problem P3 in which the data $\varphi_{\alpha}$, $\mathbf{M}$, and $\mathbf{f}$ are analytic functions of their arguments and where the impact law is chosen to be that of Moreau (4) for an arbitrary $e \in [0, 1]$. Suppose, in addition, that $\mathbf{f}$ is Lipschitz with respect to $(\mathbf{q}, \dot{\mathbf{q}})$, which ensures that problem P3 admits a unique solution for any admissible initial condition. Consider such an initial $(\mathbf{q}_0, \mathbf{v}_0)$ and the corresponding solution $\mathbf{q} \in C^0([0, T]; \mathbb{R}^d)$. Suppose that this solution $\mathbf{q}(t)$ experiments only multiple impacts that are right with respect to the scalar product $\mathbf{M}^{-1}(\mathbf{q})$:

$$ \forall t \in [0, T], \quad \forall \alpha, \beta \in J(\mathbf{q}(t)), \quad d\varphi_{\alpha, \mathbf{q}(t)} \cdot \mathbf{M}^{-1}(\mathbf{q}(t)) \cdot d\varphi_{\beta, \mathbf{q}(t)} = 0 \quad (5) $$
then the mapping:

$$ \left\{ \begin{array}{rcl} TQ & \rightarrow & C^0([0, T]; Q) \\ (\mathbf{q}_0, \mathbf{v}_0) & \mapsto & \mathbf{q}(t) \end{array} \right. $$

is continuous at $(\mathbf{q}_0, \mathbf{v}_0)$.

Physically, condition (5) (the angle between two arbitrary constraints at a multiple impact is a right angle) is satisfied for a multiple impact consisting of simultaneous collisions at different locations by two pairs of solids in a dynamical evolution of a collection of rigid bodies, whereas it is generally not satisfied in the case of three bodies colliding simultaneously.

This issue of continuous dependence of the solution to problem $P_3$ with respect to the data was further investigated by Paoli [16] who proved that, in the particular case where $e = 0$ (completely inelastic impact law), the orthogonality condition in Proposition 1 can be weakened into:

$$ \forall t \in [0, T], \quad \forall \alpha, \beta \in J(\mathbf{q}(t)), \quad d\varphi_{\alpha, \mathbf{q}(t)} \cdot \mathbf{M}^{-1}(\mathbf{q}(t)) \cdot d\varphi_{\beta, \mathbf{q}(t)} \le 0 $$

that is, in the case of the completely inelastic impact law, the angles between the constraints at a multiple impact does not need to be right, but only acute, to ensure continuous dependence of the solution on the initial condition. She also proved that this condition of orthogonality of the constraints (or acute angles, in the case of the completely inelastic impact law) ensures not only continuous dependence of the solution on the initial condition, but also on the other data of the problem, such as $M$, $f$ and $\varphi_\alpha$, in the sense of appropriate topologies.

(iii) Assuming that the impact law in problem $P_3$ is that of Moreau (4) for an arbitrary $e \in [0, 1]$. Then, Paoli [6] defines a sequence of approximants $\mathbf{q}_n \in C^0([0, T])$ ($n \in \mathbb{N} \setminus \{0\}$) by the following induction.

$$ - \mathbf{q}_n^0 = \mathbf{q}_0, \quad \mathbf{v}_n^0 = \mathbf{v}_0, $$

$$ - \forall i \in \{1, 2, \dots, n\}, $$

$$ \mathbf{Q}_n^i = \mathbf{Q}_n^{i-1} + h\mathbf{V}_n^{i-1}, $$

$$ \mathbf{V}_n^i = \mathbf{V}_n^{i-1} + h \mathbf{M}^{-1}(\mathbf{Q}_n^i) \cdot \mathbf{f}(ih, \mathbf{Q}_n^i, \mathbf{V}_n^{i-1}), $$

$$ \mathbf{V}_n^i = -e \mathbf{V}_n^{i-1} + \mathrm{Proj}_{\mathbf{M}(\mathbf{Q}_n^i)} (\mathbf{V}_n^i + e \mathbf{V}_n^{i-1}; \mathcal{V}(\mathbf{Q}_n^i)) $$

$$ - \mathbf{q}_n(t) = \mathbf{Q}_n^{i-1} + (t - (i-1)h)\mathbf{V}_n^{i-1}, \quad \forall t \in [(i-1)h, ih], $$

Now, suppose that the data $\varphi_\alpha$, $M$ and $f$ in problem $P_3$ are analytic functions of their arguments and that the impact law is chosen to be that of Moreau (4) for an arbitrary $e \in [0, 1]$. Suppose, in addition, that $f$ is Lipschitz with respect to $(\mathbf{q}, \dot{\mathbf{q}})$ and:

$$ \forall \mathbf{q}, \quad \forall \alpha \in J(\mathbf{q}) \quad d\varphi_{\alpha, \mathbf{q}} \cdot \mathbf{M}^{-1}(\mathbf{q}) \cdot d\varphi_{\beta, \mathbf{q}} = 0, \quad \text{if } e \in ]0, 1], \qquad (6) $$

$$ \forall \mathbf{q}, \quad \forall \alpha \in J(\mathbf{q}) \quad d\varphi_{\alpha, \mathbf{q}} \cdot \mathbf{M}^{-1}(\mathbf{q}) \cdot d\varphi_{\beta, \mathbf{q}} \le 0, \quad \text{if } e = 0 $$

then, Paoli proves in [6] that the sequence $\mathbf{q}_n$ converges strongly in $W^{1,1}(0, T)$ towards the unique solution $\mathbf{q}$ to problem $P_3$. Actually, her proof seems to cover also the situation where the conditions (6) are not necessarily fulfilled for all $\mathbf{q}$, but only for all the locations $\mathbf{q}(t)$ ($t \in [0, T]$) associated with the unique solution $\mathbf{q}$ of problem $P_3$.

An alternate time-stepping scheme is proposed in [5] with proof of the same results of convergence.

# 5. Frictional unilateral multibody dynamics

In most use of unilateral multibody dynamics in engineering problems, it is desirable to account for (dry) friction. Unfortunately, frictional unilateral multibody dynamics has not reached yet the same stage of completion as in the frictionless case.

## 5.1. The need for a general formulation

It is very surprising that the first attempt to obtain a general formulation of an evolution problem associated with frictional unilateral multibody dynamics was only that of Charles in his PhD thesis [14].

Formerly, the general and rather imprecise idea that was prevailing is that one should use the Lagrange equation for dynamics, complemented with the Coulomb friction law applying to real-world reaction forces and to the velocities of the material contact points. This paradigm has not led so far to any general formulation of frictional multibody dynamics. It has only resulted in two outcomes:
- in very specific cases, there is a one-to-one correspondence between generalized reaction forces and real-world reaction forces, and this paradigm can then be turned into a sound mathematical evolution problem. One such example is the Painlevé example of a rigid bar above a rigid obstacle, so that there can only be, at most, one material point at which the contact takes place. Such a precise formulation of the evolution problem can be found in Stewart [13] in the case of the Painlevé example. The motivation of Stewart was to prove the existence of a dynamical solution for the Painlevé example (which he did), since suspicion of non-existence of solution had fed a one-century-old controversy which could be entitled: ‘is the (local) Coulomb law consistent with rigid body dynamics?’;

- Numerical methods to compute discrete-time solutions in more complex cases. Examples of such methods are those of Jean [17], which was inspired by the time-stepping technique of Monteiro Marques, and of Stewart in [13]. A common feature of all these methods is that they are derived on the basis of the paradigm of the Coulomb law applying to reaction forces in the real world. But, they are not derived on the basis of a general continuous-time evolution problem. This was recognized by Stewart in [18]: ‘In many ways it is easier to write down a numerical method for rigid-body dynamics than it is to say exactly what the method is trying to compute.’

As was already pointed out when formulating the evolution problem associated with frictionless multibody dynamics, the equation of motion involves the generalized reaction force which is therefore completely determined, when the motion is calculated, whereas the precise distribution of the real-world reaction forces remains undetermined, in general. Besides, the modern point of view in mechanics stipulates, through the systematic use of the Principle of Virtual Power, that the appropriate representation of forces within any mechanical theory is the one that comes from duality with velocities. Within the dynamics of rigid bodies, the appropriate representation of forces is therefore that of Lagrange's generalized forces. Hence, it is not surprising that sticking to a formulation of the friction law in terms of the real-world reaction forces has not yet come to a successful end, with respect to the issue of obtaining a general consistent formulation of the dynamics. In the specific cases where it is possible though, such a paradigm is known to raise the following paradoxes, once again unsurprisingly.

- *Painlevé paradox.* In the case of the planar dynamics of a rigid bar in a half-plane with possible frictional contact at one extremity of the bar, there is a one-to-one correspondence between generalized reaction force and real-world reaction force. It is therefore possible to obtain a precise formulation of an evolution problem in that situation (it is to be found in [13]). However, it was recognized as early as the nineteenth century that multiple solutions can be encountered [19]. The example of the rigid bar would therefore be different in nature to that of a point particle for which the uniqueness of the solution can be proved rigorously (see the discussion in section 3). At the time of Painlevé, there was also a suspicion of possible non-existence of a solution. Later, Lecornu [20] pointed out that the existence of a solution could be recovered by allowing an impact (velocity jump and associated Dirac mass of the reaction measure), although the incoming normal velocity vanishes. This was further discussed by Moreau [21], in the framework of the motions with measure acceleration, under the name ‘tangential impact’ or ‘frictional catastrophe’. Finally, Stewart proved in [13] an existence result for the Painlevé example with completely inelastic impact law, showing that the allowance for ‘tangential impacts’ was ruling out the suspicion of non-existence in the Painlevé example. Hence, the analysis of Stewart, which was restricted to the system with one material contact point only, solved the non-existence issue in the Painlevé paradox, but let open the non-uniqueness issue.

- *Kane paradox.* The non-uniqueness part of the Painlevé paradox was sometimes claimed to be acceptable. However, another paradox connected to the use of a friction law applying to real-world reaction forces arose more recently in [22,23], which was more serious. The example is that of a double pendulum whose extremity can experiment frictional contact with a straight obstacle. The Coulomb law and an impact law with a restitution coefficient of 0.5 are adopted. In that example also, there is an exceptional one-to-one correspondence between generalized and real-world reaction forces, making it possible to formulate an evolution problem associated with the Coulomb law applying to real-world reaction forces. A motion solving that the evolution problem is exhibited, in which the kinetic energy of the system is increased during an impact, due to friction. This is an energetic inconsistency raised by the use of the Coulomb law applying to real-world reaction forces within rigid body dynamics. In the case of the Painlevé example with completely inelastic impact law, Stewart proves in [13] that the solution he constructs is dissipative, and therefore does not exhibit this energetic inconsistency. However, since no uniqueness is to be expected, there is no guarantee of energetic consistency of all the solutions.

## 5.2. The case of one unilateral constraint

In this section, we shall obtain the evolution problem governing the frictional dynamics of a discrete system obeying only one unilateral constraint of the form:

$$ \varphi(\mathbf{q}) \geq 0 $$

According to the discussion in the preceding section, we are going to resist invoking any real-world reaction force and therefore stick to generalized reaction forces and velocities.
### 5.2.1. Splitting into normal and tangential parts

At a configuration where the constraint is active, it is easy to see how to define what a tangential velocity is: it is simply a velocity that arises in a motion maintaining the contact active:

$$ \langle d\varphi_q, \dot{q} \rangle = 0 $$

However, the issue of defining what a normal velocity is turns out to be trickier, as it requires the choice of a scalar product. Even if the configuration manifold can locally be identified with $\mathbb{R}^d$ by making the choice of a parameterization, the use of the canonical scalar product of $\mathbb{R}^d$ is not an option, for, in that case, the definition of what a normal velocity is would depend on the choice of the parameterization. Actually, there is only one natural scalar product that is intrinsic, that is, it does not depend on the choice of a parameterization: the one induced by the kinetic energy whose matrix in a given parameterization is the (positive definite) mass matrix $M(q)$. We shall see later that this inevitable choice of a scalar product has the fallout of enforcing energetic consistency and therefore protecting from Kane's paradox. A normal velocity $v_n$ is therefore one that is orthogonal to a tangential velocity according to the scalar product $M(q)$:

$$ \lambda \in \mathbb{R}, \quad \mathbf{v}_n = \lambda \mathbf{M}^{-1}(\mathbf{q}) \cdot d\varphi_q = \lambda \nabla_q \varphi $$

where the scalar product $M(q)$ enables, as usual, to make the link between the gradient and the differential of the function $\varphi$. Likewise, a normal generalized reaction force $r_n$ is one that develops zero virtual power in any tangential velocity:

$$ \lambda \in \mathbb{R}, \quad \mathbf{r}_n = \lambda d\varphi_q $$

and a tangential reaction force $r_t$ is one that develops zero virtual power in any normal velocity:

$$ \mathbf{r}_t \cdot \mathbf{M}^{-1}(\mathbf{q}) \cdot d\varphi_q = 0 $$

Some insight is gained from an (n elementary) geometric point of view: velocities live in the tangent space $T_q Q$, which is endowed with an Euclidean structure by the scalar product (mass matrix) $M(q)$. Generalized forces are linear forms applying on velocities through the Virtual Power paradigm, and therefore live in the dual space $T_q^* Q$, which is endowed with a natural Euclidean structure by the dual scalar product $\mathbf{M}^{-1}(\mathbf{q})$. Hence, the natural split into tangential and normal parts reads as:

$$ \forall \mathbf{v} \in T_q Q, \quad \mathbf{v} = \mathbf{v}_t + v_n \frac{\nabla_q \varphi}{|\nabla_q \varphi|_{\mathbf{q}}}, \quad \text{with } (\mathbf{v}_t, \nabla_q \varphi)_{\mathbf{q}} = 0, \text{ and } |\nabla_q \varphi|_{\mathbf{q}} = \sqrt{d\varphi_q \cdot \mathbf{M}^{-1}(\mathbf{q}) \cdot d\varphi_q}, $$

$$ \forall \mathbf{r} \in T_{\mathbf{q}}^* Q, \quad \mathbf{r} = \mathbf{r}_t + r_n \frac{d\varphi_q}{|d\varphi_q|_{\mathbf{q}}^*}, \quad \text{with } (\mathbf{r}_t, d\varphi_q)_{\mathbf{q}}^* = 0, \text{ and } |d\varphi_q|_{\mathbf{q}}^* = \sqrt{d\varphi_q \cdot \mathbf{M}^{-1}(\mathbf{q}) \cdot d\varphi_q} $$

### 5.2.2. The generalized friction law

We are now going to rely on the splitting of generalized reaction forces and velocities into normal and tangential parts to infer the general form of a friction law expressed in terms of these quantities.

Let us first recall some classical definitions of convex analysis. Let $C$ be a nonempty closed convex subset of $\mathbb{R}^d$. Its indicatrix function $I_C$ is defined by:

$$ I_C(\mathbf{x}) = \begin{cases} 0 & \text{if } \mathbf{x} \in C \\ +\infty & \text{if } \mathbf{x} \notin C \end{cases} $$

It is proper, lower semi-continuous and convex. Its conjugate function (by the Legendre-Fenchel transform) is the support function $S_C$ of $C$:

$$ S_C(\mathbf{y}) = \sup_{\mathbf{x} \in C} \langle \mathbf{x}, \mathbf{y} \rangle $$

Introducing the closed unit disk $\mathcal{B}$ of $\mathbb{R}^2$ (with respect to the canonical Euclidean norm $|\cdot|$), it is readily checked that:

$$ S_{\mathcal{B}}(\mathbf{y}) = |\mathbf{y}| $$

The pointwise formulation (2) of the Coulomb law is therefore equivalent to anyone of the following equivalent statements:

(i) $-\mathbf{v}_t \in \partial I_{\mu r_n \mathcal{B}}[\mathbf{r}_t]$

(ii) $\mathbf{r}_t \in \partial S_{\mu r_n \mathcal{B}}[-\mathbf{v}_t]$

(iii) $\forall \hat{\mathbf{v}} \in \mathbb{R}^2, \quad (\mathbf{r}_t, \hat{\mathbf{v}} - \mathbf{v}_t) + \mu r_n (|\hat{\mathbf{v}}| - |\mathbf{v}_t|) \ge 0$
where $\partial f$ denotes the subdifferential of the proper, lower semi-continuous function $f$. Since $\partial I_C$ is the closed convex cone of all the outward normal vectors to $C$, statement (i) is also called the normality rule, when $\mu r_n$ is assumed to be given as in Coulomb’s experiments (where $r_n$ was merely the opposite of the gravity force). The equivalent statement (ii) is often referred to as the maximum dissipation principle.

To generalize the Coulomb friction law (2) for the point particle to the case of an arbitrary discrete mechanical system, we are naturally led to postulate the existence of a set $C_q^*$ in the space $T_q^* \subset T_q^*Q$ of all tangential reaction force, such that all the admissible tangential components of the generalized reaction force are characterized by the condition:

$$ \mathbf{r}_t \in r_n C_q^* $$

where the normal reaction force $r_n$ is required to be nonnegative and vanishes when the constraint is not active (that is, when $\varphi(\mathbf{q}) > 0$). As for the Coulomb friction law (2), this condition entails that whenever the normal reaction force vanishes (for example in motion episodes where the constraint is not active), then the tangential reaction also vanishes.

In general, the set $C_q^*$ will not be required to be a disk (case of isotropic friction (2) for a point particle), but only to be nonempty. It is also required to be closed and convex so that the normality rule can be adopted as a flow rule: this is nothing but assuming that for a given value of the reaction normal component, the flow rule should obey the maximum dissipation principle. With this assumption, the flow rule is formally expressed by anyone of the following equivalent statements:

$$ (i) \quad -\dot{\mathbf{q}}_t^+ \in \partial I_{r_n C_q^*} [\mathbf{r}_t] $$

$$ (ii) \quad \mathbf{r}_t \in \partial S_{r_n C_q^*} [-\dot{\mathbf{q}}_t^+] $$

$$ (iii) \quad \forall \hat{\mathbf{v}} \in \mathcal{T}_{\mathbf{q}}, \quad \langle \mathbf{r}_t, \hat{\mathbf{v}} - \dot{\mathbf{q}}_t^+ \rangle + r_n [S_{C_q^*}(-\hat{\mathbf{v}}) - S_{C_q^*}(-\dot{\mathbf{q}}_t^+)] \ge 0 $$

where the identity $S_{r_n C_q^*} = r_n S_{C_q^*}$ has been used. Here, the subdifferentials should be understood in the sense of the duality between the space $\mathcal{T}_{\mathbf{q}}$ of tangential generalized velocities and the space $\mathcal{T}_{\mathbf{q}}^*$ of tangential generalized forces.

Actually, the set $C_q^*$ can be viewed either as a subset of the space $\mathcal{T}_{\mathbf{q}}^*$ of tangential generalized forces or as a subset of the space $T_{\mathbf{q}}^*Q$ of all generalized forces. In the latter case, one has:

$$ \forall \hat{\mathbf{v}} \in T_{\mathbf{q}} Q, \quad S_{C_q^*}(\hat{\mathbf{v}}) = S_{C_q^*}(\hat{\mathbf{v}}_t) $$

since $C_q^*$ is contained in the space $\mathcal{T}_{\mathbf{q}}^*$ of tangential generalized forces. Therefore, statement (iii) is equivalent to:

$$ (iii') \quad \forall \hat{\mathbf{v}} \in T_{\mathbf{q}} Q, \quad \langle \mathbf{r}_t, \hat{\mathbf{v}} - \dot{\mathbf{q}}_t^+ \rangle + r_n [S_{C_q^*}(-\hat{\mathbf{v}}) - S_{C_q^*}(-\dot{\mathbf{q}}_t^+)] \ge 0 $$

where we recall that $T_{\mathbf{q}}Q$ stands for the space of all generalized velocities at the generalized configuration $\mathbf{q}$ (tangent space at $\mathbf{q}$ to the configuration manifold). This latter form of the (generalized) friction law is the more appropriate one in view of deriving an expression of the friction law encompassing the episodes of smooth motion and the impacts, that is, in situations where the reaction force is not pointwise defined but is only a measure (with respect to time). Having these situations in mind, the friction law has to be postulated under the form:

$$ \forall \hat{\mathbf{v}} \in C^0([0, T]; TQ), \quad \text{with } \hat{\mathbf{v}}(t) \in T_{\mathbf{q}}(t)Q \qquad \int_{[0,T]} \langle \mathbf{r}_t, \hat{\mathbf{v}} - \dot{\mathbf{q}}_t^+ \rangle + r_n [S_{C_q^*}(-\hat{\mathbf{v}}) - S_{C_q^*}(-\dot{\mathbf{q}}_t^+)] \ge 0 \quad (7) $$

As in the case of the ordinary Coulomb law where $C_q^*$ is a closed disk, the generalized friction law (7) makes sense whenever the generalized reaction force $\mathbf{r}$ is a measure and the generalized right-velocity $\dot{\mathbf{q}}_t^+$ is a function of time with bounded variation.

Equation (7) provides the general form of what should be a friction law formulated in terms of generalized reaction force, applying to both smooth episodes of motion and impacts, which is, in addition, energetically consistent. This framework cannot be avoided. Of course, the question remains to know what choice should be made in practice for the friction set $C_q^*$. It turns out that a systematic proposal of $C_q^*$ can always be made on the basis of the paradigm of a pointwise Coulomb cone containing hypothetical real-world reaction forces and, then, of the use of formulae (3) to build the cone containing the generalized reaction force $\mathbf{r}$, and therefore to deduce a proposal for $C_q^*$ (a detailed account of this systematic construction is to be found in [14]).

In the above-mentioned systematic construction for the sets $C_q^*$, it turns out that unbounded sets $C_q^*$ can appear. In that case, the function $S_{C_q^*}$ takes the value $+\infty$, that is, its domain $D(S_{C_q^*})$ is smaller than the entire space $T_{\mathbf{q}}Q \simeq \mathbb{R}^d$. Since the generalized friction law (7) entails that the tangential right-velocity $\dot{\mathbf{q}}_t^+$ must be in the closure $\bar{D}(S_{C_q^*})$, it turns out that unbounded friction sets induce forbidden values for the right-velocity $\dot{\mathbf{q}}_t^+$ at those instants belonging to Supp $\mathbf{r}$. Another significant difference of unbounded friction sets with bounded friction sets is that they allow atoms of tangential reaction
force without any atom<fnref n="1" /> in the normal force, which is impossible in the case of a bounded friction set. Typically, this is going to happen when the system **q** is arriving on the obstacle with zero normal velocity and along a forbidden value of the tangential velocity. In that case, a *tangential impact* (discontinuity of tangential velocities at an instant where the normal velocity is continuous) will be possible.

### 5.2.3. The evolution problem

Gathering all the preceding considerations, we are in measure to obtain the general formulation of the frictional dynamics of an arbitrary discrete mechanical system with only one unilateral constraint. As previously, we are given an initial condition ($\mathbf{q}_0, \mathbf{v}_0) \in TQ \simeq \mathbb{R}^d \times \mathbb{R}^d$ compatible with the obstacle:

$$ \varphi(\mathbf{q}_0) \ge 0 \quad \text{and} \quad \varphi(\mathbf{q}_0) = 0 \Rightarrow (\mathbf{v}_0, d\varphi_{\mathbf{q}_0}) \ge 0, \text{ and } \mathbf{v}_0 \in \overline{D(S_{C_{\mathbf{q}_0}^*})} $$

(the last requirement is a restriction only in the case where $C_{\mathbf{q}_0}^*$ is unbounded and aims at protecting from tangential impact at the initial instant).

**Problem P₄.** Find $\mathbf{q} \in \text{MMA}([0, T]; \mathbb{R}^d)$ and $\mathbf{r} \in \mathcal{M}([0, T]; \mathbb{R}^d)$ such that:

$$ - \mathbf{q}(0) = \mathbf{q}_0, \quad \dot{\mathbf{q}}^+(0) = \mathbf{v}_0, $$

$$ - \mathbf{M}(\mathbf{q}) \cdot \frac{D\dot{\mathbf{q}}}{Dt} = \mathbf{f}(t; \mathbf{q}, \dot{\mathbf{q}}^-) + \mathbf{r}, $$

$$ - \forall t \in [0, T], \quad \varphi(\mathbf{q}(t)) \ge 0, $$

$$ - \text{Supp} \mathbf{r} \subset \{t \in [0, T] \mid \varphi(\mathbf{q}(t)) = 0\}, $$

$$ - r_n \ge 0, $$

$$ - \forall \hat{\mathbf{v}} \in C^0([0, T]; TQ), \text{ with } \hat{\mathbf{v}}(t) \in T_{\mathbf{q}(t)} Q $$

$$ \int_{[0,T]} \{\mathbf{r}_t, \hat{\mathbf{v}} - \dot{\mathbf{q}}^+\} + r_n[S_{C_{\mathbf{q}}^*}(-\hat{\mathbf{v}}) - S_{C_{\mathbf{q}}^*}(-\dot{\mathbf{q}}^+)] \ge 0, $$

$$ - \varphi(\mathbf{q}(t)) = 0 \quad \Rightarrow \quad \dot{q}_{\mathrm{n}}^{+}(t) = -e\dot{q}_{\mathrm{n}}^{-}(t). $$

The following results about problem P₄ have been proved in [14].

- If $e \in [0, 1]$, then any solution to problem P₄ is dissipative, that is, for all $t_1, t_2 \in [0, T]$, such that $t_1 \le t_2$, one has:

$$ \frac{1}{2} \dot{\mathbf{q}}^{+}(t_2) \cdot \mathbf{M}(\mathbf{q}(t_2)) \cdot \dot{\mathbf{q}}^{+}(t_2) \le \frac{1}{2} \dot{\mathbf{q}}^{-}(t_1) \cdot \mathbf{M}(\mathbf{q}(t_1)) \cdot \dot{\mathbf{q}}^{-}(t_1) + \int_{t_1}^{t_2} \left\langle \mathbf{f}(t; \mathbf{q}(t), \dot{\mathbf{q}}^{-}(t)), \dot{\mathbf{q}}^{-}(t) \right\rangle dt \quad (8) $$

Problem P₄ is the general evolution problem raised by the frictional dynamics of an arbitrary discrete system submitted to one unilateral constraint. The impact law governs the normal component of the velocity, as in the frictionless situation. The generalized friction law rules both the impacts and the smooth episodes of motion. The above result shows that any solution is energetically consistent. Let us just recall that the example of Kane in [22,23] was that of a discrete system (double pendulum), with one unilateral constraint. However, he did not know such a formulation as our problem P₄, and used the usual Coulomb friction law applying to real-world reaction forces. As a result, some solutions were possible, in which an increase in kinetic energy during an impact was possible. The above result is therefore a demonstration of superiority of formulating the dynamics along problem P₄, in which *energetic consistency is built in*.

- The existence and the uniqueness of a solution to problem P₄ are still to be proved (under the assumption that the data are analytic). Note, however, that this existence and uniqueness result has already been proved in the case of problem P₂, which is nothing but a particular case of problem P₄, for which the friction set $C_{\mathbf{q}}^*$ reduces to a disk. In addition, a discussion around the Painlevé example is developed in [14]. When the dynamics of the Painlevé system is formulated on the basis of the usual Coulomb friction law applying to real-world reaction forces (as in [13]), multiple solutions can be exhibited for large enough friction coefficients. However, when the dynamics of the same system is formulated along problem P₄, among these explicit multiple possible dynamical evolutions, only one of them remains a solution to the corresponding problem P₄. This is a clue that the formulation of the dynamics along problem P₄ escapes not only from energetic inconsistencies as in the Kane paradox, but also from indeterminacies and multiple solutions as in the Painlevé paradox.

<footnote n="1">¹ An instant $t$ is said to be an atom of the Radon measure $\mathbf{r}$, if $\mathbf{r}(\{t\}) = \int_{\{t\}} \mathbf{r} \neq \mathbf{0}$.</footnote>
5.3. The general case

There remains to derive a systematic formulation of frictional dynamics of discrete systems, in the case where several unilateral constraints apply: $\varphi_{\alpha}(\mathbf{q}) \ge 0$ ($\alpha = 1, 2, \dots, n$). This is nothing but an extension of frictionless multibody dynamics along problem $\mathcal{P}_3$, in the spirit of problem $\mathcal{P}_4$. Such an extension is presented with details in [14].

The analysis introduces a collection $C_{\alpha,q}^*$ ($\alpha = 1, 2, \dots, n$) of friction sets, each of them being associated with one of the unilateral constraints. The corresponding generalized friction law involves the sum of each associated dissipation:

$$S_{C_{\alpha,q}^*}(-\dot{\mathbf{q}}^+)$$

A corresponding evolution problem is precisely stated and it is rigorously proved that any solution is dissipative, that is, fulfils inequality (8).

Hence, the point of view of a generalized friction law applying to generalized reaction forces (instead of the usual point of view of a friction law expressed in terms of the real-world reaction forces) yields the following pleasant fallout. It enables to formulate for the first time, a consistent abstract evolution problem associated with frictional unilateral multibody dynamics, in the most general situation. And, in addition, the energetic consistency of such a formulation is built in, meaning that any solution is necessarily dissipative.

References

[1] M. Schatzman, A class of nonlinear differential equations of second order in time, Nonlinear Anal. 2 (2) (1978) 355–373.

[2] J.-J. Moreau, Standard inelastic shocks and the dynamics of unilateral constraints, in: G. Del Piero, F. Maceri (Eds.), Unilateral Problems in Structural Analysis, Springer-Verlag, Wien, New York, 1983, pp. 173–221.

[3] L. Paoli, M. Schatzman, Mouvement à un nombre fini de degrés de liberté avec contraintes unilatérales: cas avec perte d'énergie, Math. Model. Numer. Anal. 27 (6) (1993) 673–717.

[4] M.D.P. Monteiro Marques, Differential Inclusions in Nonsmooth Mechanical Problems – Shocks and Dry Friction, Progress in Nonlinear Differential Equations and Their Applications, vol. 9, Birkhäuser, Basel, Boston, Berlin, 1993.

[5] L. Paoli, Time-stepping approximation of rigid-body dynamics with perfect unilateral constraints. I and II, Arch. Ration. Mech. Anal. 198 (2010) 457–503, 505–568.

[6] L. Paoli, A proximal-like algorithm for vibro-impact problems with a non smooth set of constraints, J. Differ. Equ. 250 (2011) 476–514.

[7] D. Percivale, Uniqueness in the elastic bounce problem, I, J. Differ. Equ. 2 (2) (1985) 206–215.

[8] D. Percivale, Uniqueness in the elastic bounce problem, II, J. Differ. Equ. 90 (1991) 304–315.

[9] P. Ballard, The dynamics of discrete mechanical systems with perfect unilateral constraints, Arch. Ration. Mech. Anal. 154 (2000) 199–274.

[10] P. Ballard, Dynamics of rigid bodies systems with unilateral or frictional constraints, in: D.Y. Gao, R. Ogden (Eds.), Advances in Mechanics and Mathematics, vol. 1, Kluwer Academic Publishers, 2002, pp. 3–88, chap. 1.

[11] P. Ballard, S. Basseville, Existence and uniqueness for dynamical unilateral contact with Coulomb friction: a model problem, Math. Model. Numer. Anal. 39 (2005) 59–77.

[12] A. Charles, P. Ballard, Existence and uniqueness of solutions to dynamical unilateral contact problems with Coulomb friction: the case of a collection of points, Math. Model. Numer. Anal. 48 (1) (2014) 1–25.

[13] D.E. Stewart, Convergence of a time-stepping scheme for rigid-body dynamics and resolution of Painlevé’s problem, Arch. Ration. Mech. Anal. 145 (1998) 215–260.

[14] A. Charles, P. Ballard, The formulation of dynamical contact problems with friction in the case of systems of rigid bodies and general discrete mechanical systems–Painlevé and Kane paradoxes revisited, Z. Angew. Math. Phys. (2016) 67–99.

[15] A. Bressan, Incompatibilità dei Teoremi di Esistenza e di Unicità del Moto per un Tipo molto Comune e Regolare di Sistemi Meccanici, Ann. Sc. Norm. Super. Pisa (3) XIV (1960) 333–348.

[16] L. Paoli, Continuous dependence on data for vibro-impact problems, Math. Models Methods Appl. Sci. 15 (1) (2005) 53–93.

[17] M. Jean, The non-smooth contact dynamics method, Comput. Methods Appl. Mech. Eng. 177 (1999) 235–257.

[18] D.E. Stewart, Rigid-body dynamics with friction and impact, SIAM Rev. 42 (1) (2000) 3–39.

[19] P. Painlevé, Sur les lois du frottement de glissement, C. r. hebd. Séances Acad. Sci. Paris 121 (1895) 112–115.

[20] L. Lecornu, Sur la loi de Coulomb, C. r. hebd. Séances Acad. Sci. Paris 140 (1905) 847–848.

[21] J.-J. Moreau, Unilateral contact and dry friction in finite freedom dynamics, in: J.-J. Moreau, P.D. Panagiotopoulos (Eds.), Non-smooth Mechanics and Applications, Springer-Verlag, Wien, New York, 1988, pp. 1–82.

[22] T.R. Kane, A dynamics puzzle, Stanf. Mech. Alumni Club Newsl. (1983) 6.

[23] T.R. Kane, D.A. Levinson, Dynamics: Theory and Applications, McGraw-Hill, New York, 1985.