# Hamilton's Principle as Variational Inequality for Mechanical Systems with Impact

Remco Leine, Ueli Aeberhard, Christoph Glocker

► To cite this version:

Remco Leine, Ueli Aeberhard, Christoph Glocker. Hamilton's Principle as Variational Inequality for Mechanical Systems with Impact. Journal of Nonlinear Science, 2009, 19 (6), pp. 633-664. 10.1007/s00332-009-9048-z. hal-01352880

HAL Id: hal-01352880

https://hal.science/hal-01352880

Submitted on 10 Aug 2016

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
# Hamilton's Principle as Variational Inequality for Mechanical Systems with Impact

R.I. Leine

U. Aeberhard

C. Glocker

**Abstract** The classical form of Hamilton's principle holds for conservative systems with perfect bilateral constraints. Several attempts have been made in literature to generalise Hamilton's principle for mechanical systems with perfect unilateral constraints involving impulsive motion. This has led to a number of different variants of Hamilton's principle, some expressed as variational inequalities. Up to now, the connection between these different principles has been missing. The aim of this paper is to put these different principles of Hamilton in a unified framework by using the concept of weak and strong extrema. The difference between weak and strong variations of the motion is explained in detail. Each type of variation leads to a variant of the principle of Hamilton in the form of a variational inequality. The conclusion of the paper is that each type of variation leads to different necessary and sufficient conditions on the impact law. The principle of Hamilton with strong variations is valid for perfect unilateral constraints with a completely elastic impact law, whereas the weak form of Hamilton's principle only requires perfect unilateral constraints and no condition on the energy.

**Keywords** Unilateral constraint · Principle of d'Alembert–Lagrange · Weierstrass–Erdmann corner conditions · Non-smooth dynamics · Contact
# 1 Introduction

Classical analytical mechanics is concerned with mechanical systems with perfect bilateral (mostly holonomic) constraints (Papastavridis 2002) and is closely related with the calculus of variations as many principles of statics and dynamics are formulated in terms of variational problems, e.g. the principle of virtual work and Hamilton's principle (see, for instance, Lanczos 1962). Unilateral constraints, which are basically inequality constraints, are completely ignored in classical analytical mechanics because inequalities are not discussed by the classical calculus of variations and, also, because unilateral constraints in dynamics lead to shocks with discontinuities in the velocity. The mathematical tools to handle unilateral problems in statics, and later also in dynamics, have only been developed since the last 4 decades. The field of non-smooth dynamics is now rapidly developing. We refer the reader to the available textbooks on this subject (Leine and Nijmeijer 2004; Leine and van de Wouw 2008; Glocker 2001; Brogliato 1999).

Variational problems involving inequality constraints can be described by variational inequalities and were first introduced by Hartman and Stampacchia (1966) to study partial differential equations. The applicability of the theory has since been expanded to include problems from mechanics, finance, optimisation, and game theory. References on variational inequalities can be found in the standard reference (Kinderlehrer and Stampacchia 1980) or in the more recent books (Goeleven et al. 2003a, 2003b; Kravchuk and Neittaanmäki 2007). Variational problems with non-convex inequality constraints are described by hemi-variational inequalities and are discussed in the seminal work of Panagiotopoulos (1993). Moreover, (hemi-)variational inequalities are related to (non-)convex super-potentials through the sub-derivative known from non-smooth analysis; see the work of Moreau (1968) and Panagiotopoulos (1981). In the following, we will simply use the term 'variational inequalities' to denote (hemi-)variational inequalities.

Various authors study the classical principle of Hamilton, which is stated in the form of a variational equality, in the context of mechanical systems with impact. In Kozlov and Treshchëv (1991), Fetecau et al. (2003) and Cirak and West (2005), it is stated that the principle of Hamilton in its classical form as variational equality is also valid for motion with completely elastic impact (see also the review of Brogliato 1999). In the proof of Kozlov and Treshchëv (1991), it becomes clear that in this case only comparison curves are considered which have an impact, i.e. the variations of the position functions are not completely free. The classical principle of Hamilton is used in Fetecau et al. (2003) and Cirak and West (2005) to develop a symplectic-momentum preserving numerical scheme for systems with impact. The principle of Hamilton in equality form used in Kozlov and Treshchëv (1991), Fetecau et al. (2003) and Cirak and West (2005) does not hold for a system with persistent contact, i.e. a contact being closed during a time lapse.

The extension of classical analytical mechanics to perfect unilateral constraints, with possible persistent contact, asks for a reformulation of the variational principles of mechanics in terms of variational inequalities. The principle of d'Alembert-Lagrange in the form of a variational inequality has been discussed in Panagiotopoulos and Glocker (1998, 2000), Goeleven et al. (1997, 1999) and May (1984a, 1984b).
Various forms of the principle of Hamilton as variational inequality can be derived from the principle of d'Alembert-Lagrange as variational inequality (May 1984a, 1984b; Panagiotopoulos and Glocker 1998, 2000). The principle of Hamilton in May (1984a, 1984b) is valid for frictionless unilateral constraints (i.e. perfect unilateral constraints) and non-impulsive motion. Hence, impacts which naturally arise due to collisions if the motion is unilaterally constrained are not considered in May (1984a, 1984b). In Panagiotopoulos and Glocker (1998), a principle of Hamilton is presented (Theorem 3) which is valid for perfect unilateral constraints with impulsive motion. Another principle of Hamilton is presented in Theorem 4 of Panagiotopoulos and Glocker (1998), which is basically identical to Proposition 4 in Panagiotopoulos and Glocker (2000), and which is valid for perfect unilateral constraints with a completely elastic impact law.

A very different approach has been taken by Tornambè (1999) and Percivale (1985, 1991) in which the action integral is modified with a Lagrangian multiplier term and the principle of Hamilton is stated as variational equality. In Percivale (1985, 1991), it is proven that a trajectory with a completely elastic impact law is a stationary point of the modified action integral. In Tornambè (1999), the inequality constraints are transformed into equality constraints by making use of so-called Valentine variables. The equality constraint on the Valentine variables is subsequently incorporated in the action integral with a Lagrangian multiplier and the stationarity of the resulting functional is shown to lead to solutions of the completely elastic impact problem. The principles of Hamilton in Tornambè (1999) and Percivale (1985, 1991) are variational equalities using a modified action integral.

The great difference between the approach used by May (1984a, 1984b) and Panagiotopoulos and Glocker (1998, 2000) and the approach used by Tornambè (1999) and Percivale (1985, 1991) is that the latter approach does not derive the principle of Hamilton from the principle of d'Alembert-Lagrange and the principle of virtual work. In the current paper, we will derive various forms of the principle of Hamilton starting from the principle of virtual work. We will therefore follow the first approach which leads to variational inequalities.

The mathematical literature on the calculus of variations distinguishes between so-called weak and strong local extrema of a functional $J(y) = \int_I f(x, y, y') dx$ (Cesar 1984; Cesari 1983; Troutman 1996). The strong norm measures the maximal function value $\max(|y(x)|)$ on the interval $I$, whereas the weak norm also takes the derivative $y'$ into account. A weak/strong local extremal is extremal with respect to comparison functions within a neighbourhood defined by the weak/strong norm. It is known that a weak local extremal fulfills the Euler-Lagrange equations almost everywhere and fulfills the first Weierstrass-Erdmann corner condition, i.e. the left and right limit of $f_{y'}$ agree on its corner points. A strong local extremal fulfills in addition the second Weierstrass-Erdmann corner condition, i.e. the left and right limit of $f - y'f_{y'}$ agree on its corner points. There is some confusion about this in literature. Some authors incorrectly state that a weak local extremal necessarily fulfills the second Weierstrass-Erdmann corner condition (see Cesar 1984 for a discussion and counter examples). Moreover, we note that the first Weierstrass-Erdmann corner condition generalises to an inclusion if it is applied to a variational inequality. Furthermore, in our knowledge, Hamilton's principle has not been put within the concepts of weak and strong local extrema.
The current paper derives weak and strong forms of the principle of Hamilton as variational inequality using the concepts of weak and strong extrema (Cesar 1984; Troutman 1996) and concepts from measure and integration theory. It therefore becomes much more clear how the existing forms of the principle of Hamilton, which have been derived in Panagiotopoulos and Glocker (1998, 2000), Tornambè (1999) and Percivale (1985, 1991), have to be understood. Moreover, we demonstrate how the concept of virtual displacements and the principle of virtual work are to be used in a non-smooth setting and that we have in fact to distinguish between weak and strong virtual displacements. The general variation of the action integral leads to a strong form of Hamilton's principle and to assumptions which are generalisations of the classical Weierstrass–Erdmann corner conditions.

The principle of Hamilton occupies a central role in classical analytical mechanics as it relates dynamics to the calculus of variations. From the practical point of view, Hamilton's principle is used in mechanics to derive (numerical) approximation methods, such as the Ritz method. The advantages of the principle of Hamilton are its concise formulation and the fact that approximations which maintain the symmetries of the Hamiltonian automatically preserve the corresponding conservation laws (Salomon 1983). Similarly, the generalisations of Hamilton's principle, discussed in this paper, can be used to derive difference schemes which can form the basis for dedicated numerical integration schemes for impulsive mechanical systems.

The paper is organised as follows. In Sect. 2, we give the reader an introduction to cones from convex analysis in functional spaces and set the notation. Variational inequalities are discussed in Sect. 3 and special attention is paid to weak and strong local extrema. Subsequently, we discuss in Sect. 4 the principle of d'Alembert-Lagrange in inequality form, which constitutes a variational inequality. The principle of Hamilton in inequality form is derived from the principle of d'Alembert-Lagrange in Sect. 5 for non-impulsive motion under the supposition of perfect geometric unilateral constraints. In Sect. 6, the general variation of the action integral is used to derive necessary and sufficient conditions for the validity of a strong form of Hamilton's principle in inequality form for impulsive motion. Concepts of measure and integration theory are used in Sect. 7 to derive a weak form of the principle of Hamilton. It turns out that only the supposition of perfect geometric unilateral constraints is necessary and sufficient for the validity of a weak form of Hamilton's principle as variational inequality. Finally, the difference between the two forms of the principle of Hamilton is illustrated in Sect. 8 by applying a Ritz-type of method on the action integral for an example system. Conclusions are drawn in Sect. 9 and the paper is closed with a discussion of the usefulness of the results and insight gained in this paper.

## 2 Cones

The analysis of inequalities requires appropriate geometric objects (Aubin and Frankowska 1990). The usual geometric objects in classical dynamics with equality constraints are manifolds and tangent spaces to manifolds. The study of inequalities,
<figure><img src="image_2.png" /><figcaption>**Fig. 1** Contingent, tangent, and normal cone on a non-convex set $K \subset \mathbb{R}^2$</figcaption></figure>

however, requires a number of cones. We first introduce these cones in the Euclidean $\mathbb{R}^n$. Let $K$ be a closed non-empty subset of the $\mathbb{R}^n$. We define the distance from $\mathbf{y}$ to $K$ as

$$d_K(\mathbf{y}) = \inf_{\mathbf{y}^* \in K} \|\mathbf{y} - \mathbf{y}^*\|, \quad (1)$$

with $\|\mathbf{y}\| = \sqrt{\sum_i y_i^2}$. We assume that $K$ is derivable (Aubin and Frankowska 1990; Rockafellar and Wets 1998) in the sense that $\lim_{t \downarrow 0} d_K(\mathbf{y} + t\xi)/t$ exists for all $\mathbf{y} \in K$ and all $\xi \in \mathbb{R}^n$. The set

$$C_K(\mathbf{y}) = \left\{ \xi \mid \lim_{t \downarrow 0} \frac{d_K(\mathbf{y} + t\xi)}{t} = 0 \right\}, \quad (2)$$

is the adjacent cone of $K$ at $\mathbf{y}$ and has the following characterisation in terms of sequences

$$\xi \in C_K(\mathbf{y}) \iff \forall t_k \downarrow 0, \exists \xi_k \to \xi \text{ such that } \mathbf{y} + t_k \xi_k \in K \forall k. \quad (3)$$

If the set $K$ is derivable, then the adjacent cone agrees with the so-called contingent cone (Aubin and Frankowska 1990). In the following, we will always assume that the set $K$ is derivable and simply speak of the contingent cone $C_K(\mathbf{y})$. Loosely speaking, the contingent cone is the cone which is obtained if one zooms in on the point $\mathbf{y}$ after a shift of the set $K$ such that $\mathbf{y}$ is at the origin (see Fig. 1). The contingent cone is closed but not necessarily convex. The tangent cone $\mathcal{T}_K(\mathbf{y})$ defined by

$$\mathcal{T}_K(\mathbf{y}) = \left\{ \xi \mid \lim_{t \downarrow 0, y^* \to_K \mathbf{y}} \frac{d_K(\mathbf{y}^* + t\xi)}{t} = 0 \right\}, \quad (4)$$

where $\to_K$ denotes the convergence in $K$, is a closed convex sub-cone of the contingent cone. We define the normal cone $\mathcal{N}_K(\mathbf{y})$ to be the polar cone of the tangent cone $\mathcal{T}_K(\mathbf{y})$

$$\mathcal{N}_K(\mathbf{y}) = \{z \mid \xi^\mathrm{T}z \le 0 \forall \xi \in \mathcal{T}_K(\mathbf{y})\}. \quad (5)$$
If the tangent cone $\mathcal{T}_K(y)$ is a sub-space of $\mathbb{R}^n$, then the normal cone $\mathcal{N}_K(y)$ is the annihilator and, therefore, a sub-space of its dual $(\mathbb{R}^n)^*$. The normal cone $\mathcal{N}_K(y)$ and tangent cone $\mathcal{T}_K(y)$ are closed convex cones. A set $K$ is called tangentially regular, if the contingent cone $\mathcal{C}_K(y)$ agrees with the tangent cone $\mathcal{T}_K(y)$. Tangential regularity excludes reentrant corners. If $K$ is a convex set, then it is also tangentially regular.

If $\mathcal{M}$ is a differentiable manifold embedded in the $\mathbb{R}^n$, then $\mathcal{T}_{\mathcal{M}}(y)$, which equals $\mathcal{C}_{\mathcal{M}}(y)$, is the tangent space to $\mathcal{M}$ at $y$. Let $\mathcal{K}$ be a closed non-empty subset of the differentiable manifold $\mathcal{M}$. The tangent cone $\mathcal{T}_{\mathcal{K}}(y)$ is a subset of the tangent space $\mathcal{T}_{\mathcal{M}}(y)$, whereas $\mathcal{N}_{\mathcal{K}}(y)$ is a subset of the cotangent space $(\mathcal{T}_{\mathcal{M}}(y))^*$.

The previously defined cones can be generalised to infinite dimensional vector spaces (see, for instance, Aubin and Frankowska 1990). Let $(X, \|\cdot\|)$ be a Banach space and $X^*$ its dual space with the duality pairing $\langle \cdot, \cdot \rangle$. Let $\mathcal{K}$ be a closed, derivable subset of the $X$ and $d_{\mathcal{K}}(y) = \inf_{y^* \in \mathcal{K}} \|y - y^*\|$ be the distance form $y$ to $\mathcal{K}$. We define the contingent cone, tangent cone, and normal cone as (see Aubin and Frankowska 1990)

$$ C_{\mathcal{K}}(y) = \left\{ \xi \in X \left| \lim_{t \downarrow 0} \frac{d_{\mathcal{K}}(y + t\xi)}{t} = 0 \right. \right\}, \quad (6) $$

$$ \mathcal{T}_{\mathcal{K}}(y) = \left\{ \xi \in X \left| \lim_{t \downarrow 0, y^* \to \mathcal{K} y} \frac{d_{\mathcal{K}}(y^* + t\xi)}{t} = 0 \right. \right\}, \quad (7) $$

$$ \mathcal{N}_{\mathcal{K}}(y) = \{z \in X^* \mid \langle \xi, z \rangle \le 0 \ \forall \xi \in \mathcal{T}_{\mathcal{K}}(y)\}. \quad (8) $$

The cones on infinite dimensional vector spaces can be related to the cones on finite dimensional space. Let the space $X = L^\infty(I, \mathbb{R}^n)$ of real-valued $n$-dimensional essentially bounded functions be defined on the compact interval $I \subset \mathbb{R}$ with the duality pairing

$$ \begin{aligned} \langle \mathbf{g}, \mathbf{f} \rangle &= \int_I (\mathbf{g}(x))^{\mathrm{T}} \mathbf{f}(x) \, dx \\ &= \sum_{i=1}^{n} \int_I g_i(x) f_i(x) \, dx, \quad \mathbf{f} \in L^{\infty}(I, \mathbb{R}^n), \quad \mathbf{g} \in L^1(I, \mathbb{R}^n), \end{aligned} \quad (9) $$

and norm $\|\mathbf{f}\|_{\infty} = \sup_{x \in I} \|f(x)\|$. Let $\mathcal{K}$ be a closed subset of $X = L^{\infty}(I, \mathbb{R}^n)$ of the form

$$ \mathcal{K} = \{ \mathbf{y} \in X \mid \mathbf{y}(x) \in K \ \forall x \in I \}, \quad (10) $$

where $K$ is a closed subset of $\mathbb{R}^n$. It is straightforward to prove the equivalence between the contingent and tangent cones in $L^\infty(I, \mathbb{R}^n)$ and $\mathbb{R}^n$, i.e.

$$ \xi \in \mathcal{C}_{\mathcal{K}}(y) \iff \xi(x) \in \mathcal{C}_K(y(x)) \quad \text{for almost all } x \in I, \quad (11) $$

$$ \xi \in \mathcal{T}_{\mathcal{K}}(y) \iff \xi(x) \in \mathcal{T}_K(y(x)) \quad \text{for almost all } x \in I, \quad (12) $$

and, using the dual pairing (9), we also have

$$ z \in \mathcal{N}_{\mathcal{K}}(y) \iff z(x) \in \mathcal{N}_K(y(x)) \quad \text{for almost all } x \in I. \quad (13) $$
### 3 Variational Inequalities

The aim of this section is to define weak and strong local extrema as variational inequalities in infinite dimensional space. First, a motivation for variational inequalities is given by showing that the sub-stationarity condition of a constrained variational minimisation problem leads to a variational inequality. Furthermore, it is explained that the sub-stationarity condition in terms of a Gâteaux derivative only considers the class of comparison functions which are in a neighbourhood to the extremal with respect to a weak norm. This leads to the concept of extremals with respect to a norm.

Let $I \subset \mathbb{R}$ be a real compact interval. Consider functions $y: I \to \mathbb{R}^n$ with the $n$ components $y_i(x) \in \mathbb{R}$, $i = 1, \dots, n$. Let $\mathcal{Y}(I, \mathbb{R}^n)$ be the class of functions $y: I \to \mathbb{R}^n$ which are absolutely continuous, and of which the derivative $y'$ is of specially bounded variation (no singular part) and, therefore, exists almost everywhere in $I$,

$$ \mathcal{Y}(I, \mathbb{R}^n) = \{y \in \mathrm{AC}(I, \mathbb{R}^n), y' \in \mathrm{SBV}(I, \mathbb{R}^n)\}, \quad (14) $$

and let $\mathcal{Y}_0(I, \mathbb{R}^n)$ be the class of functions in $\mathcal{Y}(I, \mathbb{R}^n)$ which have compact support on $I$. Functions in $\mathcal{Y}$ may be piecewise $C^1$, such that the derivative is discontinuous on a Lebesgue negligible set. Functions in $\mathcal{Y}$ may also have accumulation points for which the derivative is of bounded variation. Note that if $y \in \mathrm{AC}(I, \mathbb{R}^n)$, where $I$ is compact, then $y$ is bounded on $I$ and it therefore holds that $\mathcal{Y}(I, \mathbb{R}^n) \subset L^\infty(I, \mathbb{R}^n)$. Furthermore, let the class $\mathcal{D}(I, \mathbb{R}^n) \subset \mathcal{Y}(I, \mathbb{R}^n)$ fulfill the boundary conditions $y(x) = y_0(x)$ on the boundary of $I$, i.e.

$$ \mathcal{D}(I, \mathbb{R}^n) = \{y \in \mathcal{Y}(I, \mathbb{R}^n), y(x) = y_0(x) \ \forall x \in \mathrm{bdry} I\}. \quad (15) $$

Let $\hat{y}(\varepsilon)$ be a family of functions parameterised by $\varepsilon$ for which it holds that $\hat{y}(\varepsilon_0) = y$, i.e. $\hat{y}(\varepsilon_0, x) = y(x) \ \forall x \in I$. If $\hat{y}(\varepsilon, x)$ is locally differentiable with respect to $\varepsilon$, then the difference

$$ \hat{y}(\varepsilon, x) - \hat{y}(\varepsilon_0, x) = \hat{y}_\varepsilon(\varepsilon_0, x)(\varepsilon - \varepsilon_0) + \mathcal{O}((\varepsilon - \varepsilon_0)^2), \quad (16) $$

with $\hat{y}_\varepsilon(\varepsilon, x) = \partial\hat{y}(\varepsilon, x)/\partial\varepsilon$, gives rise to the variation $\delta y(x)$ of $y(x)$

$$ \delta y(x) = \hat{y}_\varepsilon(\varepsilon_0, x) \delta\varepsilon \quad (17) $$

with $\delta\varepsilon = \varepsilon - \varepsilon_0$. The choice of $\varepsilon_0$ is immaterial and we will take $\varepsilon_0 = 0$. Moreover, it suffices to consider $\varepsilon \ge 0$.

Let $\mathcal{K}$ be a closed subset of $\mathcal{D}(I, \mathbb{R})$, which is not necessarily convex. We seek a function $y(x)$ which is the solution of the constrained minimisation problem

$$ \min_{y \in \mathcal{K}} J(y), \quad (18) $$

where $J: \mathcal{D}(I, \mathbb{R}^n) \to \mathbb{R}$ is the functional

$$ J(y) = \int_I f(x, y, y') \, dx. \quad (19) $$
We assume that the functional $J$ is Gâteaux differentiable. We consider all functions $\mathbf{v} \in \mathcal{Y}_0(I, \mathbb{R}^n)$ for which for each sequence $t_k \downarrow 0$ there exists a sequence $\mathbf{v}_k \to \mathbf{v}$ with $\mathbf{y} + t_k \mathbf{v}_k \in \mathcal{K}$. In other words, using the characterisation (3) of the contingent cone (6), it holds that $\mathbf{v} \in C_{\mathcal{K}}(\mathbf{y})$. Because $\mathbf{y}$ minimises $J$ in $\mathcal{K}$, it holds for this sequence that $J(\mathbf{y}) \le J(\mathbf{y} + t_k \mathbf{v}_k)$. Hence, the function $\tilde{h}(t_k) = J(\mathbf{y} + t_k \mathbf{v}_k)$ attains its minimum when $t_k = 0$. Consequently, it holds that (see Kinderlehrer and Stampacchia 1980)

$$ \tilde{h}'(0) \ge 0 \qquad (20) $$

with

$$ \tilde{h}'(0) = \lim_{t_k \downarrow 0, \mathbf{v}_k \to \mathbf{v}} \frac{J(\mathbf{y} + t_k \mathbf{v}_k) - J(\mathbf{y})}{t_k} = \lim_{\varepsilon \downarrow 0} \frac{J(\mathbf{y} + \varepsilon \mathbf{v}) - J(\mathbf{y})}{\varepsilon}, \qquad (21) $$

in which we recognise the (one-sided) Gâteaux derivative

$$ dJ(\mathbf{y}; \mathbf{v}) := \lim_{\varepsilon \downarrow 0} \frac{J(\mathbf{y} + \varepsilon \mathbf{v}) - J(\mathbf{y})}{\varepsilon}. \qquad (22) $$

In this paper, we will not distinguish between the one-sided Gâteaux derivative and the Gâteaux derivative because $\varepsilon$ is defined to be non-negative. Using the Gâteaux derivative (22), we can reformulate the sub-stationarity condition as

$$ dJ(\mathbf{y}; \mathbf{v}) \ge 0 \quad \forall \mathbf{v} \in C_{\mathcal{K}}(\mathbf{y}) \cap \mathcal{Y}_0(I, \mathbb{R}^n), \qquad (23) $$

or as

$$ \int_I (f_y(x, \mathbf{y}, y') \mathbf{v} + f_{y'}(x, \mathbf{y}, y') \mathbf{v}') dx \ge 0 \quad \forall \mathbf{v} \in C_{\mathcal{K}}(\mathbf{y}) \cap \mathcal{Y}_0(I, \mathbb{R}^n), \qquad (24) $$

which is a variational inequality.

The sub-stationarity condition (23) is a necessary condition for a solution of the constrained minimisation problem (18). The sub-stationarity condition (23) is (in general) not a sufficient condition for a minimum of (18) because of two reasons:

1. The sub-stationarity condition (23) only considers the first variation which is in general not sufficient to determine a minimum.

2. The Gâteaux derivative in (23) only considers a certain class of comparison functions $\hat{y}(\varepsilon, x)$, namely those of the form $\hat{y}(\varepsilon, x) = \mathbf{y}(x) + \varepsilon \mathbf{v}(x)$, and nothing is said about the stationarity of the functional with respect to other comparison functions which do not fall within this class.

The second reason precludes the concept of a local minimum with respect to a norm (see Troutman 1996).

**Definition 1** (Local minimum with respect to a norm) In a normed vector space $(\mathcal{Y}, \| \cdot \|_{\alpha})$, a point $\mathbf{y} \in \mathcal{K} \subseteq \mathcal{Y}$ is said to provide $J$ with a local minimum on $\mathcal{K}$ in the norm $\| \cdot \|_{\alpha}$ if $\mathbf{y}$ is a minimiser for $J$ on $\mathcal{B}_{r}^{\alpha}(\mathbf{y})$ for some $r > 0$, i.e.

$$ \exists r > 0 : J(\mathbf{y}^*) \ge J(\mathbf{y}) \quad \forall \mathbf{y}^* \in \mathcal{B}_{r}^{\alpha}(\mathbf{y}), $$
where $\mathcal{B}_r^\alpha(\mathbf{y}) = \{\mathbf{y}^* \in \mathcal{K} \mid \|\mathbf{y}^* - \mathbf{y}\|_\alpha < r\}$ is an open neighbourhood of $\mathbf{y}$ with
respect to the norm $\|\cdot\|_\alpha$.

We introduce on $\mathcal{Y}(I, \mathbb{R})$ two norms:

$$
\begin{equation}
\begin{aligned}
\text{weak norm} \quad \|\mathbf{y}\|_1 &= \sum_{i=1}^{n} \max_{x \in I} |y_i(x)| + \operatorname*{ess\,sup}_{x \in I} |y_i'(x)|, \\
\text{strong norm} \quad \|\mathbf{y}\|_0 &= \sum_{i=1}^{n} \max_{x \in I} |y_i(x)|.
\end{aligned}
\tag{25}
\end{equation}
$$

The weak norm uses the essential supremum because $y_i'(x)$ is not defined for all $x \in I$.

We speak of a weak/strong local minimum if the local minimum is with respect to the weak/strong norm. The designations of the weak and the strong norm are sometimes reversed in literature, but there is a general agreement about the terminology for weak/strong local minima. It holds that $\|\mathbf{y}^* - \mathbf{y}\|_0 \le \|\mathbf{y}^* - \mathbf{y}\|_1$ and a point which is in the weak neighbourhood $\mathcal{B}_r^1(\mathbf{y})$ is therefore also in the strong neighbourhood $\mathcal{B}_r^0(\mathbf{y})$, i.e. $\mathcal{B}_r^1(\mathbf{y}) \subseteq \mathcal{B}_r^0(\mathbf{y})$. A strong local minimum minimises $J$ with respect to all neighbouring points in the strong neighbourhood $\mathcal{B}_r^0(\mathbf{y})$. A weak local minimum minimises $J$ with respect to all neighbouring points in the weak neighbourhood $\mathcal{B}_r^1(\mathbf{y})$. A strong local minimum is necessarily also a weak local minimum, because $\mathcal{B}_r^1(\mathbf{y}) \subseteq \mathcal{B}_r^0(\mathbf{y})$, but the converse is not necessarily true.

We will now give two theorems which give necessary conditions for a weak and a
strong local minimum.

**Theorem 1 (Weak local minimum)** A weak local minimum of $\min_{y \in \mathcal{K}} J(y)$ satisfies

$$
dJ(\mathbf{y}; \mathbf{v}) \ge 0 \quad \forall \mathbf{v} \in \mathcal{C}_{\mathcal{K}}(\mathbf{y}) \cap \mathcal{Y}_0(I, \mathbb{R}^n). \tag{26}
$$

*Proof* We will prove the theorem with a reductio ad absurdum. If the theorem does not hold, then **y** is a weak local minimiser of $\min_{y \in \mathcal{K}} J(y)$ for which there exists a $\mathbf{v} \in \mathcal{C}_{\mathcal{K}}(\mathbf{y})$ such that $dJ(\mathbf{y}; \mathbf{v}) < 0$. The Gâteaux derivative considers the family of comparison functions $\hat{\mathbf{y}}(\varepsilon, x) = \mathbf{y}(x) + \varepsilon\mathbf{v}(x)$. This family of comparison functions converges in the weak norm $\|\cdot\|_1$, i.e.

$$
\lim_{\varepsilon \downarrow 0} \| \hat{\mathbf{y}}(\varepsilon) - \mathbf{y} \|_1 = \sum_{i=1}^{n} \lim_{\varepsilon \downarrow 0} \max_{x \in I} |\varepsilon v_i(x)| + \lim_{\varepsilon \downarrow 0} \operatorname*{ess\,sup}_{x \in I} |\varepsilon v_i'(x)| = 0.
$$

Hence, there exists a function $\mathbf{y}^*$, in an arbitrarily small weak neighbourhood $\mathcal{B}_r^1(\mathbf{y})$
of $\mathbf{y}$, such that $J(\mathbf{y}^*) < J(\mathbf{y})$. This contradicts the assumption that $\mathbf{y}$ is a weak local
minimiser.
$\square$

We will elaborate a bit more on weak local minima and their relationships with
the Gâteaux derivative. For a weak local minimum, it holds that there exists a $r > 0$
<figure><img src="image_3.png" /><figcaption>**Fig. 2** The functions $y_i(x)$, $v_i(x)$ and $\hat{y}_i(\varepsilon, x) = y_i(x) + \varepsilon v_i(x)$</figcaption></figure>

such that $J(y^*) \ge J(y)$ for all $y^* \in \mathcal{B}_r^1(y)$; see Definition 1. We can address each $y^*$ in $\mathcal{B}_r^1(y)$ by considering families of comparison functions $\hat{y}(\varepsilon)$ which satisfy

$$ \lim_{\varepsilon \downarrow 0} \| \hat{y}(\varepsilon) - y \|_1 = 0, \qquad (27) $$

with $\hat{y}(0) = y$ and $\hat{y}(\varepsilon^*) = y^*$ for some $\varepsilon^* \ge 0$.

The weak norm $\| \cdot \|_1$ measures the point-wise difference of $\hat{y}_i(\varepsilon, x) - y_i(x)$ as well as the point-wise difference of $\hat{y}_{i,x}(\varepsilon, x) - y'_i(x)$ (whenever they exist). Consequently, the vanishing of the limit (27) implies that

$$ \lim_{\varepsilon \downarrow 0} \sup_{x \in I} |\hat{y}_{i,x}(\varepsilon, x) - y'_i(x)| = 0. \qquad (28) $$

In other words, the maximal jump height of $\hat{y}_{i,x}(\varepsilon, x) - y'_i(x)$ on $I$ has to vanish when $\varepsilon$ is approaching zero. Proposition 1 (see the Appendix) implies that all discontinuity points of $y'_i$ are discontinuity points of $\hat{y}_{i,x}(\varepsilon)$ for sufficiently small $\varepsilon$. Moreover, the jump height of those discontinuity points of $\hat{y}_{i,x}(\varepsilon, x)$ which are not discontinuity points of $y'_i(x)$ must vanish for $\varepsilon \downarrow 0$.

The family of comparison functions $\hat{y}_i(\varepsilon)$ is therefore of the form $\hat{y}_i(\varepsilon) = y_i + \varepsilon v_i + \mathcal{O}(\varepsilon^2)$, where $v_i \in \mathcal{Y}_0(I, \mathbb{R})$. This is illustrated in Fig. 2. The kink in the function $\hat{y}_i(\varepsilon, x)$ which is not a kink of $y_i(x)$ is due to a kink in $v_i(x)$ and, therefore, vanishes for $\varepsilon \downarrow 0$. If $y$ is a weak local minimiser of $J$ on $\mathcal{K}$, then it holds that $J(\hat{y}(\varepsilon)) \ge J(y)$ for all $\hat{y}(\varepsilon) \in \mathcal{B}_r^1(y)$ with sufficiently small $\varepsilon > 0$. The weak norm implies that we are only looking at comparison functions of the form $\hat{y}(\varepsilon) = y + \varepsilon v + \mathcal{O}(\varepsilon^2) \in \mathcal{K}$. If $y + \varepsilon v + \mathcal{O}(\varepsilon^2) \in \mathcal{K}$ and $\varepsilon$ is approaching 0, then $v$ must lie in the contingent cone $C_{\mathcal{K}}(y)$. Consequently, for a weak local minimum, it holds that

$$ \lim_{\varepsilon \downarrow 0} \frac{J(y + \varepsilon v) - J(y)}{\varepsilon} \ge 0 \quad \forall v \in C_{\mathcal{K}}(y) \cap \mathcal{Y}_0(I, \mathbb{R}^n), $$

in which we recognise the Gâteaux derivative $dJ(y; v)$.

Consider, for instance, the family of comparison functions $\hat{y}(\varepsilon)$ which are of the form $\hat{y}(\varepsilon, x) = y(x - \hat{x}(\varepsilon, x)) + \varepsilon v(x - \hat{x}(\varepsilon, x))$ involving a value shift $\varepsilon v$ as well
<figure><img src="image_4.png" /><figcaption>**Fig. 3** The functions $y_i(x)$, $\hat{x}(\epsilon, x)$ and $\hat{y}_i(\epsilon, x) = y_i(x - \hat{x}(\epsilon, x))$</figcaption></figure>

as a shift $\hat{x}(\epsilon)$ of the argument with $\lim_{\epsilon \downarrow 0} \hat{x}(\epsilon, x) = 0$. If this family of comparison functions converges to $\mathbf{y}$ in the weak norm $\| \cdot \|_1$, then it has to fulfill

$$ \lim_{\epsilon \downarrow 0} \sup_{x \in I} \left| y'_{i}(x - \hat{x}(\epsilon, x)) + \epsilon v'_{i}(x - \hat{x}(\epsilon, x)) - y'_{i}(x) \right| \\ = \lim_{\epsilon \downarrow 0} \sup_{x \in I} \left| y'_{i}(x - \hat{x}(\epsilon, x)) - y'_{i}(x) \right| = 0, $$

which implies that $\hat{x}(\epsilon, x_s) = 0$ if $y'_i$ is discontinuous at $x_s$. We see from this example that a family of comparison functions which involves a shift of the argument at discontinuity points can not converge in the weak norm. This is illustrated in Fig. 3. The family of comparison function $\hat{y}_i(\epsilon, x) = y_i(x - \hat{x}(\epsilon, x))$ converges in the strong norm to $y_i(x)$ but not in the weak norm, because the kinks do not appear for the same values of $x$.

Subsequently, we would like to consider a strong local minimum of the constrained minimisation problem (18). We consider therefore a family of functions $\hat{y}(\epsilon) \in \mathcal{K}$, with $\hat{y}(0) = \mathbf{y}$, which converge in the strong norm $\| \cdot \|_0$. Using

$$ h(\epsilon) = J(\hat{y}(\epsilon)) = \int_I f(x, \hat{y}(\epsilon, x), \hat{y}_x(\epsilon, x)) \, dx, \qquad (29) $$

we define the variation $\delta J$ of the functional $J$ as

$$ \delta J = h'(0)\delta\epsilon = \delta\epsilon \lim_{\epsilon \downarrow 0} \frac{J(\hat{y}(\epsilon)) - J(\mathbf{y})}{\epsilon} \qquad (30) $$

for comparison functions $\hat{y}(\epsilon)$ such that $\|\hat{y}(\epsilon) - \mathbf{y}\|_0 \to 0$ for $\epsilon \downarrow 0$.

**Theorem 2** (Strong local minimum) A strong local minimum of $\min_{y \in \mathcal{K}} J(y)$ satisfies

$$ \delta J \ge 0 \quad \forall \delta \mathbf{y} \in C_{\mathcal{K}}(\mathbf{y}). \qquad (31) $$
*Proof* If $\mathbf{y}$ is a strong local minimiser of $J$ on $\mathcal{K}$, then the function $h(\varepsilon)$ has a strong local minimum for $\varepsilon = 0$ which gives the sub-stationarity condition

$$0 \le h'(0) = \lim_{\varepsilon \downarrow 0} \frac{J(\hat{\mathbf{y}}(\varepsilon)) - J(\mathbf{y})}{\varepsilon} \quad (32)$$

in which we recognise the variation $\delta J$. The comparison function $\hat{\mathbf{y}}(\varepsilon)$ is in $\mathcal{K}$ for infinitely small values of $\varepsilon$ if $\delta \mathbf{y} = \hat{\mathbf{y}}_{\varepsilon}|_{\varepsilon=0} \delta \varepsilon$ is in the contingent cone $C_{\mathcal{K}}(\mathbf{y})$. $\square$

The properties (26) and (31) of weak and strong local minima can be viewed as sub-stationarity conditions. These properties can be used as *definition* for weak and strong local extrema of the functional $J$ on $\mathcal{K}$.

**Definition 2** (Weak local extremal) A weak local extremal $\mathbf{y}$ of $J$ on $\mathcal{K}$ is defined as a function fulfilling

$$dJ(\mathbf{y}; \mathbf{v}) \ge 0 \quad \forall \mathbf{v} \in C_{\mathcal{K}}(\mathbf{y}) \cap \mathcal{Y}_0(I, \mathbb{R}^n).$$

**Definition 3** (Strong local extremal) A strong local extremal $\mathbf{y}$ of $J$ on $\mathcal{K}$ is defined as a function fulfilling

$$\delta J \ge 0 \quad \forall \delta \mathbf{y} \in C_{\mathcal{K}}(\mathbf{y}).$$

Note that $\delta \mathbf{y} \in C_{\mathcal{K}}(\mathbf{y})$ fulfills the boundary condition $\delta \mathbf{y}(x) = \mathbf{0}$ for $x \in \text{bdry } I$, because $\mathcal{K} \subset \mathcal{D}(I, \mathbb{R})$. Weak and strong local extrema will play an important role in the following sections as the sub-stationarity conditions (26) and (31) constitute variational inequalities.

If $\mathcal{K}$ is of the form (10)

$$\mathcal{K} = \{ \mathbf{y} \in \mathcal{D}(I, \mathbb{R}^n) \mid \mathbf{y}(x) \in K \ \forall x \in I \}, \quad (33)$$

where $K$ is a closed subset of $\mathbb{R}^n$, then we can write a variational inequality such as (23) by using the notation

$$dJ(\mathbf{y}; \mathbf{v}) \ge 0 \quad \forall \mathbf{v}(x) \in C_K(\mathbf{y}(x)) \text{ almost everywhere on } I, \mathbf{v} \in \mathcal{Y}_0(I, \mathbb{R}^n), \quad (34)$$

or in the more sloppy form

$$dJ(\mathbf{y}; \mathbf{v}) \ge 0 \quad \forall \mathbf{v} \in C_K(\mathbf{y}) \text{ a.e. on } I, \mathbf{v} \in \mathcal{Y}_0(I, \mathbb{R}^n), \quad (35)$$

in which the argument $x$ has been suppressed. We will use the latter notation in the following sections.

## 4 The Principle of d'Alembert-Lagrange

In this section, we introduce the concept of a perfect constraint through the principle of d'Alembert-Lagrange. We first discuss the principle of d'Alembert-Lagrange for a
bilateral constraint. Subsequently, we discuss a unilateral constraint and define a perfect unilateral constraint through the principle of d'Alembert-Lagrange in inequality form (Goeleven et al. 1997; Panagiotopoulos and Glocker 1998).

Consider a point-mass $m$ with position $\mathbf{r}(t) \in \mathbb{R}^3$ which is confined by a bilateral geometric scleronomic constraint $g(\mathbf{r})$ with $g \in C^1(\mathbb{R}^3, \mathbb{R})$. In other words, the point-mass is constrained to the two-dimensional manifold $\mathcal{M} = \{\mathbf{r} \in \mathbb{R}^3 | g(\mathbf{r}) = 0\}$ embedded in the $\mathbb{R}^3$; see Fig. 4a. The constraint is enforced by the constraint force $\mathbf{R}$. Furthermore, the mass is subjected to an external force $\mathbf{F}$. Newton's second law for the point-mass therefore reads as

$$m\ddot{\mathbf{r}} = \mathbf{R} + \mathbf{F}. \qquad (36)$$

The principle of virtual work states that if the virtual work (Glocker 2001) of the system vanishes for all virtual displacements $\delta\mathbf{r}$, i.e.

$$\delta W = (m\ddot{\mathbf{r}} - \mathbf{F} - \mathbf{R})^T \delta\mathbf{r} = 0 \quad \forall \delta\mathbf{r}, \qquad (37)$$

then the system is in dynamic equilibrium (the Newton-Euler equations are fulfilled). The virtual displacements $\delta\mathbf{r}$ in (37) are arbitrary and can therefore be admissible as well as inadmissible with respect to the constraint $\mathbf{r} \in \mathcal{M}$, see Glocker (2001). Virtual displacements $\delta\mathbf{r}$ are admissible with respect to the constraint if they belong to the tangent space $\mathcal{T}_{\mathcal{M}}(\mathbf{r}) = \{z \in \mathbb{R}^3 | \frac{\partial g}{\partial \mathbf{r}} z = 0\}$. If the virtual work of the constraint force vanishes for all virtual displacements which are admissible with respect to the constraint, i.e.

$$\mathbf{R}^T \delta\mathbf{r} = 0 \quad \forall \delta\mathbf{r} \in \mathcal{T}_{\mathcal{M}}(\mathbf{r}), \qquad (38)$$

then we speak of a perfect bilateral constraint. A perfect constraint force of a bilateral constraint is therefore 'normal' to the constraint manifold in the sense that

$$-\mathbf{R} \in \mathcal{N}_{\mathcal{M}}(\mathbf{r}), \qquad (39)$$

where $\mathcal{N}_{\mathcal{M}}(\mathbf{r})$ is the set of all vectors which annihilate the elements of the tangent space $\mathcal{T}_{\mathcal{M}}(\mathbf{r})$. Instead of taking (38) as definition of a perfect constraint, we could also have chosen to take the normality condition (39) of the constraint force as definition. The supposition of perfect bilateral constraints in classical mechanics and, therefore, the normality of the constraint force to the constraint manifold, excludes phenomena such as friction. Combining (38) and (37) gives the variational equality

$$(m\ddot{\mathbf{r}} - \mathbf{F})^T \delta\mathbf{r} = 0 \quad \forall \delta\mathbf{r} \in \mathcal{T}_{\mathcal{M}}(\mathbf{r}), \qquad (40)$$

which is usually referred to as the classical principle of d'Alembert-Lagrange.

We now impose in addition unilateral geometric scleronomic constraints $g_i(\mathbf{r}) \ge 0$, $i = 1...m$. The position of the point-mass is now constrained to the subset $\mathcal{K} = \{\mathbf{r} \in \mathcal{M} | g_i(\mathbf{r}) \ge 0, i = 1...m\}$, which is a sub-manifold with boundary; see Fig. 4b. Virtual displacements $\delta\mathbf{r}$ are admissible with respect to the bilateral and unilateral constraints if $\delta\mathbf{r} \in \mathcal{C}_{\mathcal{K}}(\mathbf{r})$, where $\mathcal{C}_{\mathcal{K}}(\mathbf{r})$ is the contingent cone (6) on the subset $\mathcal{K}$. In the following, we will confine ourselves to the case that $\mathcal{K}$ is a tangentially regular set
<figure><img src="image_5.png" /><figcaption>Fig. 4 Constrained point-mass</figcaption></figure>

for which the contingent cone agrees with the tangent cone (7). We define a unilateral constraint as perfect if the constraint force satisfies the normality condition (39)

$$-\mathbf{R} \in \mathcal{N}_{\mathcal{K}}(\mathbf{r}), \qquad (41)$$

where $\mathcal{N}_{\mathcal{K}}(\mathbf{r})$ is the normal cone (8) to the subset $\mathcal{K}$. The tangent cone $T_{\mathcal{K}}(\mathbf{r})$ and normal cone $\mathcal{N}_{\mathcal{K}}(\mathbf{r})$ are polar in the sense that for all $\mathbf{R}$ and $\delta\mathbf{r}$ satisfying $-\mathbf{R} \in \mathcal{N}_{\mathcal{K}}(\mathbf{r})$ and $\delta\mathbf{r} \in T_{\mathcal{K}}(\mathbf{r})$ it holds that $-\mathbf{R}^T\delta\mathbf{r} \le 0$. Hence, for perfect unilateral constraints, it holds that the virtual work of the constraint force is non-negative for kinematically admissible virtual displacements, i.e.

$$\mathbf{R}^T\delta\mathbf{r} \ge 0 \quad \forall \delta\mathbf{r} \in T_{\mathcal{K}}(\mathbf{r}), \qquad (42)$$

which is sometimes referred to as Fourier's inequality (Lanczos 1962). Combining (42) and (37) gives the variational inequality

$$(m\ddot{\mathbf{r}} - \mathbf{F})^T \delta\mathbf{r} \ge 0 \quad \forall \delta\mathbf{r} \in T_{\mathcal{K}}(\mathbf{r}), \qquad (43)$$

which we will refer to as the principle of d'Alembert-Lagrange in inequality form.

## 5 The Principle of Hamilton for Non-impulsive Motion

The principle of Hamilton is derived in this section for a system with perfect bilateral as well as unilateral constraints and is shown to be a variational inequality. In this section, the motion is assumed to be non-impulsive. As a consequence, impulses are absent and velocities are absolutely continuous functions in time. The principle of Hamilton for perfect unilateral constraints with impulsive motion will be derived in Sects. 6 and 7.

Consider a system S consisting of material points, e.g. an elastic continuum body or a multi-body system consisting of rigid bodies. For convenience, we will make
the tacit assumption that each material point has a different position vector $X$ in the reference configuration at a reference time $t_0$. We can therefore address each material point in $S$ by its position vector $X \in \mathbb{R}^3$ in the reference configuration. Each material point $X$ has a position $x$ in the current configuration at time $t$. We can define a vector field $\xi$ which carries material points $X$ to places $x = \xi(X, t)$ in their current configuration. The vector field $\xi$ is called the motion of the system and induces the velocity field $\dot{\xi}$ and acceleration field $\ddot{\xi}$ defined by

$$ \dot{\xi}(X, t) = \frac{\partial \xi(X, t)}{\partial t}, \qquad \ddot{\xi}(X, t) = \frac{\partial^2 \xi(X, t)}{\partial t^2}, \qquad (44) $$

whenever they exist. The system $S$ is subjected to bilateral geometric scleronomic constraints, which constrain the motion $\xi$ to an $n$-dimensional differentiable sub-manifold $\mathcal{M}$ of the configuration space, which is called the configuration manifold. Furthermore, the system is subjected to unilateral geometric scleronomic constraints, which constrain the motion $\xi$ to a closed subset $\mathcal{K}$ of configuration manifold $\mathcal{M}$. The set of admissible motions $\mathcal{K}$ is therefore a sub-manifold with boundary. Let $dm$ be the mass element at $\xi(X, t)$. The mass-element $dm$ is subjected to external and internal forces $dF$, which are not due to the constraints (e.g. elastic forces and gravitational forces). We will assume that these forces $dF$ are potential forces, which will allow us to formulate a principle of Hamilton. Furthermore, the mass-element $dm$ is subjected to constraint forces $d\mathbf{R}$ which impose the bilateral and unilateral geometric scleronomic constraints.

We define the virtual work of the system $S$ as

$$ \delta W = \int_S \delta \xi^T (\ddot{\xi} \, dm - dF - d\mathbf{R}), \qquad (45) $$

which exists almost everywhere on the time-interval $I$. The principle of virtual work states that if the virtual work vanishes for all virtual variations $\delta\xi$, i.e.

$$ \delta W = 0 \quad \forall \delta\xi \text{ a.e. on } I, \qquad (46) $$

then the system is in dynamic equilibrium. The constraint forces $d\mathbf{R}$ are assumed to be perfect and it therefore holds that

$$ d\mathbf{R}^T \delta\xi \geq 0 \quad \forall \delta\xi \in \mathcal{T}_{\mathcal{K}}(\xi) \text{ a.e. on } I. \qquad (47) $$

A virtual displacement $\delta\xi$ which is in the tangent cone $\mathcal{T}_{\mathcal{K}}(\xi)$ is kinematically admissible with respect to the constraint $\xi \in \mathcal{K}$. The principle of d'Alembert-Lagrange in inequality form (43) for the system $S$ therefore reads as

$$ \int_S \delta \xi^T (\ddot{\xi} \, dm - dF) \geq 0 \quad \forall \delta \xi \in \mathcal{T}_{\mathcal{K}}(\xi) \text{ a.e. on } I. \qquad (48) $$

We now can define a chart with coordinates $q \in \mathbb{R}^n$ by choosing a one-to-one mapping $\varphi: U \subset \mathbb{R}^n \to \varphi(U) \subset \mathcal{M}$ which maps the coordinates $q \in U \subset \mathbb{R}^n$ in some open set $U$ to a subset of the configuration manifold $\mathcal{M}$. The motion $\xi(\cdot, t)$ is
therefore for each $t \in I$ uniquely determined (at least within some open set) by the
generalised coordinates $q$ through the mapping $\varphi$, i.e.

$$ \xi(\cdot, t) = \varphi(q(t)). \qquad (49) $$

The coordinates $q$ form a set of generalised coordinates which are minimal with
respect to the bilateral constraints. The unilateral constraints, which enforce the con-
straint $\xi \in \mathcal{K}$, induce a set $K = \varphi^{-1}(\mathcal{K})$ of admissible generalised coordinates which
we assume to be tangentially regular. The tangent cone $T_K(q)$ therefore agrees with
the contingent cone and is a subset of the chart tangent space. Let $\frac{\partial \varphi}{\partial q}$ denote the map-
ping from the chart tangent space to the tangent space $T_M(\xi)$ of the configuration
manifold. An arbitrary virtual displacement $\delta q$ on the chart tangent space causes a
variation of the motion $\delta \xi = \frac{\partial \varphi}{\partial q} \delta q \in T_M(\xi)$, which is a element of the tangent space
$T_M(\xi)$ and which is therefore kinematically admissible with respect to the bilateral
constraints. A virtual displacement $\delta q \in T_K(q)$ induces a variation of the motion
$\delta \xi = \frac{\partial \varphi}{\partial q} \delta q \in T_K(\xi)$ which is in addition admissible with respect to the unilateral
constraints.

The virtual work of the inertia forces in (48) can be rewritten as

$$ \int_S \delta \xi^T \ddot{\xi} \, dm = \frac{d}{dt} \left( \int_S \delta \xi^T \dot{\xi} \, dm \right) - \int_S \delta \dot{\xi}^T \dot{\xi} \, dm, \qquad (50) $$

or, by using the generalised coordinates $q$, as

$$ \begin{aligned} \int_S \delta \xi^T \ddot{\xi} \, dm &= \frac{d}{dt} \left( \delta q^T \int_S \left( \frac{\partial \varphi}{\partial q} \right)^T \dot{\xi} \, dm \right) - \int_S \delta \dot{\xi}^T \dot{\xi} \, dm \\ &= \frac{d}{dt} (\delta q^T p) - \delta T, \end{aligned} \qquad (51) $$

where

$$ p = \int_S \left( \frac{\partial \varphi}{\partial q} \right)^T \dot{\xi} \, dm \qquad (52) $$

is the generalised momentum and $\delta T$ is the variation of the kinetic energy

$$ T = \int_S \frac{1}{2} \dot{\xi}^T \dot{\xi} dm \dot{\xi} = \frac{1}{2} \dot{q}^T M(q) \dot{q}, \quad M(q) = \int_S \left( \frac{\partial \varphi}{\partial q} \right)^T dm \left( \frac{\partial \varphi}{\partial q} \right). \qquad (53) $$

Moreover, we introduce

$$ f = \int_S \left( \frac{\partial \varphi}{\partial q} \right)^T dF \qquad (54) $$

as the generalised force. The forces $dF$ are assumed to be potential forces and the
generalised force $f$ is therefore also a potential force $-f = (\partial V / \partial q)^T$ where $V(q)$ is
the potential energy and $\delta V = -\delta q^T f$. Substitution of (51) and (54) in the principle
of d'Alembert-Lagrange (48) gives

$$ \frac{d}{dt} (\delta q^T p) - \delta T - \delta q^T f \ge 0 \quad \forall \delta q \in T_K(q) \text{ a.e. on } I, \qquad (55) $$
which is the inequality form of the well-known Lagrange central equation (Bremer 1988; Hamel 1912; Papastavridis 2002). The Lagrange central equation (55) holds at each time-instant $t$ for which the generalised velocities $\dot{q}(t)$ and accelerations $\ddot{q}(t)$ exist. If the motion on $I$ is assumed to be non-impulsive, then the generalised positions $q(t)$ and generalised velocities $\dot{q}(t)$ are absolutely continuous on $I$ and the generalised accelerations $\ddot{q}(t)$ exist for almost all $t \in I$. Hence, we can integrate the central equation (55) over a non-impulsive time-interval $I = [t_0, t_f]$ which gives

$$[\delta \boldsymbol{q}^{\mathrm{T}} \boldsymbol{p}]_{t_0}^{t_f} - \int_{I} (\delta T + \delta \boldsymbol{q}^{\mathrm{T}} \boldsymbol{f}) \mathrm{d}t \geq 0 \quad \forall \delta \boldsymbol{q} \in \mathcal{T}_{K}(\boldsymbol{q}). \qquad (56)$$

By defining the Lagrange function $L = T - V$ we arrive at

$$[\delta \boldsymbol{q}^{\mathrm{T}} \boldsymbol{p}]_{t_0}^{t_f} - \int_{I} \delta L \mathrm{d}t \geq 0 \quad \forall \delta \boldsymbol{q} \in \mathcal{T}_{K}(\boldsymbol{q}). \qquad (57)$$

If the boundary conditions are fixed, then the variation $\delta q(t)$ vanishes at $t = t_0$ and $t = t_f$ and we are allowed to interchange the order of integration and variation such that

$$-\delta \int_{I} L \mathrm{d}t \geq 0 \quad \forall \delta \boldsymbol{q} \in \mathcal{T}_{K}(\boldsymbol{q}), \quad \boldsymbol{q}(t_0) = \boldsymbol{q}_0, \quad \boldsymbol{q}(t_f) = \boldsymbol{q}_f, \qquad (58)$$

which is the principle of Hamilton in inequality form for a non-impulsive time-interval $I = [t_0, t_f]$. From the principle of Hamilton, we can derive the Euler-Lagrange equations in inequality form by evaluating the variation in (58) as

$$-\int_{I} \left( \frac{\partial L}{\partial \boldsymbol{q}} - \frac{d}{dt} \frac{\partial L}{\partial \dot{\boldsymbol{q}}} \right) \delta \boldsymbol{q} \mathrm{d}t \geq 0 \quad \forall \delta \boldsymbol{q} \in \mathcal{T}_{K}(\boldsymbol{q}), \quad \boldsymbol{q}(t_0) = \boldsymbol{q}_0, \quad \boldsymbol{q}(t_f) = \boldsymbol{q}_f. \qquad (59)$$

Hence, for almost all $t \in I$, the variational inequality

$$-\left(\frac{\partial L}{\partial \boldsymbol{q}} - \frac{d}{dt} \frac{\partial L}{\partial \dot{\boldsymbol{q}}}\right) \delta \boldsymbol{q} \geq 0 \quad \forall \delta \boldsymbol{q} \in \mathcal{T}_{K}(\boldsymbol{q}) \qquad (60)$$

should hold, which can be cast into the form

$$(M(\boldsymbol{q})\ddot{\boldsymbol{q}} - \boldsymbol{h}(\boldsymbol{q}, \dot{\boldsymbol{q}}))^{\mathrm{T}} \delta \boldsymbol{q} \geq 0 \quad \forall \delta \boldsymbol{q} \in \mathcal{T}_{K}(\boldsymbol{q}) \quad \text{a.e. on } I, \qquad (61)$$

where the mass matrix $M(\boldsymbol{q})$ is defined by (53) and the vector

$$\boldsymbol{h}(\boldsymbol{q}, \dot{\boldsymbol{q}}) = -\left(\frac{d}{dt} \boldsymbol{M}(\boldsymbol{q})\right) \dot{\boldsymbol{q}} + \left(\frac{\partial T}{\partial \boldsymbol{q}} - \frac{\partial V}{\partial \boldsymbol{q}}\right)^{\mathrm{T}} \qquad (62)$$

contains all smooth forces. The variational inequality (61) leads to the differential inclusion

$$M(\boldsymbol{q})\ddot{\boldsymbol{q}} - \boldsymbol{h}(\boldsymbol{q}, \dot{\boldsymbol{q}}) - \boldsymbol{f}_{\mathrm{R}} = \mathbf{0}, \quad -\boldsymbol{f}_{\mathrm{R}} \in \mathcal{N}_{K}(\boldsymbol{q}) \quad \text{a.e. on } I, \qquad (63)$$

where $\boldsymbol{f}_{\mathrm{R}}$ is the generalised constraint force, which fulfills the principle of d'Alembert-Lagrange in inequality form

$$\boldsymbol{f}_{\mathrm{R}}^{\mathrm{T}} \delta \boldsymbol{q} \geq 0 \quad \forall \delta \boldsymbol{q} \in \mathcal{T}_{K}(\boldsymbol{q}). \qquad (64)$$
## 6 The Strong Principle of Hamilton

Up to now, the principle of Hamilton in inequality form (58) has been proven to hold for non-impulsive motion. Impacts and discontinuities in the velocity have therefore not been considered. Now the question rises: Under which conditions does the principle of Hamilton in inequality form (58) also hold for a time-interval with impulsive motion? We will derive necessary conditions for a strong local extremal of the action integral by making use of the general variation of a functional. This will lead to conditions which are similar to the classical Weierstrass–Erdmann corner conditions. The main idea (see Panagiotopoulos and Glocker 1998, 2000) is to split the time-interval $I$ in two non-impulsive sub-intervals $I_1$ and $I_2$ and to set up the variational conditions at the impact time $t_c$ which joins the two sub-intervals $I_1$ and $I_2$. This requires the variation of the action integral over time-intervals with variable begin or end time, i.e. a general variation of the action integral. This section is organised as follows. First, the strong variation of the position function $q(t)$ is introduced. Then the general variation of the action integral is derived. Subsequently, these results are used to treat the variation of the action integral over $I$ by splitting the interval in two non-impulsive sub-intervals which leads to generalised Weierstrass–Erdmann corner conditions. Finally, we arrive at a strong principle of Hamilton in inequality form for impulsive motion.

Let $q(t) \in D(I, \mathbb{R}^n)$ describe the motion of the system (on the chart) on the compact time-interval $I = [t_0, t_f]$. The velocity $u(t)$ is therefore of bounded variation, which means that we can define the left limit $u^-(t)$ and right limit $u^+(t)$ for each time $t$ and $u(t) = \dot{q}(t)$ for almost all $t$. We consider the family of comparison functions $\hat{q}(\varepsilon, t) = q(t - \hat{t}(\varepsilon, t)) + \varepsilon w(t - \hat{t}(\varepsilon, t))$. The function $\hat{t}(\varepsilon, t)$, being continuous and differentiable with the property $\hat{t}(0, t) = 0$, induces a virtual time-shift $\delta t(t) = \hat{t}_\varepsilon(0, t) \delta \varepsilon$. The function $w \in Y_0(I, \mathbb{R}^n)$ induces a value shift. The family of comparison functions $\hat{q}(\varepsilon, t)$ converges to $q(t)$ in the strong norm $\| \cdot \|_0$, because $\lim_{\varepsilon \downarrow 0} \hat{q}(\varepsilon, t) = \hat{q}(0, t) = q(t)$. The variation

$$ \delta q(t) = \left. \frac{\partial \hat{q}(\varepsilon, t)}{\partial \varepsilon} \right|_{\varepsilon=0} \delta \varepsilon = -u^{\pm}(t) \delta t(t) + w(t) \delta \varepsilon \quad (65) $$

is discontinuous and not defined for those time-instants for which $u(t)$ is discontinuous and not defined. Hence, the variation $\delta q(t)$ is of bounded variation and admits for each $t \in I$ a left and a right limit

$$ \delta q^{\pm}(t) = -u^{\pm}(t) \delta t(t) + w(t) \delta \varepsilon. \quad (66) $$

The variation $\delta q(t)$ should be understood as an infinitesimal difference between the comparison function $\hat{q}(\varepsilon, t)$ and the function $q(t)$ for the same value of $t$. Likewise, we can introduce the variation $\delta \overline{q}(t) \approx \hat{q}(\varepsilon, t + \delta t(t)) - q(t)$, or more explicitly using (65)

$$ \overline{\delta q}(t) = q_t(t) \delta t(t) + \hat{q}_\varepsilon(0, t) \delta \varepsilon = u(t) \delta t(t) + \delta q(t) = w(t) \delta \varepsilon, \quad (67) $$

for which also the time is varied. Clearly, the variation $\overline{\delta q}(t)$ is continuous, because $w(t)$ is continuous.
<figure><img src="image_6.png" /><figcaption>**Fig. 5** The general variation of a function $q(t)$</figcaption></figure>

We first consider a time-interval $[t_a, t_b] \subset I$ on which the motion $q(t)$ is differentiable in its interior, i.e. $u(t)$ is continuous on the open interval $(t_a, t_b)$. The general variation of the action integral

$$s(\boldsymbol{q}) = \int_{t_a}^{t_b} L \, dt \qquad (68)$$

involves a variation of the begin point $q(t_a) = \boldsymbol{q}_a$ and end point $q(t_b) = \boldsymbol{q}_b$ as well as a variable begin time $t_a$ and end time $t_b$; see Fig. 5. We introduce the function

$$h(\varepsilon) = \int_{t_a(\varepsilon)}^{t_b(\varepsilon)} L(\hat{q}(\varepsilon, t), \hat{u}(\varepsilon, t)) \, dt, \qquad (69)$$

where $\hat{q}(\varepsilon, t) = q(t - \hat{t}(\varepsilon, t)) + \varepsilon w(t - \hat{t}(\varepsilon, t))$. The general variation $\delta s$ of (68) is defined as the $\varepsilon$-derivative $\delta s = h'(0) \delta \varepsilon$, i.e.

$$\begin{aligned} \delta s &= \int_{t_a}^{t_b} \delta L \, dt + [L \, \delta t]_{t \downarrow t_a}^{t \uparrow t_b} \\ &= \int_{t_a}^{t_b} \left( \frac{\partial L}{\partial \boldsymbol{q}} \delta \boldsymbol{q} + \frac{\partial L}{\partial \boldsymbol{u}} \delta \boldsymbol{u} \right) dt + [L \, \delta t]_{t \downarrow t_a}^{t \uparrow t_b} \\ &= \int_{t_a}^{t_b} \left( \frac{\partial L}{\partial \boldsymbol{q}} - \frac{d}{dt} \frac{\partial L}{\partial \boldsymbol{u}} \right) \delta \boldsymbol{q} \, dt + \left[ \frac{\partial L}{\partial \boldsymbol{u}} \delta \boldsymbol{q} + L \, \delta t \right]_{t \downarrow t_a}^{t \uparrow t_b}. \end{aligned} \qquad (70)$$

The boundary terms in (70) are due to the variation of the begin time $t_a$ and end time $t_b$ (see Fig. 5), and are expressed with left and right limits because the velocity $\boldsymbol{u}$ and, therefore, $L$ may not exist for $t = t_a$ or $t = t_b$. Note that $\delta q(t)$ is continuous in the open interval $(t_a, t_b)$. Using the generalised momentum $\boldsymbol{p}$ (which equals (52)) and the Hamiltonian function $H$, i.e.

$$\boldsymbol{p}(\boldsymbol{q}, \boldsymbol{u}) = \left( \frac{\partial L}{\partial \boldsymbol{u}} \right)^{\mathrm{T}}, \qquad H(\boldsymbol{q}, \boldsymbol{u}) = \frac{\partial L}{\partial \boldsymbol{u}} \boldsymbol{u} - L = T + V, \qquad (71)$$
<figure><img src="image_7.png" /><figcaption>**Fig. 6** The general variation of a function $q(t)$ with kink</figcaption></figure>

we can express the boundary terms in (70) as

$$ \lim_{t \downarrow t_a} \left( \frac{\partial L}{\partial \mathbf{u}} \delta \mathbf{q} + L \delta t \right) = \mathbf{p}(\mathbf{q}_a, \mathbf{u}_a^+) \mathbf{u}^+ (\delta \mathbf{q}_a^+ + \mathbf{u}_a^+ \delta t_a) - \mathbf{H}(\mathbf{q}_a, \mathbf{u}_a^+) \delta t_a, \quad (72) $$

$$ \lim_{t \uparrow t_b} \left( \frac{\partial L}{\partial \mathbf{u}} \delta \mathbf{q} + L \delta t \right) = \mathbf{p}(\mathbf{q}_b, \mathbf{u}_b^-) \mathbf{u}^- (\delta \mathbf{q}_b^- + \mathbf{u}_b^- \delta t_b) - \mathbf{H}(\mathbf{q}_b, \mathbf{u}_b^-) \delta t_b, \quad (73) $$

in which the abbreviation $\delta \mathbf{q}^+(t_a) = \delta \mathbf{q}_a^+$, $\mathbf{u}^+(t_a) = \mathbf{u}_a^+$, etc. has been used. We recognise in (72) and (73) the variation $\overline{\delta \mathbf{q}}(t)$ of the begin and endpoint

$$ \begin{aligned} \overline{\delta \mathbf{q}_a} &= \delta \mathbf{q}_a^+ + \mathbf{u}_a^+ \delta t_a, \\ \overline{\delta \mathbf{q}_b} &= \delta \mathbf{q}_b^- + \mathbf{u}_b^- \delta t_b, \end{aligned} \quad (74) $$

see Fig. 5. The total variation (70) simplifies to

$$ \delta s = \int_{t_a}^{t_b} \left( \frac{\partial L}{\partial \mathbf{q}} - \frac{d}{dt} \frac{\partial L}{\partial \mathbf{u}} \right) \delta \mathbf{q} \, dt - \left( \mathbf{p}_a^+ \right)^\mathrm{T} \overline{\delta \mathbf{q}_a} + \left( \mathbf{p}_b^- \right)^\mathrm{T} \overline{\delta \mathbf{q}_b} + \mathbf{H}_a^+ \delta t_a - \mathbf{H}_b^- \delta t_b, \quad (75) $$

where the notation $\mathbf{p}_a^+ = \mathbf{p}(\mathbf{q}_a, \mathbf{u}_a^+)$, $\mathbf{H}_a^+ = \mathbf{H}(\mathbf{q}_a, \mathbf{u}_a^+)$, etc. has been used.

In the following, we will consider solution curves $\mathbf{q}(t)$ which may have a kink at some point in time $t_c \in I$ (see Fig. 6). The action integral, which is a Lebesgue integral over an interval $I$, can be decomposed into two differentiable parts

$$ s(\mathbf{q}) = \int_I L \, dt = \int_{I_1} L \, dt + \int_{I_2} L \, dt = s_1(\mathbf{q}) + s_2(\mathbf{q}), \quad (76) $$

where $I_1 = [t_0, t_c]$ and $I_2 = [t_c, t_f]$. The end points $t = t_0$ and $t = t_f$ are fixed and we require that the two differentiable parts of the function $\mathbf{q}(t)$ join continuously at $t = t_c$, but otherwise the point $t = t_c$ can move freely. The comparison function $\hat{\mathbf{q}}(\varepsilon, t)$ is fixed at the end points $t = t_0$ and $t = t_f$ and consists of two differentiable parts which join continuously at $t = t_c + \delta t_c$. The function $\overline{\delta \mathbf{q}}(t)$ is therefore continuous
at $t = t_c$, whereas the variation $\delta q(t)$ is discontinuous as $t = t_c$. The variation of the impact position $\overline{\delta q}_c = \overline{\delta q}(t_c)$ can be assessed from the left and from the right which gives the equality

$$ \overline{\delta q}_c = \delta \overline{q}_c^{-} + \overline{u}_c^{-} \delta t_c = \delta \overline{q}_c^{+} + \overline{u}_c^{+} \delta t_c, \quad (77) $$

with the notation $\delta q^{\pm}(t_c) = \delta \overline{q}_c^{\pm}$; see Fig. 6. It holds that $-\overline{u}_c^{-} \in \mathcal{T}_K(q(t_c))$ and $\overline{u}_c^{+} \in \mathcal{T}_K(q(t_c))$. Therefore, if $\delta t_c < 0$ and $\delta \overline{q}_c^{-} \in \mathcal{T}_K(q(t_c))$, we have $\overline{\delta q}_c = \delta \overline{q}_c^{-} + \overline{u}_c^{-} \delta t_c \in \mathcal{T}_K(q(t_c))$, whereas if $\delta t_c > 0$ and $\delta \overline{q}_c^{+} \in \mathcal{T}_K(q(t_c))$ we have $\overline{\delta q}_c = \delta \overline{q}_c^{+} + \overline{u}_c^{+} \delta t_c \in \mathcal{T}_K(q(t_c))$. Hence, if $\delta \overline{q}_c^{\pm} \in \mathcal{T}_K(q(t_c))$, then it also holds that $\overline{\delta q}_c \in \mathcal{T}_K(q(t_c))$.

We introduce the functions

$$ h_1(\varepsilon) = \int_{t_0}^{t_c(\varepsilon)} L(\hat{q}(\varepsilon, t), \hat{u}(\varepsilon, t)) dt, \quad h_2(\varepsilon) = \int_{t_c(\varepsilon)}^{t_f} L(\hat{q}(\varepsilon, t), \hat{u}(\varepsilon, t)) dt \quad (78) $$

such that $h(\varepsilon) = h_1(\varepsilon) + h_2(\varepsilon)$ and take the variation of the action integral $s$ by calculating the $\varepsilon$-derivatives $h'_1(\varepsilon)$ and $h'_2(\varepsilon)$ separately. In other words, the variation $\delta s_1$ has a variable end-point whereas $\delta s_2$ has a variable starting-point. Evaluation of $h'_1(\varepsilon)$ and $h'_2(\varepsilon)$ gives

$$ \delta s_1 = h'_1(0)\delta\varepsilon = \int_{I_1} \left( \frac{\partial L}{\partial \boldsymbol{q}} - \frac{d}{dt} \frac{\partial L}{\partial \boldsymbol{u}} \right) \delta \boldsymbol{q} \, dt + (\boldsymbol{p}_c^{-})^{\mathrm{T}} \overline{\delta \boldsymbol{q}}_c - H_c^{-} \delta t_c, \quad (79) $$

$$ \delta s_2 = h'_2(0)\delta\varepsilon = \int_{I_2} \left( \frac{\partial L}{\partial \boldsymbol{q}} - \frac{d}{dt} \frac{\partial L}{\partial \boldsymbol{u}} \right) \delta \boldsymbol{q} \, dt - (\boldsymbol{p}_c^{+})^{\mathrm{T}} \overline{\delta \boldsymbol{q}}_c + H_c^{+} \delta t_c. \quad (80) $$

Addition of $h'_1(0)$ and $h'_2(0)$ yields $h'(0)$. The variation of $s$ can be written as

$$ \delta s = \int_I \left( \frac{\partial L}{\partial \boldsymbol{q}} - \frac{d}{dt} \frac{\partial L}{\partial \boldsymbol{u}} \right) \delta \boldsymbol{q} \, dt - (\boldsymbol{p}_c^{+} - \boldsymbol{p}_c^{-})^{\mathrm{T}} \overline{\delta \boldsymbol{q}}_c + (H_c^{+} - H_c^{-}) \delta t_c, \quad (81) $$

where $\boldsymbol{p}_c^{\pm}$ and $H_c^{\pm}$ are the pre- and post-impact values of the canonical variables.

We now return to the question of finding conditions for which the principle of Hamilton in inequality form (58) also holds for a time-interval $I$ which contains an impact event. If (58) holds, in other words if $-\delta s \ge 0$ for all $\delta \boldsymbol{q} \in \mathcal{T}_K(\boldsymbol{q})$, then (81) yields

$$ -\int_I \left( \frac{\partial L}{\partial \boldsymbol{q}} - \frac{d}{dt} \frac{\partial L}{\partial \boldsymbol{u}} \right) \delta \boldsymbol{q} \, dt + (\boldsymbol{p}_c^{+} - \boldsymbol{p}_c^{-})^{\mathrm{T}} \overline{\delta \boldsymbol{q}}_c - (H_c^{+} - H_c^{-}) \delta t_c \ge 0 \quad (82) $$

for all $\delta \boldsymbol{q}(t) \in \mathcal{T}_K(\boldsymbol{q}(t))$, where $t \in I \setminus \{t_c\}$, for all $\overline{\delta \boldsymbol{q}}_c \in \mathcal{T}_K(\boldsymbol{q}(t_c))$ and for all $\delta t_c$. The first term, $-\int_I (\frac{\partial L}{\partial \boldsymbol{q}} - \frac{d}{dt} \frac{\partial L}{\partial \boldsymbol{u}}) \delta \boldsymbol{q} \, dt$, involves the non-impulsive motion and is clearly non-negative for admissible virtual displacements; see (60). This analysis suggests to make the following two suppositions about the nature of the impact:

1. The supposition of a perfect constraint remains valid during the impact. This implies that the change in momentum must be normal to the constraint in the sense that

$$ -(\boldsymbol{p}_c^{+} - \boldsymbol{p}_c^{-}) \in \mathcal{N}_K(\boldsymbol{q}(t_c)), \quad (83) $$
which is equivalent to

$$ \overline{\delta \boldsymbol{q}_c}^{\mathrm{T}} (\boldsymbol{p}_c^{+} - \boldsymbol{p}_c^{-}) \geq 0 \quad \forall \overline{\delta \boldsymbol{q}_c} \in \mathcal{T}_K (\boldsymbol{q}(t_c)). \qquad (84) $$

2. The impact is assumed to preserve the total energy in the system, i.e.

$$ H_c^+ - H_c^- = 0. \qquad (85) $$

These assumptions are for the unconstrained case ($K = \mathbb{R}^n$, $\mathcal{N}_K(\boldsymbol{q}(t_c)) = \mathbf{0}$) exactly the so-called first and the second Weierstrass–Erdmann corner conditions. Given a functional $J(y) = \int_I f(x, y, y') dx$, the first Weierstrass–Erdmann condition reads as

$$ \lim_{x \downarrow x_c} f_{y'} = \lim_{x \uparrow x_c} f_{y'}, \quad \text{or} \quad \boldsymbol{p}_c^+ = \boldsymbol{p}_c^-, \qquad (86) $$

and the second Weierstrass–Erdmann condition reads as

$$ \lim_{x \downarrow x_c} f - y' f_{y'} = \lim_{x \uparrow x_c} f - y' f_{y'}, \quad \text{or} \quad H_c^+ = H_c^-. \qquad (87) $$

The suppositions (83) and (85) lead to the following theorem (see Panagiotopoulos and Glocker 1998).

**Theorem 3** (The strong principle of Hamilton in inequality form) *Consider a conservative Lagrangian mechanical system with perfect unilateral constraints and a non-dissipative impact law. A function $\boldsymbol{q}(t) \in K$ is a motion of the system if and only if it is a strong local extremal of the action integral, i.e.*

$$ -\delta \int_I L dt \geq 0 \quad \forall \delta \boldsymbol{q} \in \mathcal{T}_K(\boldsymbol{q}) \text{ a.e. on } I, \quad \boldsymbol{q}(t_0) = \boldsymbol{q}_0, \quad \boldsymbol{q}(t_f) = \boldsymbol{q}_f, \qquad (88) $$

which is the strong principle of Hamilton in inequality form for impulsive motion.

*Proof* The supposition of a perfect unilateral constraint implies the normality condition $-(\boldsymbol{p}^{+} - \boldsymbol{p}^{-}) \in \mathcal{N}_K(\boldsymbol{q})$. The supposition of a non-dissipative impact implies the energy condition $H^{+} = H^{-}$. Hence, the conditions (83) and (85) are satisfied. If $\boldsymbol{q}(t)$ is a strong local extremal of the action integral then, by definition, $-\delta s \geq 0$ for all strong variations $\delta \boldsymbol{q} \in \mathcal{T}_K$ (88). Using (81), we obtain the variational problem (82). The variations $\delta \boldsymbol{q}(t)$, $\overline{\delta \boldsymbol{q}}_c$ and $\delta t_c$ are independent. Hence they yield two variational inequalities and one variational equality:

$$ - \int_I \left( \frac{\partial L}{\partial \boldsymbol{q}} - \frac{d}{dt} \frac{\partial L}{\partial \dot{\boldsymbol{q}}} \right) \delta \boldsymbol{q} \, dt \geq 0 \quad \forall \delta \boldsymbol{q}(t) \in \mathcal{T}_K(\boldsymbol{q}(t)), \quad \text{almost everywhere on } I, \qquad (89) $$

$$ (\boldsymbol{p}_c^+ - \boldsymbol{p}_c^-)^{\mathrm{T}} \overline{\delta \boldsymbol{q}_c} \geq 0 \quad \forall \overline{\delta \boldsymbol{q}_c} \in \mathcal{T}_K(\boldsymbol{q}(t_c)), \qquad (90) $$

$$ (H_c^+ - H_c^-) \delta t_c = 0 \quad \forall \delta t_c. \qquad (91) $$

From (89) we see that the Euler–Lagrange inequality (60) holds for almost all $t \in I$. Equation (90) requires that the unilateral constraints are perfect (in the sense of impulsive motion) and (91) that the impacts are non-dissipative. Consequently, if $\boldsymbol{q}$ is a
strong extremal of the action integral, then $q(t)$ is a motion of the system. The proof can easily be followed in the opposite direction. $\square$

## 7 A Weak Principle of Hamilton

In this section we derive a weak principle of Hamilton in inequality form for impulsive motion by directly incorporating the impulsive dynamics in the principle of virtual work. Concepts of measure and integration theory appear to be very useful in this respect.

As before, let $\xi(X, t)$ describe the motion of a mass element $dm$ with position $X$ in the reference configuration. The configuration at time $t$ is denoted by $\xi(\cdot, t)$ but we will simply write $\xi(t)$ and denote the velocity field by $v(t)$. We assume that the motion $\xi(t)$ is absolutely continuous in time and that the velocity field $v(t)$ is a function (field) of locally bounded variation without singular terms (which is sometimes called 'special bounded variation'). This implies the following:

1. At each time-instant $t$ we can define a left and right velocity field

$$ 
\begin{aligned}
v^+(t) &= \lim_{\tau \downarrow 0} \frac{\xi(t + \tau) - \xi(t)}{\tau}, \\
v^-(t) &= \lim_{\tau \uparrow 0} \frac{\xi(t + \tau) - \xi(t)}{\tau}.
\end{aligned} 
\quad (92) 
$$

2. The differential measure $d\xi$ of the configuration $\xi(t)$ contains only a density $v$ with respect to the Lebesgue measure $dt$

$$ d\xi = v \, dt. \quad (93) $$

For almost all $t$, we can define a velocity field $v(t) = \dot{\xi}(t)$.

3. The differential measure $dv$ of the velocity field $v(t)$ contains a density with respect to the Lebesgue measure $dt$ and with respect to the atomic measure $d\eta$, i.e.

$$ dv = \dot{v} \, dt + (v^+ - v^-) \, d\eta. \quad (94) $$

Let as before $dR(X, t)$ be the constraint force, associated with bilateral and unilateral geometric scleronomic constraints, and $dF(X, t)$ be the external force on a mass element $dm$ with position $\xi(X, t)$. The non-impulsive dynamics of the system is therefore described by the equation of motion

$$ dm \ddot{\xi} - dF - dR = 0. \quad (95) $$

The impulsive dynamics is described by the impact equation

$$ dm (v^+ - v^-) - dP = 0, \quad (96) $$

where $dP$ is the impulsive force field of the constraints. The introduction of the differential measure $dv$ allows us to combine the equation of motion (95) and the impact
equation (96) in a single equation

$$ dm \, dv - dF \, dt - dR \, dt - dP \, d\eta = 0, \qquad (97) $$

which is an equality of measures and which should be understood in the sense of integration. The measure $d\nu$, which has a density with respect to the atomic measure, is by definition a mapping on the space of continuous functions, i.e. $\int_I f^T d\nu$ only makes sense if $f \in C^0(I, \mathbb{R}^3)$. The equality of measures leads to a principle of virtual work in differential measures

$$ \int_S \delta\xi^T (dm \, d\nu - dF \, dt - dR \, dt - dP \, d\eta) = 0 \qquad (98) $$

for all $\delta\xi \in \mathcal{Y}_0(I, \mathbb{R}^3)$. The class $\mathcal{Y}(I, \mathbb{R}^3)$, see (14), comprises all functions on the domain $I$ which are absolutely continuous and piecewise $C^1$ as well as absolutely continuous functions with an accumulation point. We put $\hat{\xi}(\varepsilon, t) = \xi(t) + \varepsilon\omega(t)$ with $\omega \in \mathcal{Y}_0(I, \mathbb{R}^3)$. For this family of comparison functions, it holds that $\delta\xi = \hat{\xi}_\varepsilon(0) \delta\varepsilon = \omega \delta\varepsilon$ is continuous in time whereas $\delta\nu = \dot{\omega} \delta\varepsilon$ is not. Moreover, note that $\dot{\xi}(\varepsilon)$ converges for $\varepsilon \downarrow 0$ to $\xi$ in the weak norm $\| \cdot \|_1$; see (25). We now make the assumption that both the constraint force $dR$ and the constraint impulse $dP$ are perfect constraint forces/impulses, i.e.

$$ -dR \in \mathcal{N}_{\mathcal{K}}(\xi), \quad -dP \in \mathcal{N}_{\mathcal{K}}(\xi). \qquad (99) $$

The virtual work of the constraint force $dR$ and the constraint impulse $dP$ is therefore non-negative for admissible virtual displacements

$$ \delta\xi^T dR \ge 0, \quad \delta\xi^T dP \ge 0 \quad \forall \delta\xi \in \mathcal{T}_{\mathcal{K}}(\xi). \qquad (100) $$

This leads us to the principle of d'Alembert-Lagrange in inequality form for differential measures

$$ \int_S \delta\xi^T (dm \, d\nu - dF \, dt) \ge 0 \quad \forall \delta\xi = \omega \delta\varepsilon \in \mathcal{T}_{\mathcal{K}}(\xi) \qquad (101) $$

in which we explicitly write that the virtual displacements $\delta\xi$ are of the form $\omega \delta\varepsilon$. A mass element has a constant mass $dm$. Taking the differential measure-in-time of the term $\delta\xi^T dm \, \nu$ and applying the chain rule gives

$$ d(\delta\xi^T dm \, \nu) = \delta\nu^T dm \, \nu \, dt + \delta\xi^T dm \, d\nu, \qquad (102) $$

in which we used the equality $d(\delta\xi) = d\omega \delta\varepsilon = \dot{\omega} dt \delta\varepsilon = \delta\nu dt$. We therefore arrive at the variational inequality

$$ d\left(\int_S \delta\xi^T dm \, \nu\right) - \int_S \delta\nu^T dm \, \nu \, dt - \int_S \delta\xi^T dF \, dt \ge 0 \qquad (103) $$

$\forall \delta\xi = \omega \delta\varepsilon \in \mathcal{T}_{\mathcal{K}}(\xi)$. We recognise the second term $\delta T = \int_S \delta\nu^T dm \, \nu$ as being the variation of the kinetic energy $T = \int_S \frac{1}{2} \nu^T dm \, \nu$. In the same way as before, we choose
a chart $\varphi$ with generalised coordinates $q(t)$ which form a minimal set of coordinates with respect to the bilateral constraints and with which we can uniquely describe the configuration $\xi(t) = \varphi(q(t))$. Moreover, we introduce generalised velocities $u(t)$, which are assumed to be of locally bounded variation, and which are such that $d\boldsymbol{q} = u \, dt$. The family of comparison functions $\hat{\xi}(\varepsilon, t) = \xi(t) + \varepsilon\omega(t)$ is generated by the family of comparison functions $\hat{q}(\varepsilon, t) = q(t) + \varepsilon w(t)$. Hence, the variation $\delta\boldsymbol{q} = w \, \delta\varepsilon$ induces the variation $\delta\xi = \partial\varphi/\partial\boldsymbol{q} \, \delta\boldsymbol{q}$ with $\omega = \partial\varphi/\partial\boldsymbol{q} \, w$. The kinetic energy $T$ is a function of $v(q, u)$ and we can therefore write $T$ as a function $T(q, u)$. Hence, it holds that $\delta T = T_q \, \delta q + T_u \, \delta u$. The variations $\delta q(t) = w(t) \, \delta\varepsilon$ is a continuous function in time whereas $\delta u(t) = \dot{w}(t) \, \delta\varepsilon$ is discontinuous. Using the generalised momentum

$$ \boldsymbol{p}(q, u) = \int_S \left( \frac{\partial \varphi}{\partial \boldsymbol{q}} \right)^T \boldsymbol{v} \, dm \quad (104) $$

and the generalised force (54) we transform the principle of d'Alembert-Lagrange into

$$ d(\delta \boldsymbol{q}^T \boldsymbol{p}) - \delta T \, dt - \delta \boldsymbol{q}^T f \, dt \geq 0 \quad \forall \delta \boldsymbol{q} = \boldsymbol{w} \, \delta \varepsilon \in \mathcal{T}_K(\boldsymbol{q}), \quad (105) $$

which is the Lagrange central equation in differential measures. As before, we integrate over a time-interval $I = [t_0, t_f]$ and consider the generalised force $f = -V'(q)$ to be a potential force, which yields

$$ [\delta \boldsymbol{q}^T \boldsymbol{p}]_{t \uparrow t_0}^{t \downarrow t_f} - \int_I \delta L \, dt \geq 0 \quad \forall \delta \boldsymbol{q} = \boldsymbol{w} \, \delta \varepsilon \in \mathcal{T}_K(\boldsymbol{q}), \quad (106) $$

where $L = T - V$. Finally, taking fixed boundary conditions at $t_0$ and $t_f$, i.e. $\delta q(t_0) = \delta q(t_f) = 0$, we obtain

$$ -\delta \int_I L \, dt \geq 0 \quad \forall \delta \boldsymbol{q} = \boldsymbol{w} \, \delta \varepsilon \in \mathcal{T}_K(\boldsymbol{q}). \quad (107) $$

The comparison functions $\hat{q}(\varepsilon, t) = q(t) + \varepsilon w(t)$ converge in the weak norm $\| \cdot \|_1$ and lead to weak variations $\delta q$. The variation $\delta s$ of the action $s(q) = \int_I L \, dt$ in (107) reduces therefore to the Gâteaux derivative $ds(q; \delta q)$; see (22). Consequently, we arrive at a weak form of the principle of Hamilton in inequality form, i.e. the principle of Hamilton as weak local extremal of the action $s(q)$ with the weak norm $\| \cdot \|_1$.

**Theorem 4** (The weak principle of Hamilton in inequality form) *Consider a conservative Lagrangian mechanical system with perfect unilateral constraints. A function $q(t) \in K$ is a motion of the system if and only if it is a weak local extremal of the action integral $s(q) = \int_I L \, dt$, i.e.*

$$ -ds(q; \delta q) \geq 0 \quad \forall \delta q = w \, \delta \varepsilon \in \mathcal{T}_K(q), \quad q(t_0) = q_0, \quad q(t_f) = q_f, \quad (108) $$

with $w \in Y_0(I, \mathbb{R}^n)$, which is the weak principle of Hamilton in inequality form for impulsive motion.
*Proof* If $\boldsymbol{q}(t)$ is a motion of the system, or equivalently $\boldsymbol{\xi}(t) = \boldsymbol{\varphi}(\boldsymbol{q}(t))$, then it has to fulfill the equations of motion (95) for almost all $t$ and the impact equation (96) at the impact time-instants. This leads to the principle of virtual work in differential measures (98). The supposition of perfect unilateral constraint forces and impulses leads to the principle of Hamilton (107) with weak variations $\delta\boldsymbol{q}$. The proof can easily be followed in the opposite direction. $\square$

There are two interesting things to remark at this point. First of all, by making use of differential measures, we are able to treat the impulsive and non-impulsive dynamics simultaneously. This means that we do not need to split the time-interval in differentiable parts to treat an impact as has been done in Sect. 6. Secondly, as we want to use the principle of virtual work for differential measures (98), we are forced to consider comparison functions $\hat{\boldsymbol{\xi}}(\varepsilon, t) = \boldsymbol{\xi}(t) + \varepsilon\boldsymbol{\omega}(t)$ of which the variation $\delta\boldsymbol{\xi} = \boldsymbol{\omega}\delta\varepsilon$ is time-continuous. A comparison function of this form can have kinks (and, therefore, impacts), but this class of comparison functions does not include a family of curves which only varies the impact time. Accordingly, for this weak form of the principle of Hamilton, we only have to satisfy the 'generalised' first Weierstrass-Erdmann corner condition (83). Hence, we are able to prove the validity of Theorem 4 by only assuming that the unilateral geometric constraint is perfect which is expressed by (99). The supposition of energy conservation during the impact, i.e. the second Weierstrass-Erdmann corner condition $H_c^+ = H_c^-$, is therefore not a necessary condition for the validity of Theorem 4.

## 8 Example: Numerical Discretisation

The differences between the weak and strong principle of Hamilton can be illustrated by looking at their discretized forms, i.e. the action integral is discretized, similar to what is done in the Ritz method. As we will show, both principles lead to a numerical scheme for impulsive dynamics, but an impact law does (generally) not follow from the discretized principles. As we only want to convey the basic idea, the simplest example system is chosen: that of a one-dimensional free particle.

Consider a point-mass $m$ with coordinate $\boldsymbol{q}(t)$ which is subjected to the unilateral constraint $\boldsymbol{q}(t) \ge 0$. In the absence of gravitation or other external forces, the action integral reads as

$$s(\boldsymbol{q}) = \int_{t_0}^{t_f} T(u) \, dt, \quad T(u) = \frac{1}{2}mu^2, \quad u = \dot{q} \text{ a.e.} \qquad (109)$$

with given positions $\boldsymbol{q}(t_0) = q_0$ and $\boldsymbol{q}(t_f) = q_f$ at begin and end time. The evolution $\boldsymbol{q}(t)$ is approximated with approximants $q_i \approx \boldsymbol{q}(t_i)$ on a discretisation $t_0 < t_1 < t_2 < \dots < t_N < t_{N+1} = t_f$ of the time interval $[t_0, t_f]$. We can either choose this time-grid a priori, i.e. we choose a fixed grid, or the choice of $t_i$ for $i = 1, \dots, N$ has still to be determined. Correspondingly, we define approximants $u_i \approx u(t_i)$ and $T_i \approx T(u(t_i))$ for the velocity and kinetic energy at the discretisation point $t_i$:

$$u_i = \frac{q_{i+1} - q_i}{t_{i+1} - t_i}, \quad T_i = \frac{1}{2}mu_i^2 = \frac{1}{2}m \left( \frac{q_{i+1} - q_i}{t_{i+1} - t_i} \right)^2. \qquad (110)$$
Using the Euler quadrature $\int_0^l f(x) dx \approx \sum_{i=0}^N f(x_i)(x_{i+1} - x_i)$, the action integral $s(q)$ can be approximated by the sum

$$s^N(q_1, q_2, \dots, q_N) = \sum_{i=0}^{N} T_i (t_{i+1} - t_i). \qquad (111)$$

We will now consider a discretisation of the weak principle of Hamilton and the strong principle of Hamilton for this example system.

## 8.1 Weak Principle of Hamilton

If a fixed grid of time discretisation points $t_i$ is chosen, then we are only able to consider variations $\delta q_i$ of the approximants $q_i$, which can be looked upon as approximants of the variation at time $t_i$, i.e. $\delta q_i \approx \delta q(t_i)$, for fixed values of $t_i$. As time is not varied, these (approximate) variations are weak and we are able to discretise the weak principle of Hamilton in inequality form (108). The Gâteaux derivative of the action integral appearing in (108) can be approximated by

$$ds(q, \delta q) \approx ds^N(q_1, \dots, q_N; \delta q_1, \dots, \delta q_N) = \sum_{i=1}^{N} \frac{\partial s^N}{\partial q_i} \delta q_i, \qquad (112)$$

with

$$\frac{\partial s^N}{\partial q_i} = \frac{\partial T_{i-1}}{\partial q_i}(t_i - t_{i-1}) + \frac{\partial T_i}{\partial q_i}(t_{i+1} - t_i) = -m(u_i - u_{i-1}). \qquad (113)$$

The discretized weak principle of Hamilton therefore leads to the difference scheme in the form of a finite dimensional variational inequality

$$m(u_i - u_{i-1})\delta q_i \ge 0 \quad \forall \delta q_i \in \mathcal{T}_{\mathbb{R}_0^+}(q_i), \qquad (114)$$

which can be cast in the inclusion

$$m(u_i - u_{i-1}) - \lambda_i = 0, \quad -\lambda_i \in \mathcal{N}_{\mathbb{R}_0^+}(q_i), \qquad (115)$$

where $\lambda_i$ is the approximant for the percussion of the constraint during the time interval $(t_{i-1}, t_i]$. Note that the inclusion in (115) merely requires the inequality complementarity between the discrete percussion $\lambda_i$ and the position $q_i$ and does not include an impact law. The difference scheme (115) together with the boundary conditions $q_0$ and $q_{N+1} = q_f$ has therefore $N$ solutions, each with a single impact at some time discretisation point $t_j$ ($j = 1, \dots, N$),

$$q_i = \begin{cases} \frac{t_j - t_i}{t_j - t_0} q_0, & i \le j, \\ \frac{t_i - t_j}{t_f - t_j} q_f, & i > j, \end{cases} \qquad (116)$$

and one solution at constant velocity without impact $q_i = q_0 + \frac{t_i - t_0}{t_f - t_0}(q_f - q_0)$; see Fig. 7 in which all solutions are drawn. Each of these solutions are approximants of weak local extremals of the action integral.
<figure><img src="image_8.png" /><figcaption>**Fig. 7** Weak and strong solutions (in bold) of the discretized principle of Hamilton</figcaption></figure>

## 8.2 Strong Principle of Hamilton

If we consider a variable grid of time discretisation points $t_i$, then variations in time $\delta t_i$ induce the variations $\delta q(t_i) \approx u_i \delta t_i + \delta q_i$ which lead to strong variations. The strong principle of Hamilton (88) requires the sub-stationarity $-\delta s \ge 0$ of the action integral $s(q)$ with respect all strong variations $\delta q \in T_{\mathbb{R}_0^+}(q)$. The discretized form of the strong principle of Hamilton therefore reads as

$$ -\delta s^N = \sum_{i=1}^{N} \frac{\partial s^N}{\partial q_i} \delta q_i + \frac{\partial s^N}{\partial t_i} \delta t_i \ge 0 \quad \forall \delta q_i \in T_{\mathbb{R}_0^+}(q_i), \forall \delta t_i, \quad (117) $$

with independent variations $\delta q_i$ and $\delta t_i$. Hence, the partial derivatives

$$ \frac{\partial s^N}{\partial t_i} = T_i - T_{i-1} \quad (118) $$

have to vanish for $i = 1, \dots, N$, which implies the conservation of energy during each discrete time-interval. The discretized strong principle of Hamilton in equality form therefore leads to the difference scheme

$$ \begin{gathered} m(u_i - u_{i-1}) - \lambda_i = 0, \quad -\lambda_i \in \mathcal{N}_{\mathbb{R}_0^+}(q_i), \\ \frac{1}{2}mu_i^2 = \frac{1}{2}mu_{i-1}^2. \end{gathered} \quad (119) $$

The numerical scheme (119) with boundary conditions $q_0$ and $q_{N+1} = q_f$ has an impacting solution (116) with $t_j = t_c := (q_0 t_f + q_f t_0)/(q_f + q_0)$, whereas $t_i$ is free to choose for $i \ne j$, as well as the non-impacting solution $q_i = q_0 + \frac{t_i - t_0}{t_i - t_0}(q_f - q_0)$ with $t_i$ free to choose (see Fig. 7). Note that the discretized strong principle of Hamilton (119) requires the perfectness of the constraint and the conservation of energy, which implies in the one-dimensional case the completely elastic impact law $u_j = -u_{j-1}$.

From the above example, it becomes clear that the weak principle of Hamilton can be used to derive numerical integration schemes for impulsive mechanical systems on
a fixed grid. The impact law has still to be incorporated in the numerical scheme, for instance by relaxing the unilateral constraint $-\lambda_i \in \mathcal{N}_{\mathbb{R}_0^+}(q_i)$ on position level to the unilateral constraint

$$ -\lambda_i \in \mathcal{N}_{\mathbb{R}_0^+}(q_i) \left( \frac{1}{1+\varepsilon} u_i + \frac{\varepsilon}{1+\varepsilon} u_{i-1} \right), \quad (120) $$

on velocity level, which is a Newton-type of impact law with a coefficient of restitution $\varepsilon$ (Moreau 1988). More elegantly, one might try to incorporate an impact law in the weak form of Hamilton's principle by requiring the sub-stationarity of the (discretized) action integral with respect to an appropriate larger set of comparison functions, but this is still a topic of future research.

# 9 Conclusions

The aim of this paper is to clarify the existing forms of the principle of Hamilton as variational inequality, which have already been derived in Panagiotopoulos and Glocker (1998, 2000), by putting them within the context of weak and strong extrema.

In this paper, we have derived two different forms of Hamilton's principle as variational inequality. The strong form of Hamilton's principle has been derived in Sect. 6 using the general variation of the functional which leads to 'generalised' Weierstrass-Erdmann corner conditions. The first 'generalised' Weierstrass-Erdmann condition demands that the contact impulses are from the normal cone, i.e. the supposition of perfect unilateral constraints. The second Weierstrass-Erdmann condition requires that the collisions are completely elastic, i.e. there is no energy loss during the impact. This form of the principle of Hamilton takes all neighbouring functions into consideration for the sub-stationarity of the solution and is therefore a strong form of Hamilton's principle as variational inequality. A weak form of Hamilton's principle has been derived in Sect. 7 starting from the principle of virtual work in differential measures. The latter can only be posed for a special class of virtual displacements, namely those which are related to Gâteaux derivatives. This derivation therefore leads in a natural way to a weak form of Hamilton's principle. The weak form of Hamilton's principle only requires that the contact impulses are from the normal cone, which is a generalisation of the first Weierstrass-Erdmann condition. No restriction is posed in the weak form on the energy dissipation of the impact law.

The various forms of the principle of Hamilton as variational inequality which exist in literature can now be put within the context of weak and strong extrema. We conclude that Theorem 3 of Panagiotopoulos and Glocker (1998) is the weak form of Hamilton's principle, while Theorem 4 of Panagiotopoulos and Glocker (1998) (or Proposition 4 of Panagiotopoulos and Glocker 2000) is the strong form of Hamilton's principle. This insight clarifies why the various principles have different conditions on the impact law. The forms of the principle of Hamilton used by Tornambè (1999) and Percivale (1985, 1991) assume strong variations and the authors of these works therefore state that the principles are valid for completely elastic impact.

In Sect. 8, it has very briefly been shown on an example system that a direct discretisation of the action integral leads to a difference scheme, which can form the
basis of a numerical integration scheme. Future research will focus on the natural derivation of dedicated numerical schemes for the simulation of mechanical systems with frictionless unilateral contact by using a direct discretisation of the variational principles derived in this paper.

What is the practical/theoretical relevance of the results and insight gained in this paper? An obvious application has already been mentioned: the development of numerical schemes through a Ritz-type of method on the principle of Hamilton in inequality form. However, in the opinion of the authors, the relevance of the paper is more fundamental. One way to think about dynamics is in terms of variational principles. History proved that this way of thinking has been very rewarding. Variational principles form the foundation of classical analytical mechanics and have been essential for the development of modern physics, e.g. quantum mechanics. Furthermore, variational principles give the link to optimisation theory and put dynamics in an appropriate mathematical framework. Classical variational principles, however, are strictly valid for perfect bilateral constraints. For this reason, if we endeavour to develop a proper theoretical foundation for non-smooth dynamics, it is a promising step to go back to these principles and reformulate them in terms of variational inequalities.

## Appendix: Propositions

**Proposition 1** Let $f: I \to \mathbb{R}$ be a piece-wise continuous function and $\{g(\varepsilon)\}_{\varepsilon \ge 0}$ be a family of piece-wise continuous functions $g(\varepsilon): I \to \mathbb{R}$ with $g(0) = f$ almost everywhere. Denote the discontinuity points of $f$ with the set $S_f \subset I$ and the discontinuity points of $g(\varepsilon)$ with the set $S_{g(\varepsilon)} \subset I$. If it holds that

$$ \ess \sup_{x \in I} |f(x) - g(\varepsilon, x)| \xrightarrow{\varepsilon \downarrow 0} 0, \qquad (121) $$

then there exists a number $\varepsilon^* > 0$ such that

$$ \forall \varepsilon < \varepsilon^* \quad S_f \subseteq S_{g(\varepsilon)} \subset I. \qquad (122) $$

*Proof* We will prove the proposition with a reductio ad absurdum. Assume therefore that (121) holds and that there does not exist a number $\varepsilon^* > 0$ such that (122) holds. In other words, for all $\varepsilon^* > 0$, there exists $\varepsilon < \varepsilon^*$ such that $x_s$ is a discontinuity point of $f$ which is not a discontinuity point of $g(\varepsilon)$, i.e. $\forall \varepsilon^* > 0 \exists \varepsilon < \varepsilon^*$: $x_s \in S_f$ and $x_s \notin S_{g(\varepsilon)}$. Hence, the function $g(\varepsilon)$ is continuous on the open interval $(x_s - \delta, x_s + \delta)$, where $\delta(\varepsilon)$ is a function of $\varepsilon$. We will show that this assumption leads to a contradiction. According to (121), there exists $\varepsilon > 0$ for each essential upper bound $\alpha > 0$ such that

$$ \begin{cases} |f(x_s + \delta) - g(\varepsilon, x_s + \delta)| < \alpha, \\ |f(x_s - \delta) - g(\varepsilon, x_s - \delta)| < \alpha, \end{cases} \qquad (123) $$
with sufficiently small $\delta(\varepsilon)$. Because $x_s$ is a discontinuity point of $f$ it holds that

$$|f(x_s + \delta) - f(x_s - \delta)| \xrightarrow{\delta \to 0} s > 0, \qquad (124)$$

where $s$ is the jump height of the discontinuity of $f$ at $x_s$.

Take $\varepsilon^*$ to be an arbitrary strictly positive value. Let $\varepsilon_m$ and $\delta_m = \delta(\varepsilon_m)$ be such that $g(\varepsilon_m)$ is continuous on the interval $(x_s - \delta_m, x_s + \delta_m)$ and such that it holds that

$$\begin{aligned} & |f(x_s + \delta_m) - g(\varepsilon_m, x_s + \delta_m)| < \varepsilon^*, \\ & |f(x_s - \delta_m) - g(\varepsilon_m, x_s - \delta_m)| < \varepsilon^*. \end{aligned} \qquad (125)$$

Using a telescoping sum and the triangle inequality, we obtain the inequality

$$\begin{aligned} |f(x_s + \delta) - f(x_s - \delta)| &\le |f(x_s + \delta) - g(\varepsilon, x_s + \delta)| \\ &\quad + |g(\varepsilon, x_s + \delta) - g(\varepsilon, x_s - \delta)| \\ &\quad + |g(\varepsilon, x_s - \delta) - f(x_s - \delta)| \\ &< 2\varepsilon^* + |g(\varepsilon, x_s + \delta) - g(\varepsilon, x_s - \delta)| \end{aligned} \qquad (126)$$

for all $\delta < \delta_m$. Moreover, the local continuity of $g(\varepsilon)$ requires that the last term vanishes for $\delta \to 0$, whereas the left-hand side converges to the jump height $s$. Hence, we arrive at the inequality

$$s \le 2\varepsilon^*. \qquad (127)$$

Because $\varepsilon^* > 0$ is arbitrary, it follows that the jump height $s$ of $f$ is zero, which contradicts the assumption. This proves the proposition. $\square$

## References

* Aubin, J.-P., Frankowska, H.: Set-valued Analysis. Systems and Control: Foundations and Applications, vol. 2. Birkhäuser, Boston (1990)

* Bremer, H.: Dynamik und Regelung mechanischer Systeme. Teubner, Stuttgart (1988)

* Brogliato, B.: Nonsmooth Mechanics. Models, Dynamics and Control, 2 edn. Communications and Control Engineering. Springer, London (1999)

* Cesar, M.O.: Necessary conditions and sufficient conditions of weak minimum for solutions with corner points. Bol. Soc. Brasil. Mat. **15**(1-2), 109-135 (1984)

* Cesari, L.: Optimization—Theory and Applications. Applications of Mathematics, vol. 17. Springer, New York (1983)

* Cirak, F., West, M.: Decomposition contact response (DCR) for explicit finite element dynamics. Int. J. Numer. Methods Eng. **64**, 1078-1110 (2005)

* Fetecau, R.C., Marsden, J.E., Ortiz, M., West, M.: Nonsmooth Lagrangian mechanics and variational collision integrators. SIAM J. Appl. Dyn. Syst. **2**(3), 381-416 (2003)

* Glocker, Ch.: Set-Valued Force Laws, Dynamics of Non-Smooth Systems. Lecture Notes in Applied Mechanics, vol. 1. Springer, Berlin (2001)

* Goeleven, D., Panagiotopoulos, P.D., Lebau, C., Plotnikova, G.: Inequality forms of d'Alembert's principle in mechanics of systems with holonomic unilateral constraints. Z. Angew. Math. Mech. **77**(7), 483-501 (1997)

* Goeleven, D., Miettinen, M., Panagiotopoulos, P.D.: Dynamic hemivariational inequalities and their applications. J. Optim. Theory Appl. **103**(3), 567-601 (1999)
Goeleven, D., Motreanu, D., Dumont, Y., Rochdi, M.: Variational and Hemivariational Inequalities: Theory, Methods and Applications, Volume I: Unilateral Analysis and Unilateral Mechanics. Nonconvex Optimization and its Applications, vol. 69. Kluwer Academic, Dordrecht (2003a)

Goeleven, D., Motreanu, D., Dumont, Y., Rochdi, M.: Variational and Hemivariational Inequalities: Theory, Methods and Applications, Volume II: Unilateral Problems. Nonconvex Optimization and its Applications, vol. 70. Kluwer Academic, Dordrecht (2003b)

Hamel, G.: Elementare Mechanik. B.G. Teubner, Stuttgart (1912)

Hartman, P., Stampacchia, G.: On some nonlinear elliptic differential equations. Acta Math. **115**(1), 153–188 (1966)

Kinderlehrer, D., Stampacchia, G.: An Introduction to Variational Inequalities and their Applications. Pure and Applied Mathematics. Academic Press, New York (1980)

Kozlov, V.V., Treshchëv, D.V.: Billiards, a Genetic Introduction to the Dynamics of Systems with Impacts. Translations of Mathematical Monographs, vol. 89. American Mathematical Society, Providence (1991)

Kravchuk, A.S., Neittaanmäki, P.J.: Variational and Quasi-variational Inequalities in Mechanics. Solid Mechanics and its Applications, vol. 147. Springer, Dordrecht (2007)

Lanczos, C.: The Variational Principles of Mechanics. University of Toronto Press, Toronto (1962)

Leine, R.I., Nijmeijer, H.: Dynamics and Bifurcations of Non-Smooth Mechanical Systems. Lecture Notes in Applied and Computational Mechanics, vol. 18. Springer, Berlin (2004)

Leine, R.I., van de Wouw, N.: Stability and Convergence of Mechanical Systems with Unilateral Constraints. Lecture Notes in Applied and Computational Mechanics, vol. 36. Springer, Berlin (2008)

May, H.O.: Generalized variational principles and nondifferentiable potentials in analytical mechanics. J. Math. Phys. **25**(3), 491–493 (1984a)

May, H.O.: Hamilton's principle as substationarity principle. Acta Mech. **52**(3–4), 177–187 (1984b)

Moreau, J.J.: La notion de superpotentiel et les liaisons unilatérales en élastostatique. C. R. Acad. Sci., Sér. A **167** (1968), 954–957

Moreau, J.J.: Unilateral contact and dry friction in finite freedom dynamics. In: Moreau, J.J., Panagiotopoulos, P.D. (eds.) Non-Smooth Mechanics and Applications. CISM Courses and Lectures, vol. 302, pp. 1–82. Springer, Wien (1988)

Panagiotopoulos, P.D.: Non-convex superpotentials in the sense of F.H. Clarke and applications. Mech. Res. Commun. **8**(6), 335–340 (1981)

Panagiotopoulos, P.D.: Hemivariational Inequalities: Applications in Mechanics and Engineering. Springer, Berlin (1993)

Panagiotopoulos, P.D., Glocker, Ch.: Analytical mechanics: Addendum I. Inequality constraints with elastic impacts. The convex case. Z. Angew. Math. Mech. **78**(4), 219–229 (1998)

Panagiotopoulos, P.D., Glocker, Ch.: Inequality constraints with elastic impacts in deformable bodies. The convex case. Arch. Appl. Mech. **70**, 349–365 (2000)

Papastavridis, J.G.: Analytical Mechanics: A Comprehensive Treatise on the Dynamics of Constrained Systems: For Engineers, Physicists and Mathematicians. Oxford University Press, New York (2002)

Percivale, D.: Uniqueness in the elastic bounce problem I. J. Differ. Equ. **56**, 206–215 (1985)

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

Rockafellar, R.T., Wets, R.-B.: Variational Analysis. Springer, Berlin (1998)

Salomon, R.: Practical use of Hamilton's principle. J. Fluid Mech. **132**, 431–444 (1983)

Tornambè, A.: Modeling and control of impact in mechanical systems: Theory and experimental results. IEEE Trans. Autom. Control **44**(2), 294–309 (1999)

Troutman, J.L.: Variational Calculus and Optimal Control: Optimization with Elementary Convexity. Springer, New York (1996)