Understanding friction between two sliding surfaces has been a focus of research for many decades. Even under carefully controlled experimental conditions, planar friction manifests a surprising degree of variability, which makes it difficult to model and control. This paper studies the nature of that variability in a robotics context. In particular we focus on the analysis of a recent dataset of sliding frictional interactions by Yu et al.  and Bauza and Rodriguez . Our goal is to explain the nature of the observed variability (e.g., see Fig. 1) and suggest models to account for it.
Our analysis shows that a detailed understanding of the experimental process used to collect the data is key. Indeed data-driven modeling and automatic data collection are widespread techniques in robotics. In contact tasks, which involve complex and weakly observed dynamics, it is difficult to avoid data variability due to sensor noise, experiment bias, and task variations. Understanding the nature of that variability is an important step to build reliable data-driven models.
Section III describes experimental facts observed in Yu et al.’s dataset , including the dependence of the data dynamics with the initial orientation of the pushed object, the pronounced variability of the pushes, and the bias observed in the dataset. In particular we focus on sequences of 2000 repeated pushes, part of the dataset, that lead to outcomes with significant variability, as in Fig. 1. These point to three key questions that we attempt to answer in this paper:
I-a Why is pushing not behaving determiniscally?
Pushing, a task driven by the dynamics of planar friction, is a key manipulation primitive that humans and robots exploit to manipulate objects. The robotics community has developed a series of models to describe its dynamics [3, 4, 5, 6, 7] relying on a deterministic Coulomb friction law, naturally resulting in a determined pushed motion. These models have been widely used to develop control strategies [5, 8, 9, 10]. The accuracy of these models remains in question since the modelling and mechanism of friction is still an open problem [11, 12, 13, 14].
If it is natural to explain friction as a deterministic, albeit weakly observed, process, how do we explain the significant and structured variability of trajectories in Fig. 1? Although it is possible—and maybe practical for robotics—to describe that variability lumped into a stochastic process [2, 15], in this paper we show that it is also possible to explain most of the variability in simple mechanical terms. Section III-B describes the formation of the path variability for four different materials: plywood, abs, delrin, and polyurethane.
I-B What is the origin of the structure in the noise?
Yu et al.  describes the variability of friction in the dataset with respect to changes in the location, direction, time, and speed of a push. These account for a fraction of the necessary discrepancy against simple models of friction. Note however that Fig. 1 shows a marked non-Gaussian complex structure, which is more difficult to explain with mechanical arguments of simple variations of the coefficient of friction. Does the effective friction coefficient follow a complex distribution?
In this paper, we show that it is not necessary to invoke stochastic explanations, and that the multi-modality structure can be explained by a combination of anisotropic friction and bias in the data. Section III-C describes the process by which the bias in the dataset is formed by looking at the data collection as a dynamic process. This leads to stable and unstable pushing directions that affect the distribution of initial orientations of the pushed object.
I-C How can we model friction more accurately?
assume that the contact behavior follows a generalized Coulomb friction law, but that the coefficient of friction is a random variable following a Gaussian distribution. Numerical results show that this friction model can simulate some of the variability. Following the observation that different pushes yield different levels of uncertainty,Bauza and Rodriguez 
There are simpler mechanical explanations which require to account for anisotropic friction and biases in the dataset. In this paper, we propose an anisotropic friction model and simulate the data collection dynamics with it. The simulation can reproduce the bias and the direction convergence that is observed in experiment data and confirms that anisotropic friction is a key source of variability that can explain the multi-modality.
Ii Review of experimental setup
Before starting the analysis, we review the experiment setup used to collect the data.
Ii-a Experiment setup
To enable the collection of a large data-set, Yu et al.  automated a loop of pushing & re-positioning an object, with pushes carefully executed in the initial reference frame of the object. A subset of the experiments yielded continuous runs of 2000 identical pushes over four different surfaces made of plywood, abs, delrin and pu. Each experiment is composed of two phases:
Pushing phase: The robotic arm fitted with a thin cylindrical rod (see Fig. 1), pushes the object along a predetermined trajectory relative to the initial pose of the object. In the set of experiments we analyze in the paper, the pushing direction starts orthogonal to the object edge, as illustrated in Fig. 1. In the particular experiment that leads to the distribution, the pushing distance is 150mm with velocity 20mm/s.
Re-positioning phase: After a pushing phase, the robot drags the object back to the central area of the surface by inserting the thin rod in a ring attached to the object and dragging it. The dragging ring is located off-center, as shown in Fig. 1.
Ii-B Structured variability and notation
As illustrated in Fig. 2, the initial position and orientation of the object at the beginning of each push varies. The re-positioning phase brings the object close but, not exactly, to the center of the surface. This is done for simplicity in automating the experiments, but also because it naturally generates a richer dataset.
For analyzing the results, we denote Initial Object Frame (IOF) as the initial frame of the object before it is pushed. In that frame, the object starts at and ends at .
Iii Analysis of experimental data
In this section we first describe facts observed in the dataset and then show how they play a key role in the formation of structured variability. Finally, we analyze the data collection dynamics to explain biases in initial conditions.
Dependence on initial orientation. Simulators frequently assume that friction is homogeneous and isotropic. In our particular scenario, that would mean that the trajectory of a pushed object in the initial reference frame of a push (IOF) is independent of its initial position and orientation. This turns out to be far from reality.
Figures 3 and 4 illustrate it. Fig. 3 represents the output of experiments the way contact models are used in robotics (all position and direction dependence lumped in one single model). Fig. 4 unwraps the trajectories in the global reference frame which is more informative.
Trajectories starting from different initial orientations bend in different directions, which suggests a marked anisotropy. More significantly, trajectories from different initial orientation in Fig. 4 seldom cross each other, which suggests a stable deterministic friction law between object and the surface.
The common assumption that friction is uniform and isotropic leads to a dynamic system represented by the aggregated plot in the object frame Fig. 3. This shows artificial uncertainty due to ”state compression”. When we decompress the state in Fig. 4, we see that anisotropy is a likely cause. In Section III-B we will analyze this relationship.
Stochasticity vs. determinicity The density plots in Fig. 5, which show the change in orientation after a push as a function of the initial orientation, also suggest an almost deterministic relationship. The plot shows that the dynamics of frictional pushing are quite stable given the initial orientation , even after hundreds of pushes and re-positionings starting from the same initial orientation. Although there is still some noise, the plots suggest a clear functional rather than statistical relationship between them. The spread of data points on abs and delrin is likely due to wear and aging of material after hundreds of pushes. Friction shifts as experiment evolves. The material degradation for plywood and pu is not as severe since the pushes are less concentrated.
Bias in data-set Figure 6
shows histograms for the initial orientation of the object for all four materials. These data, indicating strong bias in distribution, are very far from uniformly distributed and different from each other. Onplywood, we observe one high peak at around and a small one at its opposite direction. In abs we only observe a small peak at around and a high peak at . On delrin we only observe one peak at around . We refer to these as stable directions. On the contrary, polyurethane has a more uniform spectrum with multiple peaks at an almost periodic structure.
An initial clear observation is that the dataset is highly biased toward specific initial orientations. In Section III-C, we discuss in depth how these distributions are formed, and why there are sharp peaks in some materials and not in others.
Iii-B Analysis: Formation of Structured Uncertainty – Data Collection Dynamics
In this part, we describe the structure of the variability/uncertainty of the pushing motions as a combination of three factors: a compressed state representation, bias of data-set, and anisotropic friction.
For doing so, we start by attempting to recreate the histograms of the final poses of the pushed object by:
The result is shown in Fig. 7 which is strikingly similar to the real histograms for all four materials. The only exception is an extra lump in the histogram of abs which is attributed to severe wear of the material after hundreds of pushes in the same area.
For numeric comparison, the root mean square error (RMSE) of prediction with an isotropic and homogeneous Coulomb friction and the fitted ’anisotropic law’ are shown in Table I and Tables II respectively. Table I shows the deviation of all the and in the data-set for a particular push. These numbers measure the variability of the data, and can be understood as the smallest RMSE an isotropic and homogeneous friction law can achieve at predicting it, which indicate the limitation of predictive performance by an isotropic and homogeneous friction law.
Tables II shows the RMSE and percentage errors for the prediction by the fitted ’anisotropic law’. Note that RMSE is the square root of the mean of the squared differences between predicted values and observed values. We see that the anisotropic law with the biased initial conditions yields a much better fit to the data in all prediction aspects. Note also that the prediction error is considerably smaller for materials plywood and pu.
Iii-C Analysis: Formation of Data-set Bias
Recall that the histograms in Fig. 6 show a marked bias in the initial orientation of the object, implying that the experiments happen much more frequently along certain directions. We have seen that bias plays a key role in forming the artificial multi-modal structure in the histogram of final pushed motions . In this section we study how that bias is formed.
To characterize the bias, we need to better understand the data collection dynamics by analyzing the time history of initial orientations for data collection. As illustrated in Fig. 2, the position and initial orientation at push differs from that at push . These are related through the cyclic dynamics of pushing and dragging back.
The curves in Fig. 5 show the change in orientation after the pushing phase. Now we are interested in the change in orientation after the combination of pushing and dragging phases. Let be the function that maps an initial orientation to the next initial orientation.
Figure 8 shows the evolution of as the experiment progresses, without removing the natural winding of the angle. Each material shows a slightly different behavior:
plywood: The curve shows a stair-like growth with a step height of exactly . The flatness of the step corresponds to the stable direction where experiments repeats. A step height of implies that there is only one stable direction for the initial orientation. Signs of wear on plywood are also apparent, since the ”stability” of that direction seems to decrease as the experiment progresses.
abs: The initial orientation gets ”stuck” in one stable direction where experiments repeat for about 150 times. Then it escapes and gets trapped in second stable direction from where it never escapes again. This indicates that at least two stable directions exists for the data collection dynamics on abs.
delrin: The evolution of the initial orientation is quite similar to that of abs, with a possibly less stable orientation.
pu: Unlike all other materials, we do not observe any stable orientation in polyurethane. The cumulative initial orientation is always increasing A zoomed-in view of pu in Fig. 8 shows how the initial orientation changes very regularly versus the pushing count.
It is no surprise that the stable directions in Fig. 8 correspond directly to the histogram peaks in Fig. 6. It is key then to understand how the stable directions are formed. To answer this question, we propose an anisotropic friction model in Section IV and use it to simulate the data collection process in Section V.
Iv Anisotropic friction law
The analysis in Section 5 indicates a deterministic relationship between the initial and final orientation of the object, which suggests a deterministic friction law. Anisotropic friction is a good candidate because of its dependence with the sliding direction.
Coulomb’ friction law on a point contact can be represented as a circle of radius , the limit circle . A natural generalization to include anisotropy is to define an equivalent elliptic limit circle (Fig. 9 (Left)). In the general case, the center of ellipse is offset from the origin and the two principle axes are not parallel to and axis, and are of different magnitude. We note the center of ellipse with and the two principle axes with and , rotated by from and .
We choose this model for its compactness as well as with a physical basis. Take the example of wood. Its texture represents the orientation of organized wood fibers. Although wood is usually viewed as bulk material, the micro-structure of wood fibers can yield anisotropic frictional behavior. If the wood fibers are parallel to each other, it is natural to expect two different coefficients for the two orthogonal directions following the fibers and orthogonal to them. Furthermore, while we usually think of texture as a 2D distribution, under the microscope, its micro-structures are 3D ridges. If these are not symmetric, sliding along opposite directions can produce different friction forces. This will produce a limit ellipse with an offset center.
If we denote to be the coefficient of friction in x and y direction, then and . Thus, the ellipse limit circle of anisotropic friction can be expressed as:
If we denote the non-zero sliding velocity of a point contact with anisotropic friction to be
, the maximum-power inequality, a.k.a. maximum dissipation principle, leads to a friction coefficient vector, where
The friction force is , where is the normal load. Since (2) is nonlinear on and , one essential aspect of this anisotropic friction law is that the friction force is not necessarily co-linear with the sliding velocity . The non-linearity and directional dependence of frictional force from anisotropic friction law complicates the pushed motion of the object and in consequence forms directional preferences of experiment distribution.
One of the key experimental observations in Section III is that there are directional preferences for experiments and biases in distribution of initial condition. As indicated in Section 5, anisotropic friction is a likely source of the viability. In order to reproduce the experimental phenomenon and validate the claim, we conduct numerical simulations of the data collection dynamics. For simplicity, we only carry out the simulation experiments on plywood.
V-a Dynamics of pushing and dragging
Given the generalized coordinates of the object , if we denote and
to be the mass and moment of inertia of the object and denotediag() to be its mass matrix, then the Newton-Euler equations for the object are:
where is the force and torque applied by the pusher and is the frictional load applied by the surface under the object. The interaction force in normal direction between pusher and object is modeled with a penalty function method.
V-B Friction modelling
We model the frictional forces between the pusher, object and supporting surface. If we denote and to be the tangential and normal forces between object and pusher, then the Coulomb friction law is:
As shown in Fig. 9 (right), the contact patch between object and plywood is modeled as sets of rigidly connected point contacts. To be specific, we simulate the face contact as an array of point contacts each subject to the anisotropic friction model in (1). The total frictional wrench is the sum of frictional wrenches from each point acting at the mass center of object. The normal load of each point is assumed to be of weight of the object.
V-C Parameter Identification
To identify the friction parameters for the contact between pusher and object, we analyzed the ratio of tangential and normal forces at the pushing point for one pushing experiment. Fig. 10 (Left) shows that the ratio is when the contact point pair is sliding, while it lies between and when the contact point pair is sticking. Thus Coulomb friction law can describe the friction between pusher and object approximately, with coefficient of friction .
For the contact between object and plywood surface, we identified the parameters of the elliptic limit surface via manual fitting, in Table III. Figure.10 (Right) compares the fitted limit ellipse with measured data.
V-D Simulation Results
We carry out a simulation of the push-and-drag experiments in  for 600 cycles in 6 batches, each consisting of 100 cycles with a different starting orientation: , , , , and .
Numerical results show that, independently of the initial orientation, the orientation of the pushing experiments converges to a stable direction in less than 50 cycles.
Figure 11 shows that, in contrast to an isotropic friction model, an anisotropic friction model generates significant biases in the distribution of initial orientations. Figure 12 shows the same trajectories but in the reference frame of the initial orientation of the object, recreating the plots in Fig. 3. This supports our hypothesis that anistotropic friction introduces sufficient non-linear dynamics in the data collection process to bias the collected dataset. This motivates the need for a more in depth analysis of the automated data collection process as a dynamics system itself.
We also simulate the pushed motion under an anisotropic friction law, and using the experimental distribution of initial orientations. Fig. 13 plots the dependence of with respect to the initial orientation which shows a similar structure when compared with experiment data on plywood. This indicates that anisotropic friction is key to explain the direction dependence of pushed motion.
Finally, we would like to make a reference to recent work by Zhou et al. , where they model the stochasticity in pushed trajectories by sampling the coefficient of Coulomb friction from an interval, and yielding to a similar plot. In this paper we suggest that the ’variance’ in the coefficient of friction is likely due to a rather deterministic but anisotropic friction interaction.
Vi Conclusion and Discussion
The main purpose of this paper is to understand the structured variability manifested in planar pushing dynamics  and bring to light the importance of anisotropic friction. We focus the discussion in two aspects: the uncertainty in pushing experiments and the data collection dynamics.
Anisotropic Friction and Structured Uncertainty. One of the key motivators for this paper is to better understand the nature of the uncertainty in the pushing trajectories in Fig. 3.
We attribute it to three main factors: 1) a compressed representation of the state that projects contact dynamics to the initial reference frame of the object effectively lumping variability in the environment; 2) anisotropy of friction; and 3) sharp bias of the initial orientations of the collected data. The combination of these effects explains most of the structure in the noisy plots in Fig. 3. Remaining variability can be attributed to heterogeneity and aging of the surface.
Standard deterministic models of friction dynamics with isotropic and uniform Coulomb friction imply a dynamic model of pushing that is invariant to the initial orientation of the object. But the detailed analysis indicates that anisotropy plays a key role in bending the pushing trajectories in different directions. We simulate the pushing and dragging dynamics with an anisotropic friction law and show how bias and variability is formed.
More efficient algorithms for identification of anisotropic friction parameters is an interesting topic we would like to exploit in the future. Although identification procedures for either parametric or non-parametric isotropic friction have been proposed, efficiency and accuracy are still open problems.
Micro-scale texture leads to asymmetries in friction. This could be exploited to generate useful directional behavior. We are interested in the problem of controlling anisotropic friction though embedding micro-textured patterns on contact surfaces. This opens the door to engineering friction for the purpose of robotic manipulation and locomotion.
Data collection dynamics
Practical advances in machine learning and data-driven modeling are closely tied with big-data. The availability of large and nicely balanced data-sets is increasingly key to develop high performing systems.
Data collection with real robots in real scenarios, however, is much more challenging than in applications where simulation or computational data is sufficient. We have seen that the dynamics of even a simple pushing automated data collection process can lead to significant biases in the dataset. As the system evolves, aging and wear is also a concern, which turns the dynamics that we are trying to capture a moving target.
In the case of this paper, the preference of certain initial conditions is a result of experiment design and anisotropy in the materials. Bias contributes to shape the variability of the dataset, which, if not dealt with, can result in deterioration of trained models. Hence the importance of paying attention to the data collection dynamics.
Particularly key for automating data collection in robotic manipulation is the availability of a resetting mechanism that can avoid bias. Resetting a simulation experiment is trivial, but resetting the initial conditions of a real robotic task is more challenging. Simple strategies as in , might leave the door open to bias. More complex resetting strategies that attempt to carefully control the initial conditions, might become as difficult to automate as the original problem we are trying to solve. Injecting controlled noise seems necessary.
Similar biases are present in other data collection experiments, which suggest the importance of a dynamic perspective on experimental data collection. We are interested in mathematical tools and mechanisms to track and control biases in data collection.
- Yu et al.  K.-T. Yu, M. Bauza, N. Fazeli, and A. Rodriguez, “More than a Million Ways to Be Pushed: A High-Fidelity Experimental Data Set of Planar Pushing,” in IEEE/RSJ Internatinonal Conference on Intelligent Robots and Systems (IROS), 2016.
- Bauza and Rodriguez  M. Bauza and A. Rodriguez, “A Probabilistic Data-Driven Model for Planar Pushing,” in IEEE International Conference on Robotics and Automation (ICRA), 2017.
- Mason  M. T. Mason, “Mechanics and Planning of Manipulator Pushing Operations,” The International Journal of Robotics Research, vol. 5, no. 3, pp. 53–71, 1986.
- Peshkin and Sanderson  M. A. Peshkin and A. C. Sanderson, “Motion of a pushed, sliding workpiece,” pp. 569–598, 1988.
- Lynch et al.  K. Lynch, H. Maekawa, and K. Tanie, “Manipulation And Active Sensing By Pushing Using Tactile Feedback,” Proceedings of the IEEE/RSJ International Conference on Intelligent Robots and Systems, vol. 1, pp. 416–421, 1992.
- Jia and Erdmann  Y. B. Jia and M. Erdmann, “Pose and motion from contact,” International Journal of Robotics Research, vol. 18, no. 5, pp. 466–490, 1999.
- Liu  H. Liu, “Pushing With A Physics-Based Model,” Ph.D. dissertation, Massachusetts Institute of Technology, 2011.
- Lynch and Mason  K. M. Lynch and M. T. Mason, “Stable Pushing: Mechanics, Controllability, and Planning,” The International Journal of Robotics Research, vol. 15, no. 6, pp. 533–556, 1996.
- Behrens  M. J. Behrens, “Robotic Manipulation by Pushing at a Single Point with Constant Velocity : Modeling and Techniques,” Ph.D. dissertation, University of Technology, Sydney, 2013.
- Hogan and Rodriguez  F. Hogan and A. Rodriguez, “Feedback Control of the Pusher-Slider System: A Story of Hybrid and Underactuated Contact Dynamics,” Workshop on Algorithmic Foundation of Robotics (WAFR), 2016.
- Marone  C. Marone, “The effect of loading rate on static friction and the rate of fault healing during the earthquake cycle,” Nature, vol. 391, pp. 69–72, 1998.
- Mate et al.  C. M. Mate, G. M. McClelland, R. Erlandsson, and S. Chiang, “Atomic-scale friction of a tungsten tip on a graphite surface,” Physical Review Letters, vol. 59, no. 17, pp. 1942–1945, 1987.
- Bylinskii et al.  A. Bylinskii, D. Gangloff, and V. Vuletić, “Tuning friction atom-by-atom in an ion-crystal simulator,” Science, vol. 348, no. 6239, pp. 1115 – 1118, jun 2015.
- Yamashita et al.  F. Yamashita, E. Fukuyama, K. Mizoguchi, S. Takizawa, S. Xu, and H. Kawakata, “Scale dependence of rock friction at high work rate,” Nature, vol. 528, no. 7581, pp. 254–257, dec 2015.
- Zhou et al.  J. Zhou, J. A. Bagnell, and M. T. Mason, “A Fast Stochastic Contact Model for Planar Pushing and Grasping: Theory and Experimental Validation,” in Robotics: Science and Systems (RSS), 2017.
- Meriçli et al.  T. Meriçli, M. Veloso, and H. L. Akın, “Push-manipulation of complex passive mobile objects using experimentally acquired motion models,” Autonomous Robots, vol. 38, no. 3, pp. 317–329, 2014.
- Salganicoff et al.  M. Salganicoff, G. Metta, A. Oddera, and G. Sandini, A vision-based learning method for pushing manipulation. University of Pennsylvania, 1993.
- Lau et al.  M. Lau, J. Mitani, and T. Igarashi, “Automatic learning of pushing strategy for delivery of irregular-shaped objects,” in Robotics and Automation (ICRA), 2011 IEEE International Conference on. IEEE, 2011, pp. 3733–3738.
- Walker and Salisbury  S. Walker and J. K. Salisbury, “Pushing using learned manipulation maps,” Proceedings - IEEE International Conference on Robotics and Automation, pp. 3808–3813, 2008.
- Zhou et al.  J. Zhou, R. Paolini, J. A. Bagnell, and M. T. Mason, “A convex polynomial force-motion model for planar sliding: Identification and application,” in Robotics and Automation (ICRA), 2016 IEEE International Conference on. IEEE, 2016, pp. 372–377.
- Goyal et al.  S. Goyal, A. Ruina, and J. Papadopoulos, “Planar sliding with dry friction Part 1. Limit surface and moment function,” Wear, vol. 143, no. 2, pp. 307–330, 1991.