# Modeling and control of impact in mechanical systems: theory and experimental results

Antonio Tornambè

► To cite this version:

Antonio Tornambè. Modeling and control of impact in mechanical systems: theory and experimental results. IEEE Transactions on Automatic Control, 1999, 44 (2), pp.294 - 309. 10.1109/9.746255. hal-01375341

HAL Id: hal-01375341

https://hal.science/hal-01375341

Submitted on 3 Oct 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
# Modeling and Control of Impact in Mechanical Systems: Theory and Experimental Results

Antonio Tornambè

***Abstract***— This paper considers the equations of motion of mechanical systems subject to inequality constraints, which can be obtained by looking for the stationary value of the action integral. Two different methods are used to take into account the inequality constraints in the computation of the stationary value of the action integral: the method of the Valentine variables and the method of the penalty functions. The equations of motion resulting from the application of the method of the Valentine variables, which introduces the concept of “nonsmooth” impacts, constitute the exact model of the constrained mechanical system; such a model is suitable to be employed when the impacting parts of the actual mechanical system are very stiff. The equations of motion resulting from the application of the method of the penalty functions, which introduces the concept of “smooth impacts,” constitute an approximate model of the constrained mechanical system; such a model is suitable to be employed when the impacting parts of the actual mechanical system show some flexibility. Various feedback control laws from the natural outputs and from their time derivatives are studied with reference to both models of impact; the closed-loop systems resulting from the application of the same control law to both models show pretty much the same global asymptotic stability properties. The proposed control laws are only concerned with regulation problems in the presence of possible contacts and impacts among parts of the mechanical system or with the external environment; the problem of controlling these mechanical systems along time-varying trajectories is not considered in this paper. The effectiveness of the proposed control structure has been tested experimentally with reference to a single-link robot arm, showing a valuable behavior.

***Index Terms***— Feedback stabilization, mechanical systems, smooth and nonsmooth impacts, unilateral constraints.

## I. INTRODUCTION

MODELLING and control of impact require the study of the basic physical phenomena that happen when bodies collide under the action of external forces (including the control forces) and/or due to nonzero relative velocity. Several books (see, e.g., [1]–[12]) consider in detail, with a rigorous and extensive treatment, the study of the basic physical phenomena that attend the collision of bodies.

Often, the typical sources of impact, which include clearance as between cams and followers, backlash or bearing clearances in mechanisms undergoing force or motion reversal, and mechanisms with components having large relative velocities, are neglected in modeling and control of mechanisms. The

study and the control of the impact phenomena is especially important in robotics, since all the highest and greatest stresses arise as consequence of impact, and many serious failures can be generated when impact forces are not properly recognized and taken under control; some useful short-duration effects, such as high stresses, rapid dissipation of energy, and fast acceleration and deceleration, may be achieved from low-energy sources by controlling the impact of robots operating at low force levels.

The mechanisms are usually designed so that the impact effects are reduced, by minimizing the impact velocity and the mass of the impacting bodies, and by suitably designing the mechanism so that a minimum stiffness is located in a neighborhood of the point of impact, but these design factors seem to be useful for reducing the impact effects, not for their precise control. As a matter of fact, a real control task for a mechanism implies in general several transitions between the condition of free motion and the condition of constrained motion, and therefore it implies the generation of undesired reaction forces at each sudden change from one condition to the other one. As for robotics, practical situations involving the generation of impacts are the walking robots, the manipulation with a robotic hand, and the cooperation of multiple robot arms. Numerous attempts have been made in the recent years also to properly model the impact in robotics (see, e.g., [13]–[21]), whereas the problem of controlling the impact is still open, due to the sudden change of the equations of motion that happens when the bodies involved in the impact swish sharply from a condition of noncontact to a condition of contact. The reader interested in the derivation of the Hamilton principle for mechanical systems, subject to unilateral constraints, can benefit from [22]–[24]. Interesting initial experiments have been recently found in the literature [25]–[27], but no general theory is available because of the type of the equations to be used for adequately describing the impact of a robotic manipulator with the external environment, which depends largely upon the geometry of the robot and of the external environment and upon the type of impact. A complementarity–slackness class of hybrid systems has been analyzed in [28], under a certain dissipation impact rule. A complete robotic task has been formulated in [29] by means of some complementarity conditions in case of rigid external environments and in [30] in case of compliant external environments, under the assumption that the interaction force is available for feedback during the transition phase (this will not be required in this paper). The Lyapunov second method has been extended in [31] to deal with mechanical systems
subject to unilateral constraints; in particular, a control law has been proposed for one degree-of-freedom mechanical systems, as well as for certain multidegree-of-freedom mechanical systems. The collisions of robotic manipulators are modeled and simulated in [32]. A more extensive bibliography about modeling and control of impact can be found in the recent reference [1], where 620 references can be found about the subject.

All the contributions about impact analysis and control (including this paper) suffer from several practical limitations that can limit their effectiveness. These limitations arise from the difficulty to a priori establish the approximation level to be used for the study of the various phenomena involved in the impact and to understand the manner in which the various parameters (masses, geometric dimensions, stiffness, duration of the impact) affect the actual collision. Another practical limitation is that powerful tools (such as strain gages, high-speed photography, fast velocity, and motion transducers) are needed, e.g., for precise data storing. The limitations include that the mathematics to be used for a rigorous model of the impact may be excessively complex, whereas approximate models may omit one or more important characteristics of the impact phenomenon. In spite of these limitations, we hope that the analysis that will be carried out in this paper will provide, in view of its generality and simplicity, a useful guide to the design of models and controllers of impacting systems.

The outline of the paper is as follows. Section II is devoted to the writing of the equations of motion of mechanical systems subject to inequality constraints. After some preliminary notation and results, the equations of motion in case of nonsmooth impacts are briefly introduced in Section II-A through the method of the Valentine variables, whereas those in case of smooth impacts are introduced in Section II-B through the method of the penalty functions; Section II-C gives sufficient conditions for the approximated path of motion (obtained by means of the method of the penalty functions) to tend to the actual path of motion (obtained by means of the method of the Valentine variables) as the penalization tends to infinity. The control of the equations of motion thus obtained is analyzed in Section III in case of nonsmooth impacts and in Section IV in case of smooth impacts, with reference to the simple regulation problem. Section V reports some experimental results obtained by controlling the impact of a single-link robot arm against a rigid obstacle located in the working space.

## II. EQUATIONS OF MOTION OF MECHANICAL SYSTEMS SUBJECT TO SMOOTH AND NONSMOOTH IMPACTS

Consider a finite-dimensional mechanical system, or a finite-dimensional approximation of a continuous mechanical system. Let $q(t) \in \mathbb{R}^n$, with $n \in \mathbb{Z}$, $n \ge 1$, be the vector of the generalized coordinates $q_i(t)$, $i = 1, 2, \dots, n$, which are assumed to represent uniquely the configuration at time $t \in \mathbb{R}$ of the mechanical system (in the whole); in the remainder of the paper, the generalized coordinates will be assumed to be continuous functions of time $t \in \mathbb{R}$, which, in addition, are piecewise smooth. The times $t_c \in \mathbb{R}$ at which the vector

function $q(t)$ is not differentiable (i.e., those times $t_c$ such that $\lim_{t \to t_c^-} \dot{q}(t) \neq \lim_{t \to t_c^+} \dot{q}(t)$, which correspond to corner points) will be referred to as the times at which the impacts occur (briefly, the *impact times*). The assumption that $q(t)$ is smooth between two adjacent impact times is not necessary, but it is supposed for the sake of simplicity; such a vector function (only) needs to be twice differentiable between two adjacent impact times. To simplify the notation, the symbols $\alpha(t_c^-)$ and $\alpha(t_c^+)$ will be used to denote, respectively, the values taken by the limits $\lim_{t \to t_c^-} \alpha(t)$ and $\lim_{t \to t_c^+} \alpha(t)$, when they are definite, for any function $\alpha(t)$.

Assume that some constraints are imposed on the generalized coordinates by the physical nature of the mechanical system under consideration; in particular, assume that (at each time $t \in \mathbb{R}$) the vector $q(t)$ must belong to the following admissible region of $\mathbb{R}^n$:

$$ \mathcal{A} := \{q \in \mathbb{R}^n: f_i(q) \le 0, \quad i = 1, 2, \dots, m\} \quad (1) $$

where $f_i(\cdot): \mathbb{R}^n \to \mathbb{R}, i = 1, 2, \dots, m$, are $m \in \mathbb{Z}, m \ge 1$ smooth functions of $q \in \mathbb{R}^n$. Also in this case, the assumption that the $m$ functions $f_i(q)$, $i = 1, 2, \dots, m$, are smooth is considered for the sake of simplicity; such functions (only) need to be twice differentiable with respect to all the variables at argument. Assume that the admissible region $\mathcal{A}$ is nonempty and that, in addition, there exists at least an interior point $q_0 \in \mathcal{A}$ such that all the constraints are strictly satisfied, i.e., $f_i(q_0) < 0, i = 1, 2, \dots, m$.

The vector inequality $f(q) \le 0$, with $f := [f_1 f_2 \dots f_m]^T$, will be used to denote the entire set of the $m$ inequalities $f_i(q) \le 0, i = 1, 2, \dots, m$. The manifolds of $\mathbb{R}^n$ identified by $f_i(q) = 0, i = 1, 2, \dots, m$, represent the condition of contact among parts of the mechanical system with themselves or with the external environment, which is assumed to be infinitely rigid and massive.

Let $J_i(q) := \partial f_i(q)/\partial q$ be the gradient vector of $f_i(q)$, $i = 1, 2, \dots, m$; note that if $\alpha(q)$ is a scalar function of $q$, then $\partial\alpha(q)/\partial q$ denotes the column vector having as entries the partial derivatives $(\partial\alpha(q)/\partial q_i)$, $i = 1, 2, \dots, n$. Let $q \in \mathcal{A}$ be arbitrary and let $\{j, h, \dots, k\}$ denote the (possibly, empty) set of the indexes such that $f_j(q) = 0, f_h(q) = 0, \dots, f_k(q) = 0$, and $f_i(q) \ne 0$ for all $i \notin \{j, h, \dots, k\}$ for such a $q$; in the remainder of the paper, assume that the gradient vectors $J_j(q), J_h(q), \dots, J_k(q)$ are linearly independent for such a $q$.

Let $T(q(t), \dot{q}(t))$ and $U_t(q(t))$ be the kinetic energy and the total potential energy of the mechanical system at time $t \in \mathbb{R}$, respectively. Assume that the kinetic energy $T(q(t), \dot{q}(t))$ can be expressed as follows:

$$ T(q(t), \dot{q}(t)) = \frac{1}{2} \dot{q}^T(t) B(q(t)) \dot{q}(t) \quad (2) $$

where $B(q)$ is the generalized inertia matrix, which is positive definite for all $q \in \mathbb{R}^n$. For the sake of simplicity, the entries of $B(q)$ are assumed to be smooth (although it is necessary only that they are twice differentiable) with respect to all the variables at argument.

Consider the following assumption.
**Assumption 1:** There exist two real numbers $b, \bar{b}$, with $0 < b \le \bar{b} < +\infty$, such that the following inequalities hold in the whole (i.e., for all $q \in \mathbb{R}^n$):

$$ \frac{1}{2} b \|\dot{q}\|^2 \le \frac{1}{2} \dot{q}^T B(q) \dot{q} \le \frac{1}{2} \bar{b} \|\dot{q}\|^2, \quad \forall \dot{q} \in \mathbb{R}^n \quad (3) $$

where $||\cdot||$ is the Euclidean norm of the vector $\cdot$ at argument.

Assumption 1 implies that the *limit superior* upon $\mathbb{R}^n$ of the greatest eigenvalue of $B(q)$ is finite and that the *limit inferior* upon $\mathbb{R}^n$ of the smallest eigenvalue of $B(q)$ is greater than zero.

In addition, assume throughout the paper that the total potential energy $U_t(q(t))$ can be decomposed into two parts:

$$ U_t(q(t), t) = U(q(t)) - q^T(t)E u(t) \quad (4) $$

where $U(q(t))$ is the potential energy due to the conservative forces (e.g., gravitational and/or elastic forces, if any), $u(t) \in \mathbb{R}^p$ is the vector of the control (generalized) forces, $p \in \mathbb{Z}$, with $n \ge p \ge 1$, $E \in \mathbb{R}^{n \times p}$ is a full column-rank matrix constituted by $p$ linearly independent columns of the $n$-dimensional identity matrix, suitably chosen in order to identify the actuated generalized coordinates, and $-q^T(t)E u(t)$ is the potential energy due to the control forces $u(t)$.

The equations describing the motion of the mechanical system can be carried out by considering small variations of the whole motion of the system between two times $t_1$ and $t_2$, with $t_2 > t_1$. As all the forces acting on the mechanical system have been assumed to derive from a scalar energy, the *Hamilton principle* [33] can be used to obtain the equations of motion in the closed interval $[t_1, t_2]$; the actual path of motion of the mechanical system under consideration from time $t_1$ to time $t_2$ is such that the following *action integral* has a stationary value:

$$ A = \int_{t_1}^{t_2} L_t(q(t), \dot{q}(t), t) dt \quad (5) $$

where $L_t(q(t), \dot{q}(t), t) := T(q(t), \dot{q}(t)) - U_t(q(t), t)$ subject to the following inequality constraints to be satisfied for all $t \in [t_1, t_2]$:

$$ f_i(q(t)) \le 0, \quad i = 1, 2, \dots, m. \quad (6) $$

### A. Equations of Motion of Nonsmooth Impacts Through the Valentine Variables

The *Valentine variables* are introduced so that the inequality constraints (6) are transformed into the following equality constraints:

$$ f_i(q(t)) + \gamma_i^2(t) = 0, \quad i = 1, 2, \dots, m \quad (7) $$

with $\gamma_i(t), i = 1, 2, \dots, m$, being some auxiliary nonnegative real-valued variables, which were introduced by Valentine for the solution of the problem of Lagrange with differential inequalities [34]. Since the *Valentine variables* $\gamma_i(t), i = 1, 2, \dots, m$, are taken real, the equality constraints (7) are completely equivalent to the inequality constraints (6). Since $\gamma_i^2(t) = -f_i(q(t))$, all the *Valentine variables* are continuous functions of time, which are piecewise smooth, with corner points that can occur in correspondence of the impact times.

Taking the derivative with respect to time of both sides of (7), the following differential constraints are obtained:

$$ J_i^T(q(t))\dot{q}(t) + 2\gamma_i(t)\dot{\gamma}_i(t) = 0, \quad i = 1, 2, \dots, m. \quad (8) $$

Starting from the initial conditions $q(0) \in \mathcal{A}, \gamma_i(0) = \sqrt{-f_i(q(0))}, i = 1, 2, \dots, m$, the differential constraints (8) are completely equivalent to the point constraints (7).

Then, the actual path of motion can be found by looking for the stationary value of the action integral (5), under the differential constraints (8), where the *Valentine variables* $\gamma_i(t), i = 1, 2, \dots, m$, are to be considered as nonnegative real-valued continuous functions of time, which are piecewise smooth between two adjacent impact times. In particular, introducing the *Lagrange multipliers* [35] $\lambda_i(t), i = 1, 2, \dots, m$, the actual path of motion can be determined by looking for the stationary value of the following unconstrained functional:

$$ \hat{A} = \int_{t_1}^{t_2} \hat{L} dt \quad (9) $$

where $\hat{L} := L_t(q(t), \dot{q}(t), t) + \sum_{i=1}^m \lambda_i(t)(J_i^T(q(t))\dot{q}(t) + 2\gamma_i(t)\dot{\gamma}_i(t))$.

By well-known results from functional analysis [36], the stationary value of functional (9) corresponds to the continuous piecewise smooth path of motion that is the solution (in each open interval of time that does not contain impact times) to the following $n+2m$ Euler-Lagrange equations (dependence upon time $t$ is omitted):

$$ \frac{d}{dt} \frac{\partial L}{\partial \dot{q}} - \frac{\partial L}{\partial q} + \sum_{i=1}^{m} \dot{\lambda}_i J_i(q) = E u, \quad (10a) $$

$$ 2\gamma_i \dot{\lambda}_i = 0, \quad i = 1, 2, \dots, m, \quad (10b) $$

$$ J_i^T(q)\dot{q} + 2\gamma_i\dot{\gamma}_i = 0, \quad i = 1, 2, \dots, m, \quad (10c) $$

where $L(q(t), \dot{q}(t)) := T(q(t), \dot{q}(t)) - U(q(t))$, and the impacts can occur only at the times $t_c \in \mathbb{R}$ where the following *Erdmann-Weierstrass corner conditions*, which are necessary at corner points, are satisfied:

$$ \frac{1}{2} \dot{q}^T(t_c^-) B(q(t_c^-)) \dot{q}(t_c^-) = \frac{1}{2} \dot{q}^T(t_c^+) B(q(t_c^+)) \dot{q}(t_c^+) \quad (11a) $$

$$ B(q(t_c^-)) \dot{q}(t_c^-) + \sum_{i=1}^{m} \lambda_i(t_c^-) J_i(q(t_c^-)) = B(q(t_c^-)) \dot{q}(t_c^+) + \sum_{i=1}^{m} \lambda_i(t_c^+) J_i(q(t_c^-)) \quad (11b) $$

$$ 2\gamma_i(t_c^-) \dot{\lambda}_i(t_c^-) = 2\gamma_i(t_c^-) \dot{\lambda}_i(t_c^+), \quad i = 1, 2, \dots, m. \quad (11c) $$

It is pointed out that the *Euler-Lagrange equations* (10) can also be obtained directly by balancing all the forces acting on the mechanical system, taking into account the impulsive forces due to the inequality constraints by means of some complementarity conditions such as (10b).
## B. Equations of Motion of Smooth Impacts
Through the Penalty Functions

This method is based on the observation that if the mechanical system under consideration is allowed to violate constraints (6) and, in correspondence to such a violation, the action functional (5) is strongly penalized, whereas it is not when constraints (6) are satisfied, then the solution of the problem thus modified is forced to tend to the actual path of motion of the constrained mechanical system, as the penalization tends to infinity. It is stressed that, in this case, the generalized coordinates $q(t)$ are no longer obliged to belong to the admissible region $A$.

A possible choice of the penalty function, which takes the zero value when all the constraints (6) are satisfied and high positive values when some of the constraints (6) are violated, is the following:

$$U_f(q(t)) = \sum_{i=1}^{m} U_{fi}(q(t)) \quad (12)$$

where

$$U_{fi}(q) := \begin{cases} 0, & \text{if } f_i(q) \le 0 \\ \frac{1}{2} k_i f_i^2(q), & \text{if } f_i(q) > 0 \end{cases}$$

and $k_i$ is a sufficiently high positive number, having a suitable physical dimension so that $\frac{1}{2} k_i f_i^2(q)$ has the dimension of an energy. As one can see when the $i$th constraint is violated, $U_{fi}(q(t))$ plays the role of the elastic energy due to the contact relative to the $i$th constraint, as though a linear elastic spring were located at the point of contact, with $k_i$ playing the role of the corresponding elastic constant. It is stressed that many other penalty functions can be chosen in place of (12), each one yielding a different approximate model.

Since, owing to this penalization, the total potential energy of the mechanical system is $U_t + U_f$, one can consider the action integral (5) with $L_t$ substituted by $\hat{L}_t := L_t - U_f$; assume that the path of motion corresponding to the unconstrained stationary value of the action integral (5) with $L_t$ substituted by $\hat{L}_t$ tends to the actual path of motion of the constrained mechanical system, as all the constants $k_i, i = 1, 2, \dots, m$, go to infinity (some sufficient conditions for this property to be satisfied will be given later in this section). Then, once sufficiently high values of the constants $k_i, i = 1, 2, \dots, m$, have been chosen, an approximation of the actual path of motion of the constrained mechanical system can be computed by looking for the unconstrained stationary value of (5) with $L_t$ substituted by $\hat{L}_t$. The advantage of this method, as compared with the method of the Valentine variables, is that in this case no Lagrange multipliers are to be introduced, as well as no other auxiliary variables. The main drawback of the method of the penalty functions is the choice of the values of constants $k_i, i = 1, 2, \dots, m$, which is often carried out through repeated trials, by taking greater values of $k_i$ up to the situation in which the corresponding path of motion does not change significantly.

It should be pointed out that the addition of the penalty function will cause the conditioning of the problem of computing the stationary value of the action integral to be very

large (of the order of coefficients $k_i, i = 1, 2, \dots, m$). Thus, the condition number of the problem will be proportional to coefficients $k_i, i = 1, 2, \dots, m$; as one increases such coefficients in order to obtain more accurate solutions to the original constrained problem, the rate of convergence of any algorithm used to solve in a numeric way the equations of motions thus obtained becomes extremely slow, and any control law designed on the basis of these equations of motion may strongly depend on this very high penalization.

The stationary value of (5) with $L_t$ substituted by $\hat{L}_t$ corresponds to the path of motion that is solution of the following vector Euler-Lagrange equation:

$$\frac{d}{dt} \frac{\partial L}{\partial \dot{q}} - \frac{\partial L}{\partial q} + \sum_{i=1}^{m} \eta_i J_i(q) = Eu \quad (13)$$

where $L = T - U$ and

$$\eta_i := \begin{cases} 0, & \text{if } f_i(q) \le 0, \\ k_i f_i(q), & \text{if } f_i(q) > 0. \end{cases} \quad (14)$$

The Erdmann-Weierstrass corner conditions are

$$\frac{1}{2} \dot{q}^T(t_c^-) B(q(t_c)) \dot{q}(t_c^-) = \frac{1}{2} \dot{q}^T(t_c^+) B(q(t_c)) \dot{q}(t_c^+) \quad (15a)$$

$$B(q(t_c)) \dot{q}(t_c^-) = B(q(t_c)) \dot{q}(t_c^+) \quad (15b)$$

The Erdmann-Weierstrass corner condition (15b) can be rewritten as $B(q(t_c))(\dot{q}(t_c^+) - \dot{q}(t_c^-)) = 0$, which (by the nonsingularity of $B(q(t_c))$) implies that $\dot{q}(t_c^+) = \dot{q}(t_c^-)$, namely, when the method of the penalty functions is used, no jump is possible of the generalized velocities. The impacts thus approximated are referred to as smooth, although $\dot{q}(t)$ is required to be only continuous by (15b).

By comparing (13) with (10a), $\eta_i(t)$ can be thought of as a continuous piecewise smooth approximation of the piecewise smooth $\lambda_i(t)$.

## C. Properties of the Penalty Method

The aim of this section is to give sufficient conditions under which the path of motion obtained by means of the method of the penalty functions tends to the actual path of motion of the constrained mechanical system, as the penalization tends to infinity.

For the sake of simplicity, consider the case of a single inequality constraint $f(q) \le 0$ (i.e., $m = 1$). Let $\{k_h\}, h = 1, 2, \dots$, be a sequence of real numbers such that $\lim_{h \to +\infty} k_h = +\infty$ and, for each $h \ge 1$, such that $k_h \ge 0, k_{h+1} > k_h$. Define the following functionals:

$$A(q(\cdot)) := - \int_{t_1}^{t_2} L_t(q(t), \dot{q}(t), t) dt, \quad (16a)$$

$$P_f(q) := \begin{cases} 0, & \text{if } f(q) \le 0, \\ \frac{1}{2} f^2(q), & \text{if } f(q) > 0, \end{cases} \quad (16b)$$

$$\tilde{A}(q(\cdot); k) := A(q(\cdot)) + k \int_{t_1}^{t_2} P_f(q(t)) dt, \quad (16c)$$

with $L_t(q(t), \dot{q}(t), t)$ being defined as in (5). Our aim is to find an approximation of the stationary point of $A$, subject to the inequality constraint, by looking for the stationary point of the penalized functional $\tilde{A}$. For each $h \ge 1$, assume
that the stationary point $q_h^*(\cdot)$ of $\tilde{A}(q(\cdot); k_h)$ corresponds to a minimum, as well as the stationary point $q^*(\cdot)$ of $A(q(\cdot))$ subject to the inequality constraint. A sufficient condition for $\tilde{A}(q(\cdot); k_h)$ to have a minimum is that $\tilde{A}(q(\cdot); k_h) \to +\infty$ as $\|q(\cdot)\|_{\mathcal{F}} \to +\infty$, with $\|q(\cdot)\|_{\mathcal{F}}$ being a suitable norm of function $q(\cdot)$ over $[t_1, t_2]$; similar sufficient condition can be given for $A(q(\cdot))$.

The following theorems can be easily stated and proven.

**Theorem 1:** If, for each $h \ge 1$, $q_h^*(\cdot)$ is a minimum of $\tilde{A}(q(\cdot); k_h)$, then the following inequalities are satisfied for all $h \ge 1$:

$$ \tilde{A}(q_h^*(\cdot); k_h) \le \tilde{A}(q_{h+1}^*(\cdot); k_{h+1}) \quad (17a) $$

$$ \int_{t_1}^{t_2} P_f(q_h^*(t)) dt \ge \int_{t_1}^{t_2} P_f(q_{h+1}^*(t)) dt \quad (17b) $$

$$ A(q_h^*(\cdot)) \le A(q_{h+1}^*(\cdot)). \quad (17c) $$

*Proof:* The following relations are obtained from the definition of $\tilde{A}(q(\cdot); k)$, taking into account that $P_f(q) \ge 0, \forall q \in \mathbb{R}^n$, and that $k_{h+1} > k_h$

$$ \begin{aligned} \tilde{A}(q_h^*(\cdot); k_{h+1}) &= A(q_{h+1}^*(\cdot)) + k_{h+1} \int_{t_1}^{t_2} P_f(q_{h+1}^*(t)) dt \\ &\ge A(q_{h+1}^*(\cdot)) + k_h \int_{t_1}^{t_2} P_f(q_{h+1}^*(t)) dt. \end{aligned} \quad (18) $$

Since $q_h^*(\cdot)$ is a minimum of $\tilde{A}(q(\cdot); k_h)$, one has

$$ \begin{aligned} \tilde{A}(q_h^*(\cdot); k_h) &\le \tilde{A}(q_{h+1}^*(\cdot); k_h) \\ &\le A(q_{h+1}^*(\cdot)) + k_h \int_{t_1}^{t_2} P_f(q_{h+1}^*(t)) dt. \end{aligned} \quad (19) $$

Inequality (17a) is obtained by (18) and (19).

The following inequalities can be obtained taking into account that $q_h^*(\cdot)$ is a minimum of $\tilde{A}(q(\cdot); k_h)$ and that $q_{h+1}^*(\cdot)$ is a minimum of $\tilde{A}(q(\cdot); k_{h+1})$:

$$ \begin{aligned} A(q_h^*(\cdot)) + k_h \int_{t_1}^{t_2} P_f(q_h^*(t)) dt &\le A(q_{h+1}^*(\cdot)) + k_h \int_{t_1}^{t_2} P_f(q_{h+1}^*(t)) dt \end{aligned} \quad (20a) $$

$$ \begin{aligned} A(q_{h+1}^*(\cdot)) + k_{h+1} \int_{t_1}^{t_2} P_f(q_{h+1}^*(t)) dt &\le A(q_h^*(\cdot)) + k_{h+1} \int_{t_1}^{t_2} P_f(q_h^*(t)) dt. \end{aligned} \quad (20b) $$

By adding (20a) and (20b), we obtain

$$ \begin{aligned} & (k_{h+1} - k_h) \int_{t_1}^{t_2} P_f(q_{h+1}^*(t)) dt \\ &\le (k_{h+1} - k_h) \int_{t_1}^{t_2} P_f(q_h^*(t)) dt. \end{aligned} \quad (21) $$

Inequality (17b) is obtained from (21) because $k_{h+1} - k_h > 0$.

Since

$$ \begin{aligned} A(q_{h+1}^*(\cdot)) + k_h \int_{t_1}^{t_2} P_f(q_{h+1}^*(t)) dt &\ge A(q_h^*(\cdot)) + k_h \int_{t_1}^{t_2} P_f(q_h^*(t)) dt \end{aligned} $$

from (17b) we have

$$ \begin{aligned} A(q_{h+1}^*(\cdot)) - A(q_h^*(\cdot)) &\ge k_h \left( \int_{t_1}^{t_2} P_f(q_h^*(t)) dt - \int_{t_1}^{t_2} P_f(q_{h+1}^*(t)) dt \right) \\ &\ge 0 \end{aligned} $$

thus proving (17c).

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

**Theorem 2:**

If $q^*(\cdot)$ is a minimum of $A(q(\cdot))$ subject to the inequality constraint $f(q(t)) \le 0, \forall t \in [t_1, t_2]$, and, for each $h \ge 1$, $q_h^*(\cdot)$ is a minimum of $\tilde{A}(q(\cdot); k_h)$, then the following inequalities are satisfied for all $h \ge 1$:

$$ A(q^*(\cdot)) \ge \tilde{A}(q_h^*(\cdot); k_h) \ge A(q_h^*(\cdot)). $$

*Proof:* By the definition of $q^*(\cdot)$, we have $f(q^*(t)) \le 0, \forall t \in [t_1, t_2]$, whence $\int_{t_1}^{t_2} P_f(q^*(t)) dt = 0$. This yields the following relations:

$$ \begin{aligned} A(q^*(\cdot)) &= A(q^*(\cdot)) + k_h \int_{t_1}^{t_2} P_f(q^*(t)) dt \\ &= \tilde{A}(q^*(\cdot); k_h) \\ &\ge \tilde{A}(q_h^*(\cdot); k_h) \\ &\ge A(q_h^*(\cdot)) + k_h \int_{t_1}^{t_2} P_f(q_h^*(t)) dt \\ &\ge A(q_h^*(\cdot)) \end{aligned} $$

which complete the proof.

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

The following theorem states the desired sufficient conditions of convergence.

**Theorem 3:** Suppose, for each $h \ge 1$, that $q_h^*(\cdot)$ is a minimum of $\tilde{A}(q(\cdot); k_h)$ and that $q^*(\cdot)$ is a minimum of $A(q(\cdot))$ subject to the inequality constraint $f(q(t)) \le 0, \forall t \in [t_1, t_2]$. If $A(q(\cdot))$ depends continuously on $q(\cdot)$, then any limit point of the sequence $\{q_h^*(\cdot)\}, h = 1, 2, \dots$, is a minimum of $A(q(\cdot))$ subject to the inequality constraint $f(q(t)) \le 0, \forall t \in [t_1, t_2]$.

*Proof:* Assume that the sequence $\{q_h^*(\cdot)\}, h = 1, 2, \dots$, is convergent with limit $\bar{q}(\cdot)$ as $h \to +\infty$, otherwise consider a convergent subsequence of $\{q_h^*(\cdot)\}, h = 1, 2, \dots$. By the assumed continuity of $A(q(\cdot))$, one has

$$ \lim_{h \to +\infty} A(q_h^*(\cdot)) = A(\bar{q}(\cdot)). \quad (22) $$

By Theorems 1 and 2, the sequence $\{\tilde{A}(q_h^*(\cdot); k_h)\}, h = 1, 2, \dots$, is nondecreasing and bounded above by $A(q^*(\cdot))$, whence

$$ A^* := \lim_{h \to +\infty} \tilde{A}(q_h^*(\cdot); k_h) \le A(q^*(\cdot)). \quad (23) $$
Subtracting (22) from (23), we obtain

$$ \lim_{h \to +\infty} k_h \int_{t_1}^{t_2} P_f(q_h^*(t)) dt = A^* - A(\bar{q}(\cdot)). \quad (24) $$

Since $\int_{t_1}^{t_2} P_f(q_h^*(t)) dt \ge 0, \forall h \ge 1$, and $\lim_{h \to +\infty} k_h = +\infty$, (24) yields $\lim_{h \to +\infty} \int_{t_1}^{t_2} P_f(q_h^*(t)) dt = 0$. Whence, by the continuity of $P_f(q)$ and by the property $P_f(q) \ge 0, \forall q \in \mathbb{R}^n$, we have $P_f(\bar{q}(t)) = 0, \forall t \in [t_1, t_2]$, which implies that the inequality constraint $f(q(t)) \le 0$ is satisfied over $[t_1, t_2]$ by $q(t) = \bar{q}(t)$.

Since, by Theorem 2

$$ A(\bar{q}(\cdot)) = \lim_{h \to +\infty} A(q_h^*(\cdot)) \\ \le A(q^*(\cdot)) $$

with $q^*(\cdot)$ being the minimum of $A(q(\cdot))$ subject to the inequality constraint $f(q(t)) \le 0, \forall t \in [t_1, t_2]$, the limit point $\bar{q}(\cdot)$ is also a minimum of $A(q(\cdot))$ subject to the inequality constraint $f(q(t)) \le 0, \forall t \in [t_1, t_2]$. $\square$

**Remark 1:** Theorems 1–3 can be easily emended to deal with the case when, for each $h \ge 1$, $q_h^*(\cdot)$ is a maximum of $\tilde{A}(q(\cdot); k_h)$ and $q^*(\cdot)$ is a maximum of $A(q(\cdot))$ subject to the inequality constraint $f(q(t)) \le 0, \forall t \in [t_1, t_2]$. $\square$

### III. CONTROL OF NONSMOOTH IMPACTS

Assume that the generalized coordinates that are actuated by $u(t)$ are the only coordinates available for feedback (they are the so-called natural outputs), as well as their time derivatives, $y(t) := E^T q(t)$ and $\dot{y}(t) := E^T \dot{q}(t)$. First, consider the following “derivative” control law:

$$ u(t) = -K_v \dot{y}(t) \quad (25) $$

where $K_v$ is a positive definite square matrix of dimensions $p \times p$. System (10), under the control law (25), becomes

$$ \frac{d}{dt} \frac{\partial L}{\partial \dot{q}} - \frac{\partial L}{\partial q} + \sum_{i=1}^{m} \dot{\lambda}_i J_i(q) = -E K_v E^T \dot{q}, \quad (26a) $$

$$ 2\gamma_i \dot{\lambda}_i = 0, \quad i = 1, 2, \dots, m, \quad (26b) $$

$$ J_i^T(q)\dot{q} + 2\gamma_i\dot{\gamma}_i = 0, \quad i = 1, 2, \dots, m. \quad (26c) $$

The Erdmann–Weierstrass corner conditions (11) remain unchanged.

An integral curve for system (26), under conditions (11), with initial condition $(q(0), \dot{q}(0^{-})) \in A \times \mathbb{R}^n, \gamma_i(0) = \sqrt{-f_i(q(0))}, i = 1, 2, \dots, m$, is a mapping $(s_1(\cdot; p^1, p^2), s_2(\cdot; p^1, p^2))$: $[0, T) \to A \times \mathbb{R}^n$ (with $T \in \mathbb{R}, T \ge 0$, and, possibly, $T = +\infty$) such that $(s_1(t; p^1, p^2), s_2(t; p^1, p^2))$ is differentiable with respect to time almost everywhere $(s_1(0^{-}; p^1, p^2), s_2(0^{-}; p^1, p^2)) = (p^1, p^2)$, the impact times $t_{c,i} \in [0, T), t_{c,i+1} > t_{c,i}, i = 1, 2, \dots$, at which $(s_1(t; p^1, p^2), s_2(t; p^1, p^2))$ is not differentiable with respect to time form a set $\mathcal{T}_c$ of null measure, the differential equations (26) are satisfied by $q(t), \dot{q}(t) = (s_1(t; p^1, p^2), s_2(t; p^1, p^2))$ for all $t \in [0, T) - \mathcal{T}_c$, and the Erdmann–Weierstrass

corner conditions (11) are satisfied by $(q(t_c), \dot{q}(t_c)) = (s_1(t_c; p^1, p^2), s_2(t_c; p^1, p^2))$ for all the impact times $t_c \in \mathcal{T}_c$. It is pointed out that, in general, the impact times depend on the initial condition, and the notation $t_{c,i}$ is only a shortening of $t_{c,i}(p^1, p^2)$.

Let $C$ be the greatest subset of $A \times \mathbb{R}^n$ such that for each initial condition $(q(0), \dot{q}(0^{-})) \in C, \lambda_i(0^{-}) \in \mathbb{R}, \gamma_i(0) = \sqrt{-f_i(q(0))}, i = 1, 2, \dots, m$, system (26), under conditions (11), has a unique solution upon the interval $[0, +\infty)$. As system (26), under conditions (11), is time-invariant, set $C$ is constituted by half-trajectories of system (26), under conditions (11), whence it is invariant.

**Assumption 2:** Let $C$ be nonempty and, in addition, for each initial condition $(q(0), \dot{q}(0^{-})) \in C, \lambda_i(0^{-}) \in \mathbb{R}, \gamma_i(0) = \sqrt{-f_i(q(0))}, i = 1, 2, \dots, m$, the solution of system (26), under conditions (11), is, in the variables $(q, \dot{q})$, a continuous function of the initial conditions $(q(0), \dot{q}(0^{-})) \in C$ at all times different from the impact times.

The following theorem gives sufficient conditions for the solution of system (26), under conditions (11), to exist, to be unique, and to be, in the variables $(q, \dot{q})$, a continuous function of the initial conditions at all times different from the impact times.

**Theorem 4:** Suppose that $\mathcal{B}$ is the greatest subset of $A \times \mathbb{R}^n$ such that (26), with $\dot{\lambda}_i, i = 1, 2, \dots, m$, being constant, has a unique solution $(q(t), \dot{q}(t)), t \in [0, +\infty)$, for all the initial conditions $(q(0), \dot{q}(0^{-})) \in \mathcal{B}$; in addition, suppose that such a solution is a continuous function of the initial conditions $(q(0), \dot{q}(0^{-})) \in \mathcal{B}$, at all times $t \in [0, +\infty)$. Denoting by $\partial A$ the boundary of $A$, suppose that the Erdmann–Weierstrass corner conditions (11) admit a unique solution $\dot{q}(t_c^+) = \alpha(q(t_c), \dot{q}(t_c^-))$ for all $(q(t_c), \dot{q}(t_c^-)) \in \partial A \times \mathbb{R}^n$, with $\alpha(\cdot, \cdot)$ being a suitable function. Define set $C$ as the union of $((A - \partial A) \times \mathbb{R}^n) \cap \mathcal{B}$ with the set constituted by all the points $(p^1, p^2) \in (\partial A \times \mathbb{R}^n) \cap \mathcal{B}$ such that $p^2 = \alpha(p^1, p^2)$. Let $(p^1, p^2) \in C$; if the impact times $t_{c,i}(p^1, p^2) \in [0, +\infty), t_{c,i}(p^1, p^2) < t_{c,i+1}(p^1, p^2), i = 1, 2, \dots$, form a set of null measure and $\lim_{i \to +\infty} t_{c,i}(p^1, p^2) = +\infty$, then system (26), under the Erdmann–Weierstrass corner conditions (11), has a unique solution $(q(t), \dot{q}(t)) = (s_1(t; p^1, p^2), s_2(t; p^1, p^2))$ over $[0, +\infty)$ from the initial condition $(q(0), \dot{q}(0^{-})) = (p^1, p^2)$. In addition, if the impact times $t_{c,i}(p^1, p^2) \in [0, +\infty), i = 1, 2, \dots$, are continuous functions of $(p^1, p^2)$ in $C$ and function $\alpha(\cdot, \cdot)$ is continuous in $\partial A \times \mathbb{R}^n$, then $(s_1(t; p^1, p^2), s_2(t; p^1, p^2))$ are continuous functions of $(p^1, p^2)$ in $C$ for all $t \in [0, +\infty), t \ne t_{c,i}, i = 1, 2, \dots$.

*Proof:* By assumption, system (26) has a unique solution $(\zeta^1(t), \dot{\zeta}^1(t)) \in \mathcal{B}$ starting from the initial condition $(q(0), \dot{q}(0^{-})) = (p^1, p^2)$ up to the first impact time $t_c(p^1, p^2)$; by assumption, such a solution is a continuous function of $(p^1, p^2)$. Then, $(q^1(t), \dot{q}^1(t)) = (\zeta^1(t), \dot{\zeta}^1(t))$ is the unique solution of system (26), under the Erdmann–Weierstrass corner conditions (11), over the interval $[0, t_{c,1})$, with $t_{c,1} = t_c(p^1, p^2)$ and, for each $t \in [0, t_{c,1})$, it is a continuous function of $(p^1, p^2)$.

By construction, $(\zeta^1(t_{c,1}), \dot{\zeta}^1(t_{c,1})) \in (\partial A \times \mathbb{R}^n) \cap \mathcal{B}$ and $(\zeta^1(t_{c,1}), \alpha(\zeta^1(t_{c,1}), \dot{\zeta}^1(t_{c,1}))) \in C$. By assumption,
system (26) has a unique solution $(\zeta^2(t), \dot{\zeta}^2(t)) \in \mathcal{B}$ starting from the initial condition $(q(t_{c,1}), \dot{q}(t_{c,1}^-)) = (\zeta^1(t_{c,1}), \alpha(\zeta^1(t_{c,1}), \dot{\zeta}^1(t_{c,1}^-)))$ up to the first impact time $t_c((\zeta^1(t_{c,1}), \alpha(\zeta^1(t_{c,1}), \dot{\zeta}^1(t_{c,1}^-)))$; by assumption, such a solution is a continuous function of $((\zeta^1(t_{c,1}), \alpha(\zeta^1(t_{c,1}), \dot{\zeta}^1(t_{c,1}^-)))$. Then, by the uniqueness of the solution of the Erdmann-Weierstrass corner conditions (11), the function defined by

$$ (q^2(t), \dot{q}^2(t)) = \begin{cases} (q^1(t), \dot{q}^1(t)), & t \in [0, t_{c,1}), \\ (\zeta^2(t), \dot{\zeta}^2(t)), & t \in [t_{c,1}, t_{c,2}), \end{cases} $$

is the unique solution of system (26), under the Erdmann-Weierstrass corner conditions (11), over the interval $[0, t_{c,2})$, with $t_{c,2} = t_c((\zeta^1(t_{c,1}), \alpha(\zeta^1(t_{c,1}), \dot{\zeta}^1(t_{c,1}^-)))$. In addition, by the continuity of $t_{c,1}$ on $(p^1, p^2)$ and by the continuity of function $\alpha(\cdot, \cdot)$, $(q^2(t), \dot{q}^2(t))$ is a continuous function of $(p^1, p^2)$ for all $t \in [0, t_{c,2})$, $t \neq t_{c,1}$ (note that if $t \neq t_{c,1}$, then there exists a neighborhood of $t$ that does not contain $t_{c,1}$).

Proceeding by induction, by construction, one has that $(\zeta^h(t_{c,h}), \dot{\zeta}^h(t_{c,h}^-)) \in (\partial A \times \mathbb{R}^n) \cap \mathcal{B}$ and $(\zeta^h(t_{c,h}), \alpha(\zeta^h(t_{c,h}), \dot{\zeta}^h(t_{c,h}^-))) \in \mathcal{C}$; by assumption, system (26) has a unique solution $(\zeta^{h+1}(t), \dot{\zeta}^{h+1}(t)) \in \mathcal{B}$ starting from the initial condition $(q(t_{c,h}), \dot{q}(t_{c,h}^-)) = (\zeta^h(t_{c,h}), \alpha(\zeta^h(t_{c,h}), \dot{\zeta}^1(t_{c,h}^-)))$ up to the first impact time $t_c((\zeta^h(t_{c,h}), \alpha(\zeta^h(t_{c,h}), \dot{\zeta}^h(t_{c,h}^-)))$; by assumption, such a solution is a continuous function of $((\zeta^h(t_{c,h}), \alpha(\zeta^h(t_{c,h}), \dot{\zeta}^h(t_{c,h}^-)))$. Then, by the uniqueness of the solution of the Erdmann-Weierstrass corner conditions (11), the function defined by

$$ (q^{h+1}(t), \dot{q}^{h+1}(t)) = \begin{cases} (q^h(t), \dot{q}^h(t)), & t \in [0, t_{c,h}), \\ (\zeta^{h+1}(t), \dot{\zeta}^{h+1}(t)), & t \in [t_{c,h}, t_{c,h+1}), \end{cases} $$

is the unique solution of system (26), under the Erdmann-Weierstrass corner conditions (11), over the interval $[0, t_{c,h+1})$, with $t_{c,h+1} = t_c((\zeta^h(t_{c,h}), \alpha(\zeta^h(t_{c,h}), \dot{\zeta}^h(t_{c,h}^-)))$. In addition, by the continuity of $t_{c,i}$, $i = 1, 2, \dots, h$, on $(p^1, p^2)$ and by the continuity of function $\alpha(\cdot, \cdot)$, $(q^{h+1}(t), \dot{q}^{h+1}(t))$ is a continuous function of $(p^1, p^2)$ for all $t \in [0, t_{c,h+1})$, $t \neq t_{c,i}$, $i = 1, 2, \dots, h$ (note that if $t \neq t_{c,i}$, then there exists a neighborhood of $t$ that does not contain $t_{c,i}$).

In this manner, by the assumption $\lim_{i \to +\infty} t_{c,i} = +\infty$, the solution can be extended to all of $[0, +\infty)$, also showing its continuity with respect to $(p^1, p^2)$ for all $t$ except for the impact times. $\square$

**Assumption 3:** The following algebraic equations in the unknown variables $q, \dot{\lambda}_i, \gamma_i, i = 1, 2, \dots, m$:

$$ \frac{\partial U}{\partial q} + \sum_{i=1}^{m} \dot{\lambda}_i J_i(q) = 0 \qquad (27a) $$

$$ \gamma_i \dot{\lambda}_i = 0, \quad i = 1, 2, \dots, m, \qquad (27b) $$

$$ \dot{\gamma}_i^2 + f_i(q) = 0, \quad i = 1, 2, \dots, m, \qquad (27c) $$

have a unique solution $q = q_R$, $\dot{\lambda}_i = \dot{\lambda}_{i,R}$, $\gamma_i = \gamma_{i,R} := \sqrt{-f_i(q_R)}$, $i = 1, 2, \dots, m$, with $q_R \in \mathcal{A}$ and

$\dot{\lambda}_{i,R}$ being nonnegative, $i = 1, 2, \dots, m$; assume, in addition, that $q_R$ is not an isolated point of $\mathcal{A}$ and that $(q_R, 0) \in \mathcal{C}$.

**Remark 2:** For arbitrarily fixed initial values of the Lagrange multipliers $\lambda_i(0^-)$, $i = 1, 2, \dots, m$, Assumption 3 implies that (26), under conditions (11), has a unique solution with $q, \dot{\lambda}_i$ and $\gamma_i, i = 1, 2, \dots, m$, being constant. Unfortunately, such solutions are not asymptotically stable, in the usual sense, since different initial conditions for the Lagrange multipliers $\lambda_i(0^-)$, $i = 1, 2, \dots, m$, yield different solutions that do not converge toward each other (there exists a manifold of fixed points). A possible solution to this problem could be to consider $\dot{\lambda}_i(t)$ as one of the state variables of system (26) instead of $\lambda_i(t)$. With this position, Assumption 3 implies that system (26), under conditions (11), has a unique equilibrium point $q = q_R$, $\dot{q} = 0$, $\dot{\lambda}_i = \dot{\lambda}_{i,R}$, $\gamma_i = \sqrt{-f_i(q_R)}$, $i = 1, 2, \dots, m$. Unfortunately, if $f_i(q_R) = 0$ for at least one index $i \in \{1, 2, \dots, m\}$, then such an equilibrium point cannot be stable, in the usual sense. As a matter of fact, for any real number $\varepsilon > 0$, consider the solution of (26), under conditions (11), starting from the initial condition $q(0^-) = q_R$, $\dot{q}(0^-) = \varepsilon v, \dot{\lambda}_i(0^-) = \dot{\lambda}_{i,R}, \gamma_i(0^-) = \sqrt{-f_i(q_R)}, i = 1, 2, \dots, m$, where $v \in \mathbb{R}^n, v \neq 0$, is a vector such that $J_i^T(q_R)v > 0$ for at least one index $i \in \{1, 2, \dots, m\}$ for which $f_i(q_R) = 0$ (such a $v$ can always be found). The choice of $v$ implies the existence of a jump $\dot{q}(0^+) - \dot{q}(0^-) \neq 0$ of the generalized velocities, which implies a jump $\lambda_i(0^+) - \lambda_i(0^-) \neq 0$ of at least the $i$th Lagrange multiplier; this implies an impulse (in the distribution sense) of $\dot{\lambda}_i(t)$, for any $\varepsilon > 0$, which prevents the stability of the equilibrium point. Certainly, the definition of stability could be properly modified to take into account the presence of impulses due to $\dot{\lambda}_i(t)$, but the following results will be given with respect to the usual definitions of stability $\square$ [37].

**Assumption 4:** Under Assumption 3, there exist two real numbers $a, \bar{a}$, with $0 < a \le \bar{a} < +\infty$, such that the following inequalities hold in the whole (i.e., for all $q \in \mathcal{A}$):

$$ \frac{1}{2} a \|q - q_R\|^2 \le U(q) - U(q_R) \le \frac{1}{2} \bar{a} \|q - q_R\|^2 \qquad (28) $$

where $\| \cdot \|$ is the Euclidean norm of the vector $\cdot$ at argument.

**Remark 3:** Assumption 4 implies that $q = q_R$ is a global minimum of $U(q)$ upon $\mathcal{A}$. If $q_R$ is an interior point of $\mathcal{A}$, then $(\partial U(q)/\partial q)|_{q=q_R} = 0$ and, by (27a) and by taking into account that the gradient vectors $J_j(q), J_h(q), \dots, J_k(q)$ have been assumed to be linearly independent for the indexes $j, h, \dots, k$ such that $f_j(q) = 0, f_h(q) = 0, \dots, f_k(q) = 0$, one has $\dot{\lambda}_{i,R} = 0, i = 1, 2, \dots, m$, whereas if $q_R$ belongs to the boundary $\partial \mathcal{A}$ of $\mathcal{A}$, then $(\partial U(q)/\partial q)|_{q=q_R}$ may be different from zero, as well as the corresponding $\dot{\lambda}_{i,R}, i = 1, 2, \dots, m$. It is stressed that Assumption 4 implies the positive definiteness of $U(q) - U(q_R)$ about $q_R$ upon $\mathcal{A}$, but not necessarily upon the whole $\mathbb{R}^n$; in particular, if the latter property were true, then $(\partial U(q)/\partial q)|_{q=q_R}$ would be equal to zero, with the corresponding $\dot{\lambda}_{i,R}, i = 1, 2, \dots, m$, being zero, as well. $\square$

It is important to stress that the Lyapunov stability of a particular class of mechanical systems subject to inequality constraints has been studied in [38] when the potential energy has a minimum in correspondence of the equilibrium
point; sufficient conditions are given in [38, Th. 1] for these mechanical systems to be Lyapunov stable in the simplified case of equilibrium points located upon the surface of contact, with the relative reaction force being zero, under the additional assumption that during an impact there is no loss of energy. The assumption of zero exchanged force at the point of contact (which is not made in this paper, as well as the one regarding the loss of energy at the impact times) renders difficult the application of this result in any task involving exchange of forces among parts of the mechanical systems or with the external environment. A technique for the verification of the sufficient conditions given in [38] for such mechanical systems to be Lyapunov stable has been proposed in [39].

**Theorem 5:** Under Assumptions 1–4, if the mechanical system under consideration is fully actuated (i.e., if $p = n$ and $E = I$), then the following properties are satisfied.

1) For each real number $\varepsilon > 0$, there exists a real number $\delta > 0$ such that for every initial condition $(q(0), \dot{q}(0^-)) \in C, \lambda_i(0^-) \in \mathbb{R}, \gamma_i(0) = \sqrt{-f_i(q(0))}, i = 1, 2, \dots, m$, satisfying $||q(0) - q_R|| < \delta$ and $||\dot{q}(0^-)|| < \delta$, the corresponding solution $q(t), \dot{q}(t)$ of (26a), under conditions (11), is such that $||q(t) - q_R|| < \varepsilon$ and $||\dot{q}(t^-)|| < \varepsilon, ||\dot{q}(t^+)|| < \varepsilon$ for all times $t \ge 0$.

2) For each real number $\delta > 0$, with $\delta$ being arbitrarily large, and for any initial condition $(q(0), \dot{q}(0^-)) \in C, \lambda_i(0^-) \in \mathbb{R}, \gamma_i(0) = \sqrt{-f_i(q(0))}, i = 1, 2, \dots, m$, satisfying $||q(0) - q_R|| < \delta$ and $||\dot{q}(0^-)|| < \delta$, the corresponding solution $q(t), \dot{q}(t), \lambda_i(t), \gamma_i(t) = \sqrt{-f_i(q(t))}, i = 1, 2, \dots, m$, of (26), under conditions (11), is such that

$$ \lim_{t \to +\infty} ||q(t) - q_R|| = 0 \qquad (29a) $$

$$ \lim_{t \to +\infty} ||\dot{q}(t^-)|| = 0 \qquad (29b) $$

$$ \lim_{t \to +\infty} ||\dot{q}(t^+)|| = 0 \qquad (29c) $$

$$ \lim_{t \to +\infty} ||\dot{\lambda}_i(t^-) - \dot{\lambda}_{i,R}|| = 0, \quad i = 1, 2, \dots, m, \quad (29d) $$

$$ \lim_{t \to +\infty} ||\dot{\lambda}_i(t^+) - \dot{\lambda}_{i,R}|| = 0, \quad i = 1, 2, \dots, m. \quad (29e) $$

**Proof:** First of all, note that all the solutions of (26), under conditions (11), in the variables $q(t)$ and $\dot{q}(t)$ are independent of the initial conditions $\lambda_i(0^-), i = 1, 2, \dots, m$, chosen for the Lagrange multipliers. In the following, denote by $(s_1(t; p^1, p^2), s_2(t; p^1, p^2))$ the solution of (26), under conditions (11), in the variables $(q, \dot{q})$, starting from the initial conditions $(q(0), \dot{q}(0^-)) = (p^1, p^2), \gamma_i(0) = \sqrt{-f_i(p^1)}, i = 1, 2, \dots, m$, with $(p^1, p^2) \in C$; sometimes, the shortening $(q(t), \dot{q}(t))$ will be used instead of $(s_1(t; q(0), \dot{q}(0^-)), s_2(t; q(0), \dot{q}(0^-)))$. It is pointed out that at the impact times $t_c, \dot{q}(t) = s_2(t; q(0), \dot{q}(0^-))$ is not defined, and only the limits $\dot{q}(t_c^-)$ and $\dot{q}(t_c^+)$ have meaning. By Assumption 2, at each time $t$ different from an impact time, $s_1(t; p^1, p^2)$ and $s_2(t; p^1, p^2)$ are continuous functions of $(p^1, p^2) \in C$. It is stressed again that, by (26c) and Assumption 2, if $q(0) \in A$ and $\gamma_i(0) = \sqrt{-f_i(q(0))}, i = 1, 2, \dots, m$, then $q(t) \in A$ for all times $t \ge 0$, for all the values of the other initial conditions $\dot{q}(0^-) \in \mathbb{R}^n, \lambda_i(0^-) \in \mathbb{R}, i = 1, 2, \dots, m$, for which the solution exists and is unique.

In the remainder of the proof, even though not explicitly mentioned, it will be assumed that, for any $q(0) \in A$, the initial conditions of the Valentine variables are taken such that $\gamma_i(0) = \sqrt{-f_i(q(0))}, i = 1, 2, \dots, m$.

Consider

$$ V(q, \dot{q}) := T(q, \dot{q}) + U(q) - U(q_R) \quad (30) $$

which is a positive definite function of $(q, \dot{q}) \in C$ about $(q_R, 0)$, by Assumptions 1, 3, and 4 (note that $V$ is only a positive semidefinite function of $q, \dot{q}, \lambda_i, i = 1, 2, \dots, m$, about $q_R, 0, \lambda_{i,R}, i = 1, 2, \dots, m$). It is stressed that the Erdmann-Weierstrass corner condition (11a) implies that $V(q(t), \dot{q}(t))$ is a continuous, piecewise smooth function of $t$. The total time derivative of $V(q(t), \dot{q}(t))$ along the solutions of (26), under conditions (11), is

$$ \dot{V}(q(t), \dot{q}(t^-)) = -\dot{q}^T(t^-)K_v \dot{q}(t^-) + \sum_{i=1}^{m} \dot{\lambda}_i(t^-) J_i^T(q(t)) \cdot \dot{q}(t^-) \quad (31a) $$

$$ \dot{V}(q(t), \dot{q}(t^+)) = -\dot{q}^T(t^+) K_v \dot{q}(t^+) + \sum_{i=1}^{m} \dot{\lambda}_i(t^+) J_i^T(q(t)) \cdot \dot{q}(t^+) \quad (31b) $$

Multiplying both sides of (26c), rewritten with $t = \tau$, by $\dot{\lambda}_i(\tau)$ and taking the limit of the resulting expression for $\tau \to t^-$ and $\tau \to t^+$, respectively, by (26b) and (11c), it is easy to see that $\dot{\lambda}_i(t^-) J_i^T(q(t)) \dot{q}(t^-) = 0$ (and, respectively, $\dot{\lambda}_i(t^+) J_i^T(q(t)) \dot{q}(t^+) = 0$), for all times $t$ (possibly coincident with the impact times) and for all $i \in \{1, 2, \dots, m\}$, whence that (31) can be recast as follows for all times $t$:

$$ \dot{V}(q(t), \dot{q}(t^-)) = -\dot{q}^T(t^-) K_v \dot{q}(t^-) \quad (32a) $$

$$ \dot{V}(q(t), \dot{q}(t^+)) = -\dot{q}^T(t^+) K_v \dot{q}(t^+) \quad (32b) $$

which are semidefinite negative functions of $(q(t), \dot{q}(t^-)), (q(t), \dot{q}(t^+)) \in C$ about $(q_R, 0)$.

For each real number $\delta > 0$, define $\beta(\delta)$ as the limit superior of $V(q, \dot{q})$ for $(q, \dot{q}) \in C$ such that $||q - q_R|| \le \delta$ and $||\dot{q}|| \le \delta$, which is finite by the continuity of $V(q, \dot{q})$ with respect to $(q, \dot{q}) \in C$ and by the compactness of the domain chosen in $C$. For each real number $\varepsilon > 0$, there exists a real number $\delta > 0$ such that $\beta(\delta) < \frac{1}{2} \min(a, b)\varepsilon^2$. This $\delta$ can always be found because $\frac{1}{2} \min(a, b)\varepsilon^2$ is positive and $\beta(\delta)$ tends to zero as $\delta$ goes to zero. Suppose that $(q(0), \dot{q}(0^-)) \in C, ||q(0) - q_R|| \le \delta$ and $||\dot{q}(0^-)|| \le \delta$; then, $V(q(0), \dot{q}(0^-)) \le \beta(\delta) < \frac{1}{2} \min(a, b)\varepsilon^2$. But, by (32), since $\dot{V}(q(t), \dot{q}(t^-)) \le 0, \dot{V}(q(t), \dot{q}(t^+)) \le 0$ for all times $t \ge 0$, by the continuity of $V(q(t), \dot{q}(t))$ with respect to time $t$, it follows that

$$ V(q(t), \dot{q}(t)) \le V(q(0), \dot{q}(0^-)) < \frac{1}{2} \min(a, b)\varepsilon^2, \quad \forall t \in \mathbb{R}, t \ge 0. \quad (33) $$

Now, taking into account Assumptions 1, 3, and 4

$$ V(q(t), \dot{q}(t^-)) \ge \frac{1}{2} a ||q(t) - q_R||^2 + \frac{1}{2} b ||\dot{q}(t^-)||^2, \quad \forall t \in \mathbb{R}, t \ge 0 \quad (34a) $$

$$ V(q(t), \dot{q}(t^+)) \ge \frac{1}{2} a ||q(t) - q_R||^2 + \frac{1}{2} b ||\dot{q}(t^+)||^2, \quad \forall t \in \mathbb{R}, t \ge 0. \quad (34b) $$
Inequalities (33) and (34) together imply that

$$
\begin{align*}
\|q(t) - q_R\|^2 &\le \varepsilon^2, && \forall t \in \mathbb{R}, t \ge 0 \\
\|\dot{q}(t^-)\|^2 &\le \varepsilon^2, && \forall t \in \mathbb{R}, t \ge 0 \\
\|\dot{q}(t^+)\|^2 &\le \varepsilon^2, && \forall t \in \mathbb{R}, t \ge 0
\end{align*}
$$

which, taking the square root, prove item 1) of the theorem.
For each real number $c>0$, consider the set

$$
\Omega_V(c) := \{(q, \dot{q}) \in C : V(q, \dot{q}) \le c\}
$$

and define the level set $L_V(c)$ as the connected component of $\Omega_V(c)$ containing the point $(q_R, 0)$ of $C$. By Assumptions 1, 3, and 4, $V(q, \dot{q})$ is a radially unbounded function of $(q, \dot{q})$, whence set $\Omega_V(c)$ is bounded, as well as $L_V(c)$, for any real number $c>0$. Since the time behavior of $(q(t), \dot{q}(t))$ is independent of the initial conditions chosen for the Lagrange multipliers, by a reasoning similar to the one used for the proof of item 1) of the theorem, it is easy to see that if $(q(0), \dot{q}(0^-)) \in L_V(c)$, then $(q(t), \dot{q}(t^-))$, $(q(t), \dot{q}(t^+)) \in L_V(c)$ for all times $t \ge 0$. For each $(q(0), \dot{q}(0^-)) \in L_V(c)$, a point $(p^1, p^2) \in C$ is a limit point of $(q(t), \dot{q}(t))$ if there exists a sequence $\{t_i\}_{i \in \mathbb{N}}$ of times $t_i$ (possibly coincident with the impact times) such that

$$
\begin{align*}
\lim_{i \to +\infty} t_i &= +\infty \\
\lim_{i \to +\infty} q(t_i) &= p^1 \\
\lim_{i \to +\infty} \dot{q}(t_i^-) &= p^2.
\end{align*}
$$

For each $(q(0), \dot{q}(0^{-})) \in L_V(c)$, the set constituted by all the limit points of $(q(t), \dot{q}(t))$ is called the limit set of $(q(t), \dot{q}(t))$ and is denoted hereafter by $\Lambda(q(0), \dot{q}(0^{-}))$. Now, some properties of $\Lambda(q(0), \dot{q}(0^{-}))$ are stated and proven.

Since, for each $(q(0), \dot{q}(0^{-})) \in L_V(c)$, the corresponding solution $(q(t), \dot{q}(t))$ belongs to $L_V(c)$ (whence it is bounded) for all times $t \ge 0$, then the corresponding limit set $\Lambda(q(0), \dot{q}(0^{-}))$ is: 1) nonempty; 2) bounded; and 3) closed, as shown in the following.

As for 1), since $(q(t), \dot{q}(t))$ is bounded, then for any sequence $\{t_i\}_{i \in \mathbb{N}}$ of times $t_i$ (possibly coincident with the impact times) approaching $+\infty$ as $i$ goes to infinity, the corresponding sequence $\{(q(t_i), \dot{q}(t_i^-))\}_{i \in \mathbb{N}}$ is bounded, whence it contains a convergent subsequence. By the definition of $\Lambda(q(0), \dot{q}(0^{-}))$, the limit of this convergent subsequence must belong to $\Lambda(q(0), \dot{q}(0^{-}))$, which is, therefore, nonempty.

As for 2), the boundedness of $\Lambda(q(0), \dot{q}(0^{-}))$ follows trivially by the boundedness of $(q(t), \dot{q}(t))$.

As for 3), let $\{(p_i^1, p_i^2)\}_{i \in \mathbb{N}}$ be a sequence of points of $\Lambda(q(0), \dot{q}(0^{-}))$ converging to $(p^1, p^2) \in C$; for $\Lambda(q(0), \dot{q}(0^{-}))$ to be closed, it must be shown that $(p^1, p^2) \in \Lambda(q(0), \dot{q}(0^{-}))$. Let the real numbers $\varepsilon > 0$ and $T < +\infty$ be arbitrary; a time $t \ge T$ (possibly coincident with an impact time) must be found such that the following relations hold for such a $t$:

$$
\|p^1 - q(t)\| < \varepsilon, \quad \|p^2 - \dot{q}(t^-)\| < \varepsilon. \quad (35)
$$

First, choose an index $i$ such that

$$
\|p^1 - p_i^1\| < \frac{\varepsilon}{2}, \quad \|p^2 - p_i^2\| < \frac{\varepsilon}{2}. \qquad (36)
$$

Such an index exists because $(p_i^1, p_i^2)$ tends to $(p^1, p^2)$ as index $i$ goes to infinity. Secondly, choose a time $t \ge T$ such that

$$
\|p_i^1 - q(t)\| < \frac{\varepsilon}{2}, \quad \|p_i^2 - \dot{q}(t^-)\| < \frac{\varepsilon}{2}. \quad (37)
$$

Such a time exists because $(p_i^1, p_i^2) \in \Lambda(q(0), \dot{q}(0^-))$. Inequalities (36) and (37) together yield (35), namely the closure of $\Lambda(q(0), \dot{q}(0^-))$.

Now, define the distance $d((r^1, r^2), \Lambda(q(0), \dot{q}(0^{-})))$ between a point $(r^1, r^2) \in C$ and the nonempty, closed, and bounded set $\Lambda(q(0), \dot{q}(0^{-}))$ as follows:

$$
\begin{align}
d((r^1, r^2), \Lambda(q(0), \dot{q}(0^{-}))) & \nonumber \\
&:= \min_{(p^1, p^2) \in \Lambda(q(0), \dot{q}(0^{-}))} \left( \max \left( \|r^1 - p^1\|, \|r^2 - p^2\| \right) \right) \tag{38}
\end{align}
$$

where $\| \cdot \|$ is the Euclidean norm of the vector $\cdot$ at argument; such a minimum exists (whence $d((r^1, r^2), \Lambda(q(0), \dot{q}(0^{-})))$ is well defined) because $\Lambda(q(0), \dot{q}(0^{-}))$ is nonempty, closed, and bounded, by the continuity of the Euclidean norm.

For each $(q(0), \dot{q}(0^{-})) \in L_V(c)$, taking into account that $(q(t), \dot{q}(t))$ is bounded, then

$$
\lim_{t \to +\infty} d((q(t), \dot{q}(t^{-})), \Lambda(q(0), \dot{q}(0^{-)})) = 0. \quad (39)
$$

As a matter of fact, if (39) is false, then there exists a real number $\varepsilon > 0$ and a sequence $\{t_i\}_{i \in \mathbb{N}}$ of times $t_i$ (possibly coincident with the impact times) approaching infinity as index $i$ goes to infinity, such that

$$
d((q(t_i), \dot{q}(t_i^-)), \Lambda(q(0), \dot{q}(0^{-}))) \ge \varepsilon, \quad \forall i \in \mathbb{N}. \quad (40)
$$

However, since the sequence $\{(q(t_i), \dot{q}(t_i^-))\}_{i \in \mathbb{N}}$ is bounded, it contains a convergent subsequence, whose limit must belong to $\Lambda(q(0), \dot{q}(0^{-}))$ by definition, which is contradiction of (40). As a consequence, (39) holds.

Finally, for each $(q(0), \dot{q}(0^{-})) \in L_V(c)$, taking into account that $(q(t), \dot{q}(t))$ is bounded, for any $(\dot{q}(0), \dot{q}(0^{-})) \in \Lambda(q(0), \dot{q}(0^{-}))$, then the corresponding solution $(\dot{q}(t), \dot{q}(t)) = (s_1(t; \dot{q}(0), \dot{q}(0^{-})), s_2(t; \dot{q}(0), \dot{q}(0^{-)}))$ satisfies $(\dot{q}(t), \dot{q}(t^{-})), (\dot{q}(t), \dot{q}(t^{+})) \in \Lambda(q(0), \dot{q}(0^{-}))$ for all times $t \ge 0$. As a matter of fact, since $(\dot{q}(0), \dot{q}(0^{-})) \in \Lambda(q(0), \dot{q}(0^{-}))$, there exists a sequence $\{t_i\}_{i \in \mathbb{N}}$ of times $t_i$ (possibly coincident with the impact times) approaching infinity as index $i$ goes to infinity, such that

$$
\begin{align*}
\lim_{i \to +\infty} s_1(t_i; q(0), \dot{q}(0^{-})) &= \dot{q}(0) \\
\lim_{i \to +\infty} s_2(t_i; q(0), \dot{q}(0^{-})) &= \dot{q}(0^{-}).
\end{align*}
$$

Then, by the continuity of $s_1(t; p^1, p^2)$, $s_1(t; p^1, p^2)$ with respect to $(p^1, p^2) \in C$, for an arbitrary time $\tau \ge 0$ (different
from an impact time)

$$
\begin{align*}
\lim_{i \to +\infty} & s_1(\tau; s_1(t_i; q(0), \dot{q}(0^{-})), s_2(t_i; q(0), \dot{q}(0^{-}))) \\
&= s_1(\tau; \hat{q}(0), \dot{\hat{q}}(0^{-})) \\
&= \hat{q}(\tau), \\
\lim_{i \to +\infty} & s_2(\tau; s_1(t_i; q(0), \dot{q}(0^{-})), s_2(t_i; q(0), \dot{q}(0^{-}))) \\
&= s_2(\tau; \hat{q}(0), \dot{\hat{q}}(0^{-})) \\
&= \dot{\hat{q}}(\tau).
\end{align*}
$$

Since $s_1, s_2$ are independent of the initial time and of the
initial conditions chosen for the Lagrange multipliers, one has
(by taking the limit for $\tau \to t^-$)

$$
\begin{align*}
& s_1(t; s_1(t_i; q(0), \dot{q}(0^{-})), s_2(t_i; q(0), \dot{q}(0^{-}))) \\
&= s_1(t + t_i; q(0), \dot{q}(0^{-})) \\
&= q(t + t_i), \\
& s_2(t^{-}; s_1(t_i; q(0), \dot{q}(0^{-})), s_2(t_i; q(0), \dot{q}(0^{-}))) \\
&= s_2((t + t_i)^{-}; q(0), \dot{q}(0^{-})) \\
&= \dot{q}((t + t_i)^{-}).
\end{align*}
$$

Such relationships show the existence of a sequence {$\tau_i$}$_{i \in \mathbb{N}}$ of times $\tau_i := t+t_i$ (possibly coincident with the impact times) approaching infinity as index $i$ goes to infinity, such that

$$
\begin{align*}
\lim_{i \to +\infty} q(\tau_i) &= \hat{q}(t) \\
\lim_{i \to +\infty} \dot{q}(\tau_i^-) &= \dot{q}(t^-)
\end{align*}
$$

which show that $(\hat{q}(t), \dot{\hat{q}}(t^-)) \in \Lambda(q(0), \dot{q}(0^-))$ for all times $t \ge 0$; since $\lim_{\tau \to t^+} \dot{q}(\tau^-) = \dot{q}(t^+)$, the property $(\hat{q}(t), \dot{\hat{q}}(t^-)) \in \Lambda(q(0), \dot{q}(0^-))$ and the fact that $\Lambda(q(0), \dot{q}(0^-))$ is closed imply $(\hat{q}(t), \dot{\hat{q}}(t^+)) \in \Lambda(q(0), \dot{q}(0^-))$ as to be proven.

Consider now what happens to function $V(q(t), \dot{q}(t))$. Since $\dot{V}(q(t), \dot{q}(t^-)) \le 0$ and $\dot{V}(q(t), \dot{q}(t^+)) \le 0$ for all times $t \ge 0$ (i.e., since $\dot{V}(q(t), \dot{q}(t^-))$ and $\dot{V}(q(t), \dot{q}(t^+))$ have constant nonpositive sign for all times), the nonnegative function $V(q(t), \dot{q}(t))$ is monotonic and has a definite limit as time $t$ goes to infinity; in addition, $\dot{V}(q(t), \dot{q}(t^-))$ and $\dot{V}(q(t), \dot{q}(t^+))$ tend to zero as time $t$ goes to infinity, at least almost everywhere. For each $(q(0), \dot{q}(0^-)) \in L_V(c)$, let $(p^1, p^2)$ be a point of $\Lambda(q(0), \dot{q}(0^-))$. Then, by definition, there exists a sequence $\{t_i\}_{i \in \mathbb{N}}$ of times $t_i$ (possibly coincident with the impact times) approaching infinity as index $i$ goes to infinity, such that

$$
\begin{align*}
\lim_{i \to +\infty} q(t_i) &= p^1 \\
\lim_{i \to +\infty} \dot{q}(t_i^-) &= p^2 \\
\lim_{i \to +\infty} \dot{V}(q(t_i), \dot{q}(t_i^-)) &= 0.
\end{align*}
$$

By the continuity of $\dot{V}(q, \dot{q})$ with respect to $(q, \dot{q}) \in C$,
one has $\dot{V}(p^1, p^2) = 0$; the arbitrariness of $(p^1, p^2) \in$
$\Lambda(q(0), \dot{q}(0^-))$ shows that $\dot{V}(q, \dot{q}) = 0$ along $\Lambda(q(0), \dot{q}(0^-))$.
Since, by Assumptions 1-4, $q(t) = q_R$, $\dot{q}(t) = 0$, $\lambda_i(t) =$
$\lambda_{i,R}t + \lambda_i(0^-)$, $\gamma_i(t) = \sqrt{-f_i(q_R)}$, $i = 1, 2, \dots, m$, $t \in$
$\mathbb{R}$, $t \ge 0$, are (for arbitrary $\lambda_i(0^-) \in \mathbb{R}$, $i = 1, 2, \dots, m$) the
only solutions of (26), under conditions (11), for which $\dot{V}(t) =$

0 identically, then for any $(q(0), \dot{q}(0^{-})) \in L_V(c)$ and for
any $\lambda_i(0^{-}) \in \mathbb{R}, i = 1, 2, \dots, m$, the corresponding solution
$q(t), \dot{q}(t), \lambda_i(t), \gamma_i(t) = \sqrt{-f_i(q(t))}, i = 1, 2, \dots, m, t \in$
$\mathbb{R}, t \ge 0$, of (26), under conditions (11), is, by the previous
reasoning and by (11a), such that (29) holds. The proof
of item 2) of the theorem is completed by observing that
$(q(0), \dot{q}(0^{-})) \in C, \|q(0) - q_R\| \le \delta$ and $\|\dot{q}(0^{-})\| \le \delta$ imply
$(q(0), \dot{q}(0^{-})) \in L_V(c)$, if one takes $c \ge \frac{1}{2}(\bar{a} + \bar{b})\delta^2$. $\square$

With a little abuse of terminology, property 1) of Theorem 5 states the “stability” of the solution of (26), under conditions (11), with reference only to the components $q(t)$ and $\dot{q}(t)$, whereas property 2) of Theorem 5 states the “global attractiveness” of such a solution as expressed by (29) (note that (29d) and (29e) do not imply that $\|\lambda_i(t) - \dot{\lambda}_{i,R}t\| \to 0$ as $t \to +\infty$).

Remark 4: Assumptions 3 and 4 can be weakened by requiring that the algebraic equations (27) have an isolated (instead of unique) solution $q = q_R$, $\dot{\lambda}_i = \dot{\lambda}_{i,R}$, $\gamma_i = \sqrt{-f_i(q_R)}$, $i = 1, 2, \dots, m$, with $q_R \in \mathcal{A}$ and $\dot{\lambda}_{i,R}$ being nonnegative, $i = 1, 2, \dots, m$, and that inequalities (28) hold for all $q \in \mathcal{E} \subset \mathcal{A}$, with $\mathcal{E} = \mathcal{D} \cap \mathcal{A}$ and $\mathcal{D} \subset \mathbb{R}^n$ being a sufficiently small neighborhood of $q = q_R$ (instead of the whole $\mathcal{A}$); Assumption 4 can be weakened further by using some other comparison functions $\varphi(r)$ instead of $r^2$. Then, Theorem 5 still holds under the Assumptions 3 and 4 thus weakened if, in item 1) of the theorem, the phrase “for each real number $\delta > 0$, with $\delta$ being arbitrarily large,” is replaced by the phrase “there exists a sufficiently small real number $\delta > 0.” The proof of the theorem thus modified is pretty much the same as the proof of Theorem 5, and is omitted for the sake of brevity. $\square$

Remark 5: The mechanical system so far considered has been assumed to be subject to conservative forces (the forces that can be derived from a potential energy), control forces, and reaction forces due to the inequality constraints. In reality, all the mechanical systems have always inherent damping. Nevertheless, all the previous analysis is still valid. As a matter of fact, the presence of internal damping will, in general, improve the performance of the derivative control law (25), with respect to stability. In all the cases in which the damping is to be considered, it is possible to add to the Euler–Lagrange equation (10a) a dissipation term $\partial R(\dot{q})/\partial \dot{q}$, obtaining

$$
\frac{d}{dt} \frac{\partial L}{\partial \dot{q}} - \frac{\partial L}{\partial q} + \sum_{i=1}^{m} \dot{\lambda}_i J_i(q) + \frac{\partial R(\dot{q})}{\partial \dot{q}} = Eu
$$

where $R(\dot{q})$ is the Rayleigh dissipation function [33], which
is a positive semidefinite function of $\dot{q}$. The corresponding
(26a) becomes (for fully actuated mechanical systems, i.e., for
$p = n$ and $E = I$)

$$
\frac{d}{dt} \frac{\partial L}{\partial \dot{q}} - \frac{\partial L}{\partial q} + \sum_{i=1}^{m} \dot{\lambda}_i J_i(q) = - \frac{\partial R(\dot{q})}{\partial \dot{q}} - K_v \dot{q}. \quad (41)
$$

If the dissipation term $\partial R(\dot{q})/\partial \dot{q}$ is proportional to the gener-
alized velocity $\dot{q}$ (i.e., if $(\partial R(\dot{q})/\partial \dot{q}) = D\dot{q}$, with $D$ being a
$n$-dimensional positive semidefinite square matrix), then (41)
is equal to (26a) with the positive definite square matrix $K_v$
substituted by the positive definite square matrix $K_v + D$.
Whence, Theorem 5 still holds in case of such an internal dissipation. □

**Remark 6:** Theorem 5 gives a tool for the control of nonsmooth impacts in mechanical systems, by the derivative control law (25), provided that such mechanical systems are fully actuated (i.e., when $p = n$ and $E = I$). Such an assumption is sufficient to show that $q(t) = q_R, \dot{q}(t) = 0, \dot{\lambda}_i(t) = \dot{\lambda}_{i,R}, \gamma_i(t) = \sqrt{-f_i(q_R)}, i = 1, 2, \dots, m$, is the only solution of system (26), under conditions (11), such that $\dot{q}^T(t)EK_vE^T\dot{q}(t) = 0$ identically, but of course it is not necessary. Theorem 5 still holds for the general case $p \le n$ if the phrase “if the mechanical system under consideration is fully actuated (i.e., if $p = n$ and $E = I$),” is replaced by the phrase “if system (26), under conditions (11), has $q(t) = q_R, \dot{q}(t) = 0, \dot{\lambda}_i(t) = \dot{\lambda}_{i,R}, \gamma_i(t) = \sqrt{-f_i(q_R)}, i = 1, 2, \dots, m$, as the only solution such that $\dot{q}^T(t)EK_vE^T\dot{q}(t) = 0$ identically.” The proof of the theorem thus modified is exactly the same as the proof of Theorem 5. □

**Remark 7:** It is stressed that Assumption 4 requires that $q = q_R$ is a minimum (at least local, in its weakened version) of $U(q)$ upon $A$, which, in general, may not be true. Nevertheless, Theorem 5 can be applied for the control of nonsmooth impacts in general mechanical systems, provided that a preliminary feedback control law from the natural outputs is applied to the mechanical system under consideration. Suppose that $p = n$; in such a case, the natural outputs $y$ coincide with the generalized coordinates $q$. Let $\tilde{U}(q)$ be a function such that Assumptions 2–4 hold with $U(q)$ replaced by $\tilde{U}(q)$. Then, consider the following preliminary feedback control law:

$$u = \frac{\partial U(q)}{\partial q} - \frac{\partial \tilde{U}(q)}{\partial q} + \tilde{u} \quad (42)$$

where $\tilde{u}$ is the vector of the new control forces. Taking into account that $L = T - U$, (10a), under control law (42), can be rewritten as follows (remember that the assumption $p = n$ implies $E = I$):

$$\frac{d}{dt} \frac{\partial \tilde{L}}{\partial \dot{q}} - \frac{\partial \tilde{L}}{\partial q} + \sum_{i=1}^{m} \dot{\lambda}_i J_i(q) = E \tilde{u} \quad (43)$$

where $\tilde{L}:=T-\tilde{U}$, whereas (10b) and (10c) remain unchanged. Equations (43), (10b), and (10c), under conditions (11), can be regarded as the Euler–Lagrange equations of a new mechanical system having $T$ as kinetic energy, $\tilde{U}$ as potential energy, and $\tilde{u}$ as the vector of the control forces, subject to the inequality constraints $f(q) \le 0$. Since Assumption 1–4 hold for such a mechanical system, Theorem 5 can be applied to prove the “stability” of $(q(t), \dot{q}(t))$ as stated in item 1), and the “global attractivity” of the solution $q(t) = q_R, \dot{q}(t) = 0, \dot{\lambda}_i(t) = \dot{\lambda}_{i,R}, \gamma_i(t) = \sqrt{-f_i(q_R)}, i = 1, 2, \dots, m$, as stated in item 2), when the further derivative feedback control law $\tilde{u} = -K_v \dot{q}(t)$, with $K_v$ being definite positive, is used.

If the mechanical system is not fully actuated (i.e., if $p < n$ and $E \neq I$), then the term $(\partial U(q)/\partial q)$ cannot be exactly compensated by feedback as in (42). In such a case, instead of (42), consider a linear feedback from the natural outputs

$$u = -K_p E^T q + \tilde{u} \quad (44)$$

where $K_p$ is a positive definite square matrix and $\tilde{u}$ is the vector of the new control forces. Also in this case, (10a), under control law (44), can be rewritten as in (43) with $\tilde{L}:=T-\tilde{U}$ and $\tilde{U}:=U + \frac{1}{2}q^T E K_p E^T q$. Since $E \neq I$ when $p < n$, the square matrix $EK_p E^T$ is only semidefinite positive, and Assumptions 3 and 4 are not necessarily true, with $U(q)$ being replaced by $\tilde{U}(q)$. However, if they hold and Remark 6 can be applied, then a properly modified version of Theorem 5 can be applied also in this case to prove the “stability” of $(q(t), \dot{q}(t))$ as stated at item 1), and the “global attractivity” of the solution $q(t) = q_R, \dot{q}(t) = 0, \dot{\lambda}_i(t) = \dot{\lambda}_{i,R}, \gamma_i(t) = \sqrt{-f_i(q_R)}, i = 1, 2, \dots, m$, as stated at item 2), when the further derivative feedback control law $\tilde{u} = -K_v E^T \dot{q}(t)$, with $K_v$ being definite positive, is used. □

**Remark 8:** In the simplest case of two bodies colliding along the common line connecting their centers of mass, the coefficient of restitution was defined by Newton as the ratio between the relative velocity of the two impacting bodies after the impact time and the relative velocity of the two impacting bodies just prior the impact time. However, for more general mechanical systems subject to (possible, multiple) impacts, it is most convenient to regard the coefficients of restitution (each one for each possible type of impact) as energy-loss functions, since all the impacts are basically processes involving energy exchange and energy transformation. In particular, if $T(q(t_c), \dot{q}(t_c^-))$ and $T(q(t_c), \dot{q}(t_c^+))$, with $T(q(t_c), \dot{q}(t_c^-)) \neq 0$ and $T(q(t_c), \dot{q}(t_c^+)) \le T(q(t_c), \dot{q}(t_c^-))$, are the kinetic energies immediately before and after the impact time $t_c$ due to the sole $i$th constraint, the coefficient of restitution $c_i(t_c)$ can be defined as follows:

$$c_i^2(t_c) := \frac{T(q(t_c), \dot{q}(t_c^+))}{T(q(t_c), \dot{q}(t_c^-))}$$

By its definition, each coefficient of restitution is nonnegative and less than or equal to one. When there is a collision between parts of a mechanical system, a portion of the original kinetic energy is converted in strain energy within the impacting parts. Subsequently, some fraction of the strain energy is reconverted back into the kinetic energy of the impacting parts, whereas the remainder of the energy is trapped within the mechanical system in the form of exciting various modes of vibration and/or is dissipated as energy of plastic deformation. The coefficient of restitution is commonly assumed to be constant and equal to one (as in the previous part of this paper), although its magnitude is dependent upon the geometry of the bodies involved in the impact, upon the presence or absence of slip at the point of contact, upon the duration of the impact, as well as upon some basic material properties of the bodies, such as Young modulus, mass density, elastic limits, etc. For instance, as for compact bodies such as spheres, the portion of kinetic energy lost because of excitation of various modes of vibration is quite small; for such bodies the coefficient of restitution is primarily controlled by the plastic deformation about the point of contact and by the friction within the impacting bodies. As both such phenomena are, in first approximation, linearly dependent on the impact relative velocity, the coefficient of restitution corresponding to the $i$th constraint can be taken as follows (for sufficiently small values
of $||\dot{q}(t^-)||$ for all times $t$ (possibly coincident with the impact times):

$$c_i(t) = \begin{cases} 1, & \text{if } f_i(q(t)) < 0, \\ 1, & \text{if } f_i(q(t)) = 0 \text{ and } J_i^T(q(t))\dot{q}(t^-) \le 0, \\ 1 - c_i J_i^T(q(t))\dot{q}(t^-), & \text{if } f_i(q(t)) = 0 \text{ and } J_i^T(q(t))\dot{q}(t^-) > 0, \end{cases} \quad (45)$$

where $c_i$ is a nonnegative real number. Then, the Erdmann-Weierstrass corner condition (11a) must be modified properly to take into account the coefficients of restitution

$$\frac{1}{2} \dot{q}^T(t_c^-) B(q(t_c)) \dot{q}(t_c^-) \prod_{i=1}^{m} c_i^2(t_c) = \frac{1}{2} \dot{q}^T(t_c^+) B(q(t_c)) \dot{q}(t_c^+). \quad (46)$$

Although (46) does not imply any longer the continuity of $V(q(t), \dot{q}(t))$ with respect to all times $t$, since $\prod_{i=1}^m c_i^2(t_c) \le 1$ by (45) for small values of $||\dot{q}(t^-)||$, the condition $V(q(t_c), \dot{q}(t_c^+)) \le V(q(t_c), \dot{q}(t_c^-))$ [owing to (45) and (46)], together with (32), implies that inequalities (33) are satisfied, whence Theorem 5 still holds locally if the Erdmann-Weierstrass corner condition (11a) is substituted by (46), i.e., if the coefficients of restitution are taken into account. Coefficients of restitution more general than (45) are to be used to deal with the global version of Theorem 5. $\square$

**Remark 9:** It is stressed that Theorem 4 (which gives sufficient conditions for the solution of system (26), under the Erdmann-Weierstrass corner conditions (11), to exist and to be continuously dependent on the initial conditions) cannot be applied when some of the coefficients of restitution are less than one. In particular, it seems that Theorem 4 cannot be easily amended to cover this case, for two different reasons: 1) the presence of one or more coefficients of restitution less than one could imply a finite accumulation point of the impact times (i.e., $\lim_{i \to +\infty} t_{c,i}(p^1, p^2) \neq +\infty$) and 2) even if the modified Erdmann-Weierstrass corner conditions (46), (11b), and (11c) admit a unique solution $\dot{q}(t_c^+) = \alpha(q(t_c), \dot{q}(t_c^-))$ for all $(q(t_c), \dot{q}(t_c^-)) \in \partial A \times \mathbb{R}^n$, function $\alpha(\cdot, \cdot)$ could be not continuous in $\partial A \times \mathbb{R}^n$. For these reasons, when there is a loss of kinetic energy at the impact times, other results guaranteeing existence and continuity of the solution of (26), under the modified Erdmann-Weierstrass corner conditions (46), (11b), and (11c), should be considered; the reader can benefit from the results given in [40]–[43]. $\square$

#### IV. CONTROL OF SMOOTH IMPACTS

Assume again that the natural outputs are the only coordinates available for feedback, as well as their time derivatives, and consider, first, the “derivative” control law (25). System (13), under the control law (25), becomes

$$\frac{d}{dt} \frac{\partial L}{\partial \dot{q}} - \frac{\partial L}{\partial q} + \sum_{i=1}^{m} \eta_i J_i(q) = -EK_v E^T \dot{q} \quad (47)$$

with $\eta_i$ given by (14).

**Assumption 5:** For each initial condition $q(0) \in \mathbb{R}^n$, $\dot{q}(0) \in \mathbb{R}^n$, system (47) has a unique solution upon the interval $[0, +\infty)$, which is a continuous function of $(q(0), \dot{q}(0))$.

Reference [44, Ths. 2.4.25 and 2.4.57] can be used for checking if Assumption 5 holds.

**Assumption 6:** The following algebraic equation in the unknown variable $q$:

$$\frac{\partial U}{\partial q} + \sum_{i=1}^{m} \eta_i J_i(q) = 0 \quad (48)$$

with $\eta_i$ given by (14), has a unique solution $q = q_R$.

**Assumption 7:** Under Assumption 6, there exist two real numbers $a, \bar{a}$, with $0 < a \le \bar{a} < +\infty$, such that the following inequalities hold in the whole (i.e., for all $q \in \mathbb{R}^n$):

$$\frac{1}{2} \bar{a} ||q - q_R||^2 \le U(q) - U(q_R) + U_f(q) - U_f(q_R) \\ \le \frac{1}{2} \bar{a} ||q - q_R||^2 \quad (49)$$

where $||\cdot||$ is the Euclidean norm of the vector $\cdot$ at argument.

**Remark 10:** Assumption 7 implies that $q = q_R$ is a global minimum of $U(q) + U_f(q)$ upon $\mathbb{R}^n$. If $q_R$ is such that $f(q_R) \le 0$, then $U_f(q_R) = 0$ and $q = q_R$ is a global minimum of the sole $U(q)$ upon $\mathbb{R}^n$, whence $(\partial U(q)/\partial q)|_{q=q_R} = 0$; this, by (14) and (48), implies that the $\eta_{i,R}$ given by (14), with $q$ being substituted by $q_R$, is equal to zero ($i = 1, 2, \dots, m$), whereas if $q_R$ is such that $f_i(q_R) > 0$ for at least one index $i$, then $(\partial U(q)/\partial q)|_{q=q_R}$ may be different from zero, as well as the corresponding $\eta_{i,R}$, $i = 1, 2, \dots, m$. $\square$

**Theorem 6:** Under Assumptions 1, 5, 6, and 7, if the mechanical system under consideration is fully actuated (i.e., if $p = n$ and $E = I$), then the solution $(q(t) = q_R, \dot{q}(t) = 0)$ of (47) is globally asymptotically stable.

*Proof:* Consider

$$V(q, \dot{q}) := T(q, \dot{q}) + U(q) - U(q_R) + U_f(q) - U_f(q_R) \quad (50)$$

which is a positive definite, radially unbounded, function of $(q, \dot{q}) \in \mathbb{R}^n \times \mathbb{R}^n$ about $(q_R, 0)$, by Assumptions 1, 6, and 7. The total time derivative of $V(q, \dot{q})$ along the solutions of (47) is (remember that $E = I$)

$$\dot{V}(q, \dot{q}) = -\dot{q}^T K_v \dot{q}. \quad (51)$$

Define $\mathcal{R} := \{(q, \dot{q}) \in \mathbb{R}^n \times \mathbb{R}^n: \dot{V} = 0\}$; then, by (51), it is easy to see that $\mathcal{R}$ does not contain any half-trajectory of (47) other than the trivial trajectory $(q(t) = q_R, \dot{q}(t) = 0)$. Then, by [44, Th. 5.3.79], the theorem is proven. $\square$

**Remark 11:** Assumptions 6 and 7 can be weakened by requiring that the algebraic equation (48) has an isolated (instead of unique) solution $q = q_R$, with $q_R \in \mathbb{R}^n$, and that inequalities (49) hold for all $q \in \mathcal{E}$, with $\mathcal{E}$ being a sufficiently small neighborhood of $q = q_R$ (instead of the whole $\mathbb{R}^n$); Assumption 7 can be weakened further, by using some other comparison functions $\varphi(r)$ instead of $r^2$. Then, Theorem 6 still holds under Assumptions 6 and 7 thus weakened, but stating only the local asymptotic stability of the trivial solution $(q(t) = q_R, \dot{q}(t) = 0)$. The proof of the theorem thus modified can be done by applying [44, Th. 5.3.77]. $\square$
**Remark 12:** A remark similar to Remark 5 can be stated to deal with damping terms, when the method of the penalty functions is used. Also in this case, the presence of internal damping will, in general, improve the performance of the derivative control law (25), with respect to stability. In particular, as for linear damping terms, since the dissipation in the contact period will be different from the dissipation in the noncontact period, the Rayleigh dissipation function can be taken as $R(\dot{q}) = \sum_{i=1}^{m} R_i(\dot{q})$, where

$$R_i(\dot{q}) := \begin{cases} \frac{1}{2} \dot{q}^T D_i \dot{q}, & \text{if } f_i(q) < 0, \\ \frac{1}{2} \dot{q}^T \hat{D}_i \dot{q}, & \text{if } f_i(q) \ge 0, \end{cases}$$

with $D_i, \hat{D}_i, i = 1, 2, \dots, m$, being positive semidefinite matrices. The dissipation matrices $\hat{D}_i, i = 1, 2, \dots, m$, relative to the contact period, can be used to characterize the loss of kinetic energy during the impact (the same role as the coefficients of restitution). □

**Remark 13:** Theorem 6 can be applied to fully actuated mechanical systems as $p = n$ and $E = I$ imply that $(q(t) = q_R, \dot{q}(t) = 0)$ is the only solution of (47) such that $\dot{q}^T(t)EK_vE^T\dot{q}(t) = 0$ identically. Theorem 6 still holds for the general case $p \le n$ if the phrase “if the mechanical system under consideration is fully actuated (i.e., if $p = n$ and $E = I$),” is replaced by the phrase “if system (47) has $(q(t) = q_R, \dot{q}(t) = 0)$ as the only solution such that $\dot{q}^T(t)EK_vE^T\dot{q}(t) = 0$ identically.” The proof of the theorem thus modified is exactly the same as the proof of Theorem 6.□

**Remark 14:** Theorem 6 can be used for the control of smooth impacts in mechanical systems, by the derivative control law (25). However, it is stressed that Assumption 7 requires that $q = q_R$ is a minimum (at least local, in its weakened version) of $\tilde{U}(q) + U_f(q)$ upon $\mathbb{R}^n$, which, in general, may not be true. Nevertheless, Theorem 6 can be applied for the control of smooth impacts in general mechanical systems, provided that a preliminary feedback control law from the natural outputs is applied to the mechanical system under consideration. Suppose that $p = n$; in such a case, the natural outputs $y$ coincide with the generalized coordinates $q$. Let $\tilde{U}(q)$ be a function such that Assumptions 5–7 hold with $U(q)$ replaced by $\tilde{U}(q)$. Then, consider the following preliminary feedback control law:

$$u = \frac{\partial U(q)}{\partial q} - \frac{\partial \tilde{U}(q)}{\partial q} + \tilde{u} \quad (52)$$

where $\tilde{u}$ is the vector of the new control forces. Taking into account that $L = T - U$, (13), under control law (52), can be rewritten as follows (remember that the assumption $p = n$ implies $E = I$):

$$\frac{d}{dt} \frac{\partial \tilde{L}}{\partial \dot{q}} - \frac{\partial \tilde{L}}{\partial q} + \sum_{i=1}^{m} \eta_i J_i(q) = E \tilde{u} \quad (53)$$

where $\eta_i$ is given by (14) and $\tilde{L}:=T - \tilde{U}$. Equation (53) can be regarded as the Euler–Lagrange equation of a new unconstrained mechanical system having $T$ as kinetic energy, $\tilde{U}+U_f$ as potential energy, and $\tilde{u}$ as the vector of the control forces. Since Assumption 1, 5, 6, and 7 hold for such a mechanical system, Theorem 6 can be applied to prove the

global asymptotic stability of $(q(t) = q_R, \dot{q}(t) = 0)$, when the further derivative feedback control law $\tilde{u} = -K_v \dot{q}(t)$, with $K_v$ being definite positive, is used.

If the mechanical system is not fully actuated (i.e., if $p < n$ and $E \neq I$), then the term $(\partial U(q)/\partial q)$ cannot be exactly compensated by feedback as in (52). In such a case, instead of (52), consider a linear feedback from the natural outputs

$$u = -K_p E^T q + \tilde{u} \quad (54)$$

where $K_p$ is a positive definite square matrix and $\tilde{u}$ is the vector of the new control forces. Also in this case, (13), under control law (54), can be rewritten as in (53) with $\tilde{L}:=T - \tilde{U}$ and $\tilde{U}:=U + \frac{1}{2}q^T E K_p E^T q$. Since $E \neq I$ when $p < n$, the square matrix $EK_p E^T$ is only semidefinite positive, and Assumptions 6 and 7 are not necessarily true, with $U(q)$ being replaced by $\tilde{U}(q)$. However, if they hold and Remark 13 can be applied, then a properly modified version of Theorem 6 can be applied also in this case to prove the global asymptotic stability of $(q(t) = q_R, \dot{q}(t) = 0)$, when the further derivative feedback control law $\tilde{u} = -K_v E^T \dot{q}(t)$, with $K_v$ being definite positive, is used. □

## V. EXPERIMENTAL RESULTS ABOUT THE IMPACT CONTROL OF A SINGLE-LINK ROBOT ARM

The mechanical system under consideration is constituted by a beam, which is situated in a plane, where an inertial frame $(x, y)$ is defined and is constrained by a hinge to rotate about one of its extremities, whereas the other extremity is completely free: the frame $(x, y)$ is defined so that its origin $(0, 0)$ coincides with the hinge. The control input is an external torque $u(t)$ exerted at the hinge (the system is fully actuated); the angular position $q(t)$ of the hinge is the natural output of this mechanical system. For the sake of simplicity, the vector of gravity is assumed to be perpendicular to the motion plane, so that the effects of the gravity force can be neglected. An infinitely rigid and massive obstacle is located in the plane of motion at a point $(x_0, y_0)$ so that the beam is constrained to satisfy the inequalities $2\pi + q_0 \ge q(t) \ge q_0$, with $q_0 := \arctan(y_0/x_0)$, for all times $t \ge 0$; assume, without loss of generality, that $q_0 > 0$. The admissible region is then

$$A := \{q \in \mathbb{R}: f_1(q) := q - q_0 - 2\pi \le 0, f_2(q) := q_0 - q \le 0\}$$

which is nonempty. The two constraints cannot be simultaneously satisfied with the equality signs, and their gradient vectors $J_1(q)$ and $J_2(q)$ are constant and equal to 1 and -1, respectively. The kinetic energy is given by $T = \frac{1}{2} I_r \dot{q}^2(t)$, where $I_r$ is the inertia of the beam: Assumption 1 holds with $\underline{b} = \bar{b} = I_r$. By the assumptions, the potential energy due to the conservative forces is equal to zero; then, the total potential energy is given by $U_t = -q(t)u(t)$. The Euler–Lagrange equation (10a) becomes

$$I_r \ddot{q}(t) + \dot{\lambda}_1(t) - \dot{\lambda}_2(t) = u(t) \quad (55)$$

where $\lambda_1(t), \lambda_2(t)$ are the Lagrange multipliers; similarly, the Euler–Lagrange equation (13) becomes

$$I_r \ddot{q}(t) + \eta_1(t) - \eta_2(t) = u(t) \quad (56)$$
<figure><img src="image_4.png" /><figcaption>Fig. 1. Measured angular positions.</figcaption></figure>

where

$$ 
\begin{aligned}
\eta_1 &:= \begin{cases} 0, & q - q_0 - 2\pi \le 0, \\ k_1(q - q_0 - 2\pi), & q - q_0 - 2\pi > 0, \end{cases} \\
\eta_2 &:= \begin{cases} 0, & q_0 - q \le 0, \\ k_2(q_0 - q), & q_0 - q > 0, \end{cases}
\end{aligned} 
\quad (57) $$

and $k_1, k_2 \gg 0$ are real numbers. Consider the following control law for both systems (55) and (56):

$$ u(t) = -K_p q(t) - K_v \dot{q}(t) \quad (58) $$

with $K_p, K_v > 0$; systems (55) and (56), under control law (58), can be recast as in (43) and (53), respectively, with $\tilde{U}(q) = \frac{1}{2} K_p q^2$

$$ I_r \ddot{q}(t) + \dot{\lambda}_1(t) - \dot{\lambda}_2(t) + K_p q(t) = -K_v \dot{q}(t) \quad (59a) $$

$$ I_r \ddot{q}(t) + \eta_1(t) - \eta_2(t) + K_p q(t) = -K_v \dot{q}(t). \quad (59b) $$

It is easy to see that Assumption 2 holds for system (59a), under the relevant Erdmann-Weierstrass corner conditions, whereas Assumption 5 holds for system (59b). As for system (59a), the algebraic equations (27), with $U$ being replaced by $\tilde{U}$, are

$$ 
\begin{aligned}
\dot{\lambda}_1 - \dot{\lambda}_2 + K_p q &= 0 \\
\gamma_1 \dot{\lambda}_1 &= 0 \\
\gamma_2 \dot{\lambda}_2 &= 0 \\
\gamma_1^2 + q - q_0 - 2\pi &= 0 \\
\gamma_2^2 + q_0 - q &= 0
\end{aligned} 
$$

which have the unique solutions $q_R = q_0, \dot{\lambda}_{R,1} = 0, \dot{\lambda}_{R,2} = K_p q_0, \gamma_{R,1} = \sqrt{2\pi}, \gamma_{R,2} = 0$, with $q_R \in A$ and $\dot{\lambda}_{R,1}, \dot{\lambda}_{R,2}$ being nonnegative; Assumption 3 holds for system (59a). Note that the positive value of $\dot{\lambda}_{R,2} = K_p q_0$ can be arbitrarily fixed

by choosing properly the value of $K_p$. As for system (59b), the algebraic equation (48), with $U$ being replaced by $\tilde{U}$, is

$$ \eta_1 - \eta_2 + K_p q = 0 $$

with $\eta_1, \eta_2$ given by (57), which has the unique solution $q_R = (k/(k+K_p))q_0, \eta_{R,1} = 0, \eta_{R,2} = (kK_p/(k+K_p))q_0$ (it is easy to see that $(k/(k+K_p))q_0 \to q_0$ and $\eta_{R,2} \to \dot{\lambda}_{R,2}$ as $k \to +\infty$); Assumption 6 holds for system (59b). Finally, as for system (59a), it is easy to see that $\tilde{U}(q) - \tilde{U}(q_0)$ is a global positive definite function about $q = q_0$ upon $A$ (a different comparison function $\varphi(r)$ must be used in this case instead of $r^2$), whereas as for system (59b), it is easy to see that $\tilde{U}(q) - \tilde{U}((k/(k+K_p))q_0) + \tilde{U}_f(q) - \tilde{U}_f((k/(k+K_p))q_0)$ with $\tilde{U}_f(q) = \tilde{U}_{f1}(q) + \tilde{U}_{f2}(q)$ given by

$$ 
\begin{aligned}
\tilde{U}_{f1}(q) &= \begin{cases} 0, & q - q_0 - 2\pi \le 0 \\ \frac{1}{2} k_1 (q - q_0 - 2\pi)^2, & q - q_0 - 2\pi > 0 \end{cases} \\
\tilde{U}_{f2}(q) &= \begin{cases} 0, & q_0 - q \le 0 \\ \frac{1}{2} k_2 (q_0 - q)^2, & q_0 - q > 0 \end{cases}
\end{aligned} 
$$

is a global positive definite function about $q = (k/(k+K_p))q_0$ (a different comparison function $\varphi(r)$ must be used also in this case instead of $r^2$). Therefore, Theorems 5 and 6 can be applied, showing the desired stability property.

The effectiveness of the control law (58) has been tested experimentally. The mechanical system used for the experiments has been developed at the Robotics and Industrial Automation Laboratory of the University of Rome and consists of a link made of hardened steel (having dimensions of 600 × 40 × 2 mm), clamped on an aluminum flange secured to the rotor of a direct drive dc motor. The measuring device consists of a set of local bending sensors, electrical strain-gauges located along the link at suitable positions, which are used for monitoring the possible deformations and forces due to impacts and contacts, and of an incremental optical encoder,
<figure><img src="image_5.png" /><figcaption>Fig. 2. Measured deformations.</figcaption></figure>

which is used to measure the joint displacement and velocity, as needed by the controller. The measurement and control system is completed with signal transducers, signal amplifiers, and an I/O board installed in a PC computer, equipped with signal processing and numerical computation software. The local deformation is measured at two points: in each of the two curvature measuring points, four strain-gauges have been located, connected in a Wheatstone bridge configuration. The forces due to impacts and contacts are measured at the tip of the beam, by the local deformation measured by additional strain-gauges. Because of the light intrinsic damping of the experimental system that has been used (direct drive motors have very light intrinsic damping), the assumption of absence of friction is well reproduced by the experimental environment. Two different experiments have been carried out. In Case a, the motion plane has been kept free of obstacles (i.e., the beam is completely free to rotate), whereas in Case b, an obstacle has been located in the motion plane at a certain position $(x_0, y_0)$ such that $q_0 = \arctan(y_0/x_0)$. Fig. 1 reports the time histories of the angular position of the motor in both Cases a and b in the closed interval $[0, 2]$. As can be seen, in Case a (the unconstrained one), the angular position $q(t)$ of the motor asymptotically goes to zero, where in Case b (the constrained one), the first impact between the beam and the obstacle occurs at a time $t_c \approx 0.25$ s and, after a sequence of other impacts one close to the other, the beam reaches the steady-state configuration in contact with the obstacle; the first impact has been sufficiently strong so to produce a certain deformation of the beam, and the motor has been able to overtake its steady-state limit. Fig. 2 shows the deformations recorded by strain-gauges, connected in a Wheatstone bridge configuration, at distances of 9.5 and 25.5 cm from the motor; in Case a, the deformations asymptotically go to zero, whereas, in Case b, the deformations asymptotically reach some steady-state values, due to the steady-state contact between the beam and the obstacle.

## VI. CONCLUSION

The main contribution of this paper is to have shown that the classical feedback control laws used for unconstrained mechanical systems can be used effectively for mechanical systems subject to inequality constraints, thus obtaining pretty much the same results, independently of the fact that smooth or nonsmooth impacts are considered. Experimental results have shown the effectiveness of these control laws. Future work will regard the tracking problem for mechanical systems subject to inequality constraints, as well as the regulation problem with a deadbeat transient response.

## REFERENCES

[1] B. Brogliato, *Nonsmooth Impact Mechanics*. London, U.K.: Springer-Verlag, 1996.

[2] J. A. Zukas, *Impact Dynamics*. New York: Wiley, 1982.

[3] R. J. Wasley, *Stress Wave Propagation in Solids, an Introduction*. New York: Marcel Dekker, 1973.

[4] W. Goldsmith, *Impact*. London, U.K.: E. Arnold, 1960.

[5] J. S. Rhinehart and J. Pearson, *Behavior of Metals Under Impulsive Loads*. Cleveland, OH: American Soc. Testing and Materials, 1954.

[6] J. D. Achenbach, *Wave Propagation in Elastic Solids*. New York: American Elsevier, 1975.

[7] J. Wittenberg, *Dynamics of Systems of Rigid Bodies*. B. G. Teubner, 1977.

[8] W. Johnson, *Impact Strength of Materials*. New York: Crane, Russak, 1972.

[9] U. S. Lindholm, *Mechanical Behavior of Materials Under Dynamic Loads*. New York: Springer-Verlag, 1968.

[10] N. J. Huffington, *Behavior of Materials Under Dynamic Loading*. New York: Amer. Soc. Mechanical Eng., 1965.

[11] R. Kinslolo, *High Velocity Impact Phenomenon*. New York: Academic, 1970.

[12] H. Kosly, *Stress Waves in Solids*. Fair Lawn, NJ: Oxford Univ. Press, 1953.

[13] W. B. Gevarter, "Basic relations for control of flexible vehicles," *AIAA J.*, vol. 8, pp. 666-672, Apr. 1972.

[14] B. Chapnik, G. Heppler, and J. Aplevich, "Modeling impact on a one-link flexible robotic arm," *IEEE Trans. Robot. Automat.*, vol. 7, pp. 479-488, Aug. 1991.
[15] A. Yigit, A. Ulsoy, and R. Scott, “Dynamics of a radially rotating beam with impact,” *ASME J. Vibration and Acoustics*, vol. 112, pp. 515–525, Jan. 1990.

[16] Q. F. Wei, P. S. Krishnaprasad, and W. P. Dayawansa, “Modeling of impact on a flexible beam,” in *Proc 32nd Conf. Decision and Control*, San Antonio, TX, vol. 2, pp. 1377–1382, Dec. 1993.

[17] B. Brogliato, P. Orhant, and R. Lozano, “On the transition phase in robotics—Part 1: Impact models and dynamics,” in *IFAC Symp. Robot Control*, Capri, Italy, Sept. 1994.

[18] ——, “On the transition phase in robotics—Part 2: Control,” in *IFAC Symp. Robot Control*, Capri, Italy, Sept. 1994.

[19] Y. F. Zheng and H. Hemami, “Mathematical modeling of a robot collision with its environment,” *J. Robot. Syst.*, vol. 2, no. 3, pp. 289–307, 1985.

[20] Y. Hurmuzlu and D. B. Marghitu, “Rigid body collision of planar kinematic chains with multiple contact points,” *Int. J. Robot. Res.*, vol. 13, pp. 82–92, Feb. 1994.

[21] Y. Hurmuzlu and T. H. Chang, “Rigid body collisions of a special class of planar kinematic chain,” *IEEE Trans. Syst., Man, and Cybern.*, vol. 22, pp. 964–971, Sept.–Oct. 1992.

[22] V. V. Kozlov and D. V. Treshchev, *Biliards: A Genetic Introduction to the Dynamics of Systems with Impacts*. Providence, RI: Amer. Math. Soc., 1991.

[23] G. Buttazzo and D. Percivale, “On the approximation of the elastic bounce problem on riemann manifolds,” *J. Differential Equations*, vol. 47, pp. 227–275, 1981.

[24] D. Percivale, “Uniqueness in the elastic bounce problem,” *J. Differential Equations*, vol. 56, pp. 206–215, 1985.

[25] R. Volpe and P. Khosla, “A theoretical and experimental investigation of impact control for manipulators,” *Int. J. Robotics Res.*, vol. 12, no. 4, pp. 351–365, 1993.

[26] K. Youcef-Toumi and D. A. Gutz, “Impact and force control: Modeling and experiments,” *ASME J. Dynamic Syst., Measurements, and Contr.*, vol. 116, pp. 89–98, 1994.

[27] T. Tarn, Y. Wu, N. Xi, and A. Isidori, “Force regulation and contact transition control,” *IEEE Contr. Syst.*, pp. 32–40, Feb. 1996.

[28] A. van der Schaft and J. M. Schumacher, “The complementarity-slakness class of hybrid systems,” *Math. Contr., Signal, and Syst.*, vol. 9, pp. 266–301, 1996.

[29] H. P. Huang and N. H. MacClamroch, “Time optimal control for a robotic contour following problem,” *IEEE J. Robot. Automat.*, vol. 4, pp. 140–149, 1988.

[30] J. K. Mills and D. M. Lokhorst, “Control of robotic manipulators during general task execution: A discontinuous control approach,” *Int. J. Robot. Res.*, vol. 12, no. 2, pp. 144–163, 1993.

[31] B. Brogliato, S. I. Niculescu, and P. Orhant, “On the control of finite-dimensional mechanical systems with unilateral constraints,” *IEEE Trans. Automat. Contr.*, vol. 42, no. 2, 1997.

[32] J. K. Mills and C. Nguyen, “Robotic manipulator collisions: Modeling and simulation,” *Trans. ASME J. Dyn. Syst. Meas. and Contr.*, vol. 114, no. 4, 1993.

[33] H. Goldstein, *Classical Mechanics*. Reading, MA: Addison-Wesley, 1980.

[34] F. A. Valentine, *The Problem of Lagrange with Differential Inequalities as Added Side Conditions*. Chicago, IL: Univ. Chicago Press, 1937.

[35] D. G. Luenberger, *Linear and Nonlinear Programming*. Reading, MA: Addison-Wesley, 1989.

[36] G. A. Bliss, *Lectures on the Calculus of Variation*. Chicago, IL: Univ. Chicago Press, 1946.

[37] W. Hahn, *Stability of Motion*. Berlin, Germany: Springer-Verlag, 1967.

[38] A. Ivanov, “On the stability of systems with nonretaining constraints,” *J. Applied Math. and Mechanics*, vol. 48, no. 5, 1984.

[39] L. B. Rapoport, “Stability of the equilibrium of a system with single-sided constraints and the sign-definiteness of a pencil of quadratic forms in a cone,” *J. Applied Math. and Mechanics*, vol. 56, no. 4, pp. 501–507, 1992.

[40] M. Carriero and E. Pascoli, “Il problema del rimbalzo unidimensionale e sue approssimazioni con penalizzazioni non convesse,” *Rendiconti di Matematica*, vol. 13, no. 1, pp. 541–553, 1980.

[41] L. Paoli and M. Schatzman, “Mouvement a un nombre fini de degrees de liberté avec contraintes unilaterales: Cas avec perte d'energie,” *Math. Modeling and Numerical Anal.*, vol. 27, no. 26, pp. 673–717, 1993.

[42] M. D. P. Monteiro-Marques, *Differential Inclusion in Nonsmooth Mechanical Problems: Shock and Dry Friction*. Boston, MA: Birkhauser, 1993.

[43] W. W. Schmaldeke, “Optimal control theory for nonlinear vector differential equations containing measures,” *SIAM J. Contr.*, vol. 3, no. 2, pp. 231–280, 1965.

[44] M. Vidyasagar, *Nonlinear Systems Analysis*. Englewood Cliffs, NJ: Prentice-Hall, 1993.

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

Antonio Tornambè received the “Dottore” degree in electrical engineering, summa cum laude, from Università’ di Roma La Sapienza in 1985.

He was a Researcher at Fondazione Ugo Bordoni, from 1985 to 1990, Research Associate at Università degli Studi di Roma “Tor Vergata,” from 1991 to 1992, Associate Professor of Industrial Automation at Politecnico di Torino, from 1992 to 1993, Associate Professor of Control Theory and of Operation Research at Università degli Studi di Roma Tre, and a Visiting Professor of System Theory at Università degli Studi di Siena. Currently he is a Professor at the Terza Università di Roma. He is the Editor of *Modeling and Control of Mechanisms and Robots* (Singapore, World Scientific, 1996) and *Modeling and Control of Mechanical Systems* (London, U.K.: Imperial College Press, 1997) and the author of *Discrete-Event System Theory: An Introduction*” (Singapore, World Scientific Publishing, 1995), *Mathematical Methods for System Theory* (Singapore, World Scientific Publishing, 1998), and *Modeling and Control of Mechanical Systems Subject to Constraints* (Boston, MA: Birkhauser, 1998). He is also the author of about 200 papers, published in international journals and presented at international conferences. His research interests include robotics and control theory.

Dr. Tornambè was the General Chairman of the 6th IEEE Mediterranean Conference, Alghero, Italy, June 1998. He has been an Associate Editor of the Conference Editorial Board of the IEEE Control System Society since 1995.