I Introduction
Modern robots are fast and strong, and, in some situations, their capabilities eclipse those of humans. However, when these robots interact with their environment, whether by manipulating objects or traversing over uneven surfaces, they do so with far less skill than their human counterparts. The speed and accuracy that robots display in static or structured settings disappears when faced with less structured tasks. Critical challenges facing the field lie in modeling, planning, and control of robots in these complex, multicontact settings, particularly for locomotion [38] and manipulation [18].
Rigidbody models of dynamics and contact (see Stewart [33] or Brogliato [4] for an overview) are widely used in robotics, as they can lead to far more tractable methods than approaches which explicitly attempt to capture the stiff interaction between objects. These approaches have also led to complementaritybased simulation schemes, such as [32, 1, 17, 13] and others. Recent research, using complementarity models, has also been conducted into multicontact optimal planning [22, 20, 23] and control [24, 12]. Similar applications have been seen for manipulation (e.g. [29]), including quasistatic approaches [7, 11]. When impacts occur, rigidbody models approximate the event as an instantaneous change in velocity due to an impulsive force.
The approaches above, now deeply ingrained within the robotics community, universally assume that it is possible to determine a single potential postimpact velocity, even during simultaneous multicontact. However, as observed in [37, 31, 36] and others, including recent analysis of robot locomotion [25], the resolution of simultaneous impacts is dependent upon the sequence in which they are resolved. Simulation schemes to this problem (e.g. [31, 36, 17, 19, 14]
) focus on generation of a single solution via a heuristic (symmetry
[17], potential energy [36], etc.). However, for many practical applications in robotics, it is not possible to create a model detailed enough to reliably disambiguate between the multiple potential solutions. Furthermore, even given such detail, this lack of uniqueness often represents an extreme sensitivity to initial conditions: slight perturbations in the initial state of the system might lead to different impact sequences.As the motivating examples in IIIA will demonstrate, simultaneous impacts are not limited to unlikely, pathological events but are, in fact, regular occurrences in robotics and require careful analysis. From the perspective of planning, learning and control, it is critical to understand the role of this nonuniqueness (alternatively, extreme sensitivity), as some of the broad challenges in executing dynamic, multicontact motion likely arise from these issues. For example, methods which use a simulator to learn or plan a motion may, unwittingly, be planning for an ambiguous, therefore unstable, outcome due to multicontact. Furthermore, as the set of these ambiguous outcomes is often nonconvex, it is insufficient to try to capture this sensitivity via simple models of uncertainty.
Many methods have been proposed for modeling single, frictional impacts (e.g. [26, 34, 6, 3], and others) along with recent datadriven models [9, 15], experimental validation [10], and efforts to translate multicontact simulated motions to real robots [35]. However, comparable results for simultaneous impacts have been limited to restricted models; addressing the complexity of Coulomb friction has been particularly difficult. Seghete and Murphey [30]
presented a model where solutions were guaranteed to exist, but assumed that impacts are frictionless and that contact normal vectors are linearly independent.
Burden et al. [5] studied discontinuous vector fields, with strong results and applications to robot impacts, but are similarly restricted to frictionless contact. Johnson et al. [16] treated a limited form of friction, but assumed that the effects of contact forces are inertially decoupled. For a quasistatic model, thus without impact, Halm and Posa [11] guaranteed existence of solutions for multicontact motion.This work extends Routh’s original graphical model [26]
to address simultaneous, inelastic impacts by permitting impulses to occur in arbitrary sequences. As a result, the model produces a setvalued map that captures the lack of uniqueness inherent in the problem. We believe this is the appropriate description for robotic planning and control, as motions that present as nonunique will, for physical systems, display extreme sensitivity to any errors in estimation or control. In contrast with prior literature, the presented model captures a broad class of frictional systems. In
III, we describe the model and a number of its theoretical properties and in IV we prove the key result that the impact model is guaranteed to terminate. To the best of the authors’ knowledge, this work presents the first known proof of existence of solutions for simultaneous frictional impact.Ii Background
We now introduce notation for and study the limiting behaviors of the frictional impact dynamics of rigid multibody systems. Denote the interior, closure, and convex hull of a set as , , and . We identify the norm and unit direction of a vector as and , respectively. We define the open radius ball in as . We denote as the vectors with strictly positive entries and definite a function to be positive definite if it is strictly positive on . For a singlevalued function and a setvalued function , we denote the image of under and as and respectively.
Iia Functional Analysis
The results herein are broadly derived from measure theory and functional analysis; for a thorough background, see Rudin [27, 28]. For a set , unless otherwise specified we equip with the standard Euclidean metric and norm, and integrals on are taken with respect to the Lebesgue measure. We denote the valued Lebesgue spaces and Sobolev space of degree one on a compact interval as and respectively, and particularly note that is a Hilbert space. The total time derivative of an absolutely continuous functions is taken in the Lebesgue sense and defined almost everywhere () such that
(1) 
Convergence of sequences of functions will be frequently used. Convergence of a sequence of functions to almost everywhere, uniformly, and with respect to the norm of a Hilbert space are denoted , , and respectively. Furthermore, weak convergence in is denoted as . It is well know that for any compact interval , is compactly embedded in , leading to a key result for the derivations in this work:
[Rellich] Let be a bounded sequence in for some compact interval . if (weakly), then (strongly).
IiB Differential Inclusions
The dynamics of many robots can be captured accurately with a system of ordinary differential equations (ODEs)
which relates , the state of the robot (typically a generalized notion of its position and velocity), to , a set of inputs (such as motor torques) that can be manipulated. However, the dynamics of rigid bodies under frictional contact present two additional complexities that this formulation cannot capture. Impacts between bodies induce instantaneous jumps in velocity that in general cannot described by an ODE (nonsmooth behaviors). Additionally, when contact occurs at many points, multiple frictional forces that obey Coulomb’s laws of friction may exist (nonunique behaviors). It is therefore useful to define an object that is more permissive than ODEs by allowing for the derivative at each state to lie in a set of possible values(2) 
As the derivative map associated with Coulomb friction may not be differentiable or even continuous in , conditions for a function to be considered a solution to the differential inclusion (2) are weakened from those of an ODE: For a compact interval , is a solution to the differential inclusion if is absolutely continuous and on . Denote the set of such solutions as . Solutions to initial value problems for (2) are defined in a similar manner: For a compact interval , denote the set of functions that are solutions to the initial value problem as . For example, consider the differential inclusion
(3) 
where is the setvalued unit direction function
(4) 
For any compact interval , the initial value problem admits the unique solution
(5) 
We note that is differentiable except at , thus no ODE would admit as a solution. In general, nonemptiness, regularity, and closure of are highly dependent on the structure of ; fortunately, frictional dynamics admit an upper semicontinuous (u.s.c.) structure that moderates the behavior of their solution sets: A function with values closed in is upper semicontinuous if with , , and , we have . [Aubin and Cellina [2]] Let and be a compact interval. If is uniformly bounded; u.s.c.; and closed, convex, and nonempty at all , is u.s.c. in . Furthermore as well as are nonempty and closed under uniform convergence.
Intuitively, a map is u.s.c. if its value at each is not significantly smaller than its value at any near . for example obeys all requirements of Proposition 5. We have already seen that is closed, nonempty, and convex. It is easy to see that if , then with . An illustration of this system as well as the function can be found in Figure 2.
IiC Frictional Impact Dynamics
The dynamics of robotic systems can be modeled as a system of rigid bodies experiencing frictional contact at up to points (for a thorough introduction, see [33] and [4]). The state of such a system can be represented by generalized coordinates and velocities . The continuous evolution is governed manipulator equations
(6) 
where is the generalized inertial matrix; encompasses Coriolis and gravitational forces; projects the velocity onto the contact normals; and projects onto to the contact tangents of the frictional contacts. We identify the behavior with a set of contacts , and identify each contact with it’s related vectors: row of and rows and of , denoted as and , respectively. Denote the collection of potential contact sets as , thus . We furthermore denote , the collection of sets of contacts of which are frictional. The worldframe contact normal and frictional forces and are constrained to lie within the Coulomb friction cone ; that is, for all and ,
(7)  
(8) 
where and are identified similarly to and and is the friction coefficient for the th contact. Additionally, we denote the lumped terms
(9)  
(10)  
(11) 
Intuitively, is the set of actively penetrating velocities, where impact is guaranteed to occur. are separating velocities, where no impact can occur. Note that , and velocities in this set may require impacts, as in Painlevé’s Paradox [33].
In this work, we focus on inelastic impulsive impacts. During such an impact, the velocity changes instantaneously. Letting represent an impulse, the relationship between pre and postimpact velocities, and , is
Coulomb friction poses challenges in computing , as an impact may cause stickslip transitions or the direction of slip to change. For a single contact , Routh [26] proposed a graphical method describing a path in velocity space (equivalently impulse space), from pre to postimpact velocities, where the path must differentially satisfy Coulomb friction. To briefly summarize Routh’s technique,

Increase the normal impulse with some slope .

Increment the tangential impulse with slope , satisfying to Coulomb friction, identical to (8) and depending on , the velocity after net impulse .

Terminate when the normal contact velocity vanishes^{1}^{1}1To permit resolutions to Painlevé’s Paradox, terminate only when consistency no longer requires an instantaneous change in velocity., , and take .
To later proceed to the multicontact case, we observe that this process could be modeled as an upper semicontinuous differential inclusion:
(12) 
where is equal to the net increment in velocity due to the “force” applied in steps 1) and 2) of Routh’s method. Since is constant during an impact, we will apply the transformation , leaving
(13) 
where we retain the use of for ease of notation. For any , we can associate a set of forces such that
(14) 
Note that for a frictionless contact (), this simplifies to
(15) 
A diagram depicting the resolution of a potential planar impacts is shown in Figure 2. Solutions may transition between sliding and sticking, and the direction of slip may even reverse as a result of each impact. While the path is piecewise linear in the planar case, this is not true in three dimensions.
From this point forward, in a slight abuse of notation, we will take to be the “time” during the resolution of an impact event. We will also define the net impulse over a subinterval of an impact resolution:
(16) 
Implicit in Routh’s method is an assumption that the terminal condition in step 1) will be eventually be reached by any valid choice of increment on ; if it is possible to get “stuck” with , then Routh’s method would be illdefined and not predict a post impact state. It is easy to see that this does not happen in the frictionless case, as has constant positive derivative . The frictional case requires more careful treatment. Intuitively, the added effect of the frictional impulse will be to dissipate kinetic energy quickly. One may conclude that termination happens eventually as zero velocity is a valid postimpact state:
Proof.
Let be a matrix with columns that constitute an orthogonal basis of . By equivalence of norms there exists such that
(17) 
Pick . Let . Assume for .
(18)  
(19)  
(20)  
(21) 
on and thus . Therefore . ∎
The implications of Lemma IIC is that a priori, one can determine a proportional to the preimpact velocity such that any solution to the differential inclusion (12) on can be used to construct the postimpact velocity . We will see, however, that the extension of this methodology to multiple concurrent impacts is nontrivial, and that the physicals systems associated with these models often exhibit a high degree of indeterminacy.
Iii Simultaneous Impact Model
Iiia Motivating Examples
We include, as motivation, two common robotics examples that exhibit simultaneous impacts: one related to legged locomotion and the other to manipulation. Both examples, depending on initial conditions and model properties, can exhibit nonuniqueness. Before describing our model in full detail, we present these examples by considering the outcome of applying Routh’s method to a single contact at a time.
IiiA1 Rimless Wheel
The rimless wheel is a commonly used description of simple robotic walking [8]. Here, we will analyze the case where two feet contact the ground. This can occur if the robot were to fall on two feet, simultaneously, or when one foot is in sustained ground contact and the other impacts the ground. Note that this example is not limited to a legged robot with locked hip and knee joints; see [25] for a thorough analysis of similar legged examples.
For a simple example, illustrated in Fig. 3, we assume that both feet strike the ground vertically, with friction sufficient to sustain sticking. In this case, existing simulation schemes ([32, 1] and others) predict that equal impulses are generated on both feet, brining the robot to rest immediately. However, as illustrated in the figure, if the contacts are sequenced one at a time, other postimpact states are possible where one leg separates from the ground. For other configurations of this problem, nonunique solutions exist spanning sticking, sliding, and separation all for a single initial condition.
IiiA2 Nonprehensile Pushing
In this second example, motivated by nonprehensile pushing of an object, we take a boxlike object (Fig. 4) to have one corner sliding along a surface before impacting a frictionless second surface. Here, the impact on the right wall causes a secondary, frictional impact against the lower wall.
If the first impact is taken to termination before activating the contact on the right wall, the solution in Fig. 3(b) is discovered. Here, the bottom contact is separating and the right contact is sliding upward. Instead, if the impact switches prior to termination, shown in Fig. 3(c), a slightly different solution emerges. This example illustrates that, in simple cases, reminiscent of common robotics applications, subtly different nonunique solutions can emerge from multiple contacts.
IiiB Model Construction
As postimpact velocity is sensitive to the ordering of individual impact resolutions, if we would like to predict as many reasonable postimpact velocities as possible, we must use as relaxed of a notion of impact resolutions as possible. A similar model, without theoretical results or a detailed understanding, was proposed by Posa et al. [24] where it proved useful for stability analysis of robots undergoing simultaneous impact. We consider a formulation in which at any given instant during the resolution process, the impacts are allowed to concurrently resolve an any relative rate:

Monotonically increase the normal impulse on each nonseparating contact at rate such that
(22) 
Increment the tangential impulse for each frictional at rate such that

Terminate when .
We can understand the constraint (22) on as choosing a net force that comes from a convex combination of the forces that Routh’s method might select for any of the individual contacts . As in the single contact case, we might instead think of the selection of a as picking an element of a set of admissible values for . As before, we construct a u.s.c. differential inclusion to capture this behavior:
(23) 
(24) 
We also define a “total impulse” over an interval as identically to (16). Similar to (14), one can extract from a solution such that
(25) 
IiiC Properties
The construction of (24) appears quite similar to that of the single contact systems (13) and (15)—in fact, it is equivalent when is a singleton. We now detail several properties of the multicontact system that are useful for the analysis of its solution set.
IiiC1 Existence and Closure
For any , is trivially closed, uniformly bounded, and convex as it is constructed from the convex hull of a set of bounded vectors. Therefore, in light of Proposition 5, we obtain the following: For all , velocities , and compact intervals , and are nonempty and closed under uniform convergence.
IiiC2 Homogeneity
As each and are conic, and therefore are positively homogeneous in . That is to say, , . Positive homogeneity induces a similar property on : [Solution Homogeneity] For all , , and compact intervals , if , .
IiiC3 Equivalent Minimal Coordinate Systems
In light of (25), for a set of contacts , we have that for all solutions . As in Lemma IIC, it will be useful to analyze the the evolution of a minimalcoordinate representation of ’s projection onto . Let be a matrix with columns that constitute an orthogonal basis of . Therefore is an orthogonal projector onto and
(26)  
(27) 
Therefore, by defining a new set of contacts with equal size to such that , we have that
(28)  
(29)  
(30) 
We denote the collection of contact sets of this size that comply with the full rank condition (30) as
(31) 
Note that does not imply that all of the rows of are linearly independent; for systems with many contacts, will have many more rows than columns, and in this case implies that for every element of the velocity , there exists a contact that can perturb it.
IiiC4 Energy Dissipation
A basic behvaior of inelastic impacts is that they dissipate kinetic energy . We now examine the dissipative properties of the model, which function both as a sanity check on its physical realism and as a device to prove critical theoretical properties. On inspection of (7)–(8) and (14), must be nonincreasing, and furthermore, unless is stationary, it will decrease by a nonzero amount: [Dissipation] Let and . Then is nonincreasing. Let , and let be a compact interval. If and , constant implies constant.
Proof.
See Appendix A. ∎
One might then wonder if is strictly decreasing on . This would obviously be untrue if there existed a with as would be an element of for all compact . We will denote the collection of contacts that do not have this property as
(32) 
Critically, covers most situations in robotics, including grasping and locomotion, with the notable exception being jamming against an immovable surface. Sumsofsquares programming [21], a form of convex optimization, can be used to certify membership in .
Iv Finite Time Termination
While solutions to the underlying differential inclusion are guaranteed to exist in the multicontact model, we have yet to show that they are guaranteed to terminate, as in Routh’s singlecontact method. There have been termination proofs for other simultaneous impact models (e.g. [30]), but only under nontrivial assumptions (no friction; limited number of contacts; smaller space of solutions; etc.). We show what we understand to be the most permissive guarantee of termination: For any preimpact velocity for a contact set , The differential inclusion (24) will resolve the impact by some time proportional to .
We will show that this claim is true as a consequence of kinetic energy decreasing fast enough to force termination—a significant expansion of the claim in Theorem IIIC4. Even though must always decrease, Theorem IIIC4 does not forbid . In fact, it is not possible to create an instantaneous bound . For example, consider 2 frictionless, axisaligned contacts such that . For every , we can pick a velocity
(33) 
(34) 
and arrive at . However as we take , converges to to boundary of and thus will only be ably to sustain a small for a small amount of time before terminating the impact. It remains possible that the aggregate energy dissipation over an interval of fixed nonzero length can be bounded away from zero. We establish a rigorous characterization of this quality by defining dissipativity:
[dissipativity] For a positive definite function , the system is said to be dissipative, if for all , for all s.t. , if , . Denote the collection of contact sets with this property as
(35) 
Unsurprisingly, if on and decreases at a known negative rate, we can show that any trajectory of the multicontact system will exit at a time linearly bounded in : [Bounded Exit] Let be positive definite and let be dissipative. Then , , .
Proof.
It is trivial to show that any contact set that complies with the strong assumption of dissipativity is an element of , as otherwise and could be constant (i.e. ). Far more useful is that we will show Theorem IV arises from the converse: that every exhibits dissipativity. This is particularly surprising for systems with not full rank, as could be large, yet the projection of onto could be arbitrary small, which in turn allows to be small. We observe that the rank of does not effect whether or not , as all solutions will fall into two categories: either is large, or the related minimal coordinate system will exit very quickly: ,
Proof.
See Appendix B. ∎
Finally, we prove the primary claim of this work. Intuitively, if there exists that is not dissipative, then one could construct a sequence of convergent solutions to that dissipate arbitrarily small amounts of energy. Therefore their limit, also a solution to , as the solution set is closed, dissipates no energy—leading to a contradiction with Theorem IIIC4. This argument will be used in an inductive manner, incrementing the size of the contact sets: [Dissipation Inductive Step] Assume for all with or and . Then .
Proof.
Suppose not. Then by Theorem IV, there is a set of contacts , , and a corresponding sequence of solutions , , all starting with velocity magnitude () and never exiting . We must also have that each dissipates less energy than the last: . As is uniformly bounded, are equicontinuous and bounded in . By Theorem IIA, we may assume that such that . Therefore for all and by Theorem IIIC4 is constant. As , by Theorem IIIC4, is not an element of (i.e., ). As is full rank, . Let be the corresponding force vector for each .
Case 1: One contact has strictly deactivated (, ). But then as , by taking a subsequence starting from sufficiently high we may assume that never activates (), and therefore at least one of the other contacts is always active (). But then only the forces from determine , and thus . As the shrinking the contact set shrinks the set of possible forces to apply (), and contains only contacts. But then, by assumption, for some , is dissipative. But . Contradiction!
Case 2: At least one contact always slides . Let be the set of contacts that slide for velocity . Then as is upper semicontinuous, ; that is to say, convergence of the velocity to implies convergence of the direction of sliding on each contact in . Therefore WLOG by taking a subsequence starting from sufficiently high we may assume sufficiently close to and associated new contacts such that
(36)  
(37)  
(38) 
for . Denote and . (37) and (38) in conjunction imply that, for velocities that stay close to , each sliding frictional contact pushes mostly in one direction. Furthermore, the associated frictional force can be generated by three frictionless contacts tilted away from the sliding direction () which never deactivate (). An illustration for this construction is available in Figure 5. As has strictly fewer frictional contacts than and is not dissipative (), by assumption we must have that . By definition of there must exist some penetrating velocity such that is a permissible net force. We therefore must be able to find individual contact forces with and for each contact such that and . As no combination of the original contacts can create zero net force alone, one of the must strictly activate (). By construction of and and the assumption of Case 2, we have , and thus for each and for each . Thus . But then . Contradiction! ∎
We are now ready to prove the main result of this section.
V Conclusion
Nonunique behavior is a pervasive complexity that is present in both realworld robotic systems and common models capturing frictional impacts between rigid bodies—and thus accurate incorporation of such phenomena is an essential component of robust planning, control, and estimation algorithms. Our model presents a stateoftheart theoretical foundation for the capture of this behavior, because despite the high versatility of allowing impacts to resolve at arbitrary relative rates, it is guaranteed to terminate in finite time under far more modest conditions than shown for previous models.
The logical progression from these theoretical results is to develop a numerical scheme to generate postimpact velocities. Constructing approximate solutions to the differential inclusion poses significant challenges associated with discontinuities in ; tools from timestepping schemes (e.g. [1, 32]) may circumvent this issue. Another strategy is to precompute a formula for the entire postimpact set as a function of . Sumsofsquares programming presents potential for construction of an outer approximation.
Future generalizations of the model include elastic impacts using Poisson restitution; resolution of Painlevé’s Paradox; and a full rigid body dynamics model that has continuous solutions through impact.
a Proof of Theorem IiiC4
Let with nonconstant. Let be the associated vector of force variables. As is continuous, we may select such that , is nonconstant on . Let
Comments
There are no comments yet.