Model Identification and Control of a Low-Cost Wheeled Mobile Robot Using Differentiable Physics

With the availability of affordable 3D printers and micro-controllers such as Arduino [1], Beaglebone Black [2], and Raspberry Pi [3], light-weight high-performance computing platforms such as Intel’s Next Unit of Computing (NUC) [4] and NVIDIA’s Jetson Nano [5], and programmable RGB-D cameras, such as Intel’s Realsense [6], there is renewed interest in building low-cost robots for various tasks [7]. Motivated by these developments, the long-term goal of the present work is to develop affordable mobile robots that can be easily assembled using off-the-shelf components and 3D-printed parts. The assembled affordable robots will be used for exploration and scene understanding in unstructured environments. They can also be augmented with affordable robotic arms and hands for manipulating objects. Our ultimate objective is to remove the economic barrier to entry that has so far limited research in robotics to a relatively small number of groups that can afford expensive robot hardware. To that end, we have designed our own wheeled mobile robot for exploration and scene understanding, illustrated in Figure 1. The first contribution of this work is then the hardware and software design of the proposed robot.

High-end robots can easily be controlled using software tools provided by manufacturers. Their physical properties, such as inertial and frictional parameters, are also precisely measured, which eliminates the need for further calibrations. Affordable robots assembled and fabricated in-house are significantly more difficult to control due to uncertainties in the manufacturing process, which result in differences in size, weight and inertia according to the manufacturing technique. Due to this uncertainty, hand-crafting precise and shared models for these robots is challenging. For example, the wheels of the robot in Figure 1 cannot be precisely modeled manually because of their complex structure and uncertain material properties. Moreover, the frictions between the wheels and the terrain vary largely when the robot is deployed on an unknown non-uniform terrain. Statistical learning tools such as Gaussian processes and neural networks have been largely used in the literature to deal with this uncertainty and to learn dynamic and kinematic models directly from data. While such methods have the advantage of being less brittle than classical analytical models, they typically require large amounts of training data collected from each individual robot and for every type of terrain.

In this work, we propose a hybrid data-driven approach that combines the versatility of machine learning techniques with the data-efficiency of physics models derived from first principles. The main component of the proposed approach is a self-tuning differentiable physics simulator of the designed robot. The proposed simulator takes as inputs the robot’s pose and generalized velocity at a given time, a sequence of control signals, and returns a trajectory of predicted future poses and velocities. After executing the sequence of controls on the real robot, the resulting ground-truth trajectory of the robot is recorded and systematically compared to the predicted one. The difference, known as the reality gap, between the predicted and the ground-truth trajectories is then used to automatically identify the unknown coefficients of friction between each of the robot’s wheels and the present terrain. Since the identification process must happen on the fly and in real time, black-box optimization tools cannot be effectively used. Instead, we show how to compute analytically the derivatives of the reality gap with respect to each unknown coefficient of friction, and how to use these derivatives to identify the coefficients by following the gradient-descent technique. A key novelty of our approach is the integration of a differentiable forward kinematic model with a neural-network dynamic model. The kinematic model represents the part of the system that can be modeled analytically in a relatively easy manner, while the dynamics neural network is used for modeling the more complex relation between the frictions and the velocities. But since both parts of the system are differentiable, the gradient of the simulator’s output is back-propagated all the way to the coefficients of friction and used to update them.

The time and data efficiency of the proposed technique are demonstrated through two series of experiments that we have performed with the robot illustrated in Figure 1. The first set of experiments consists in executing different control signals with the robot, and recording the resulting trajectories. Our technique is then used to identify the friction coefficients of each individual wheel, and to predict future trajectories accordingly. The second set of experiments consists in using the identified parameters in a model-predictive control loop to select control signals that allow the robot to track predefined trajectories. The proposed gradient-based technique is shown to be more efficient computationally than black-box optimization methods, and more accurate than a neural network trained using the same small amount of data.

The problem of learning dynamic and kinematic models of skid-steered robots has been explored in several past works. Vehicle model identification by integrated prediction error minimization was proposed in [8]. Rather than calibrate the system differential equation directly for unknown parameters, the approach proposed in [8] calibrates its first integral. However, the dynamical model of the robot is approximated linearly using numerical first-order derivatives, in contrast with our approach that computes the gradient of the error function analytically. A relatively similar approach was used in [9] for calibrating a kinematic wheel-ground contact model for slip prediction. The present work builds on the kinematic model for feedback control of an omnidirectional wheeled mobile robot proposed in [10].

A learning-based Model Predictive Control (MPC) was used in [11] to control a mobile robot in challenging outdoor environments. The proposed model uses a simple a priori vehicle model and a learned disturbance model, which is modeled as Gaussian Process and learned from local data. An MPC technique was also applied to the autonomous racing problem in simulation in [12]. The proposed system identification technique consists in decomposing the dynamics as the sum of a known deterministic function and noisy residual that is learned from data by using the least mean square technique. The approach proposed in the present work shares some similarities with the approach presented in [13], wherein the dynamics equations of motion are used to analytically compute the mass of a robotic manipulator. A neural network was also used in a prior work to calibrate the wheel–terrain interaction frictional term of skid-steered dynamic model [14]. An online estimation method that identifies key terrain parameters using on-board robot sensors is presented in [15]. A simplified form of classical terramechanics equations was used along with a linear-least squares method for terrain parameters estimation. Unlike in the present work that considers a full body simulation, [15] considered only a model of a rigid wheel on deformable terrains. A dynamic model is also presented in [16] for omnidirectional wheeled mobile robots, including surface slip. However, the friction coefficients in [16] were experimentally measured, unlike in the present work where the friction terms are automatically tuned from data by using the gradient of the distance between simulated trajectories and the observed ones.

Classical system identification builds a dynamics model by minimizing the difference between the model’s output signals and real-world response data for the same input signals [17, 18]. Parametric rigid body dynamics models have also been combined with non-parametric model learning for approximating the inverse dynamics of a robot [19].

There has been a recent surge of interest in developing natively differentiable physics engines. For example, [20] used the Theano framework to develop a physics engine that can be used to differentiate control parameters in robotics applications. The same framework can be altered to differentiate model parameters instead. This engine is implemented for both CPU and GPU, and it has been shown how such an engine speeds up the optimization process for finding optimal policies. A combination of a learned and a differentiable simulator was used to predict action effects on planar objects [21]. A differentiable framework for integrating fluid dynamics with deep networks was also used to learn fluid parameters from data, perform liquid control tasks, and learn policies to manipulate liquids [22]. Differentiable physics simulations were also used for manipulation planning and tool use [23]. Recently, it has been observed that a standard physical simulation formulated as a Linear Complementary Problem (LCP) is also differentiable and can be implemented in PyTorch [24]. In [25], a differentiable contact model was used to allow for optimization of several locomotion and manipulation tasks. The present work is another step toward the adoption of self-tuning differentiable physics engines as both a data-efficient and time-efficient tool for learning and control in robotics.

The design of our mobile robot consists of the hardware assembly, and the software drivers that provide control signals for actuation. These components are described below.

Consider a simplified cylindrical model for the wheel. Let $J$ be the scalar component of its $3\times 3$ (diagonal) inertia tensor matrix about the rotational axis of the motor shaft (we discard the components in the plane orthogonal to this axis, as the wheel is not allowed to rotate in this plane). Then, its equation of motion can be written as:

where $T_{s}$ is the motor stall torque, $\omega$ is the wheel’s angular velocity, $\omega_{s}$ is the desired angular velocity (specified by the PWM signal), $\mu$ is the coefficient of friction, $g$ is the acceleration due to gravity, $R$ is the wheel radius, and $M$ is the mass of the robot. We assume that the weight of the robot is balanced uniformly by all four wheels (thus, the term $Mg/4$). The first term on the right hand side was derived in [27]. Equation (1) can be analytically integrated to obtain:

Since we are interested in large time spans for mobile robot navigation, only the steady state terms in equation (2) are important. Thus, we discard the transient terms to arrive at:

which gives us an expression for the wheel’s angular velocity as a function of the ground friction forces. We assume that the friction coefficient $\mu_{j}$ for each wheel $j$ is different. The chassis of our mobile robot is similar to that of the Uranus robot, that was developed in the Robotics Institute at Carnegie Mellon University [10]. The following expression was derived in [10] for the linear and angular velocity of the Uranus robot, using the wheel angular velocities as input:

where $l_{ab}=l_{a}+l_{b}$ and $l_{a},l_{b}$ are as defined in Figure 3, $v_{x},v_{y}$ are the linear velocity components of the robot along the $X$ and $Y$ axes, and $\omega_{z}$ is the angular velocity of the robot about the $Z$ axis. Equations (3), (4) define a motion model for our mobile robot including the effects of friction.

Let $(x^{i},y^{i},\theta^{i})$ be the generalized position of the mobile robot at time $t^{i}$, and $(v_{x}^{i},v_{y}^{i},\omega_{z}^{i})$ be its generalized velocity. Then, its predicted state at time $t^{i+1}$ can be computed as:

where $\Delta t=t^{i+1}-t^{i}$. The loss function computes the simulation-reality gap, defined as the divergence of a predicted robot’s trajectory from an observed ground-truth one. The loss is defined as follows,

where $(x^{k}_{gt},y^{k}_{gt})$ are the ground-truth position values at time-step $k$, estimated using Apriltags [28] and an overhead camera. $(x^{k},y^{k})$ are the positions predicted in simulation by using Equation 5 and the same sequence of control signals as the ones provided to the real robot. The video is recorded using the overhead camera at $60$fps, whereas the PWM signal is sent to the motor at a 3-4$\times$ higher frequency. Thus, there are several simulation time steps between two consecutive ground truth position values.

In practice, we observed that the predicted state can often “lag behind” the ground truth values, even when the robot is exactly following the overall path, leading to a high loss value. Thus, we fit a spline curve to the ground truth values and compute another loss function that uses the point $(x^{k}_{sp},y^{k}_{sp})$ closest to this curve from the predicted state at time-step $k$:

Note that the loss $L_{sp}$ alone is not sufficient, as it does not penalize the simulated robot for not moving at all from its starting position. Thus, we use a weighted linear combination of $L_{gt}$ and $L_{sp}$ as our actual loss function, as given below:

where we set $w_{1}=0.8,w_{2}=0.2,M=4,g=9.8,r=0.03$, and $T_{s}=0.6$. The pseudocode for loss computation is shown in Algorithm 1. Vector quantities are shown in bold. $B$ is defined in Equation 4, and $\mathbf{\mu}$ is a vector of the friction coefficients of the different wheels. Line 4 uses component-wise vector multiplication. $T$ is the total number of time-steps, and $t^{i}$ is the time when the PWM signal was applied to the motor at time-step $i$. The function IsGroundTruthSample checks if a ground truth Apriltag estimate exists at time $t^{i}$, and if so, then the function GroundTruthIndex returns the index of that estimate.

To minimize the loss function in equation (8) with respect to unknown parameters $\mu_{j}$, corresponding to the friction values for the four wheels, we derive analytical expressions for the gradient of the loss with respect to each variable $\mu_{j}$:

Let $(\Delta x^{k}_{gt},\Delta y^{k}_{gt})=(x^{k}-x^{k}_{gt},y^{k}-y^{k}_{gt})$ and $d^{k}$ denote the length $\left(|\Delta x^{k}|^{2}+|\Delta y^{k}|^{2}\right)^{1/2}$. Then the second term in equation (9) can be expanded using the chain rule as follows:

The derivatives in equation (10) can be computed using equations (4), (5) at the higher frequency of the PWM signal sent to the motors, used for predicting the next state, as:

where $b_{rs}$ is the $(r,s)$ entry in the matrix $B$, as defined in equation (4). The derivative $\partial\omega^{i}/\partial\mu_{j}$ can be computed using equation (3). The expression for $\partial L_{sp}/\partial\mu_{j}$ in equation (9) can be derived similarly. Note that, strictly speaking, the closest point $(x^{k}_{sp},y^{k}_{sp})$ on the spline curve is a function of the point $(x^{k},y^{k})$. However, we have empirically found that estimating its derivative is not necessary and can be ignored, when computing the term $\partial L_{sp}/\partial\mu_{j}$, for $j\in\{1\ldots 4\}$.

As shown in Section VII, using the gradients derived in Section VI to minimize the loss function in equation (8) yields friction parameters that give good agreements with the experimental trajectory that is estimated using Apriltags. However, we observed that the computed friction parameters can differ in values for two different trajectories (with different control signals). This implies that the friction $\mu_{j}$ for each wheel is not a constant, but a function $\mu_{j}(\omega_{s})$ of the applied control signal. Thus, we first generate a sequence of trajectories with fixed control signals, and estimate friction parameters for each of them by separately minimizing the loss function in equation (8) using gradient descent. We then train a small neural network with 4 input nodes, 16 hidden nodes, and 4 output nodes. The input to the neural network are the applied control signals to the wheels, and the output are the friction parameters estimated using gradient descent. As shown in Section VII, a sequence of only $8$ input trajectories is enough to obtain reasonable predictions from the neural network, and allowed us to autonomously drive the mobile robot along an “eight curve” with high precision.

The pseudocode for gradient computation is shown in Algorithm 2. Matrix $J$ stores the accumulated gradients for the generalized position, computed in equations (11)-(13). The notation $J[\{0,1\},:]$ refers to the first two rows of $J$.

We remotely controlled the robot for 6-8 seconds and collected 8 different trajectories, with 2 samples for each trajectory that were estimated using Apriltags [28], as shown in Figure 5 in red and blue. The control signals applied to the wheels were constant for the entirety of each trajectory, but different per trajectory. This was done for simplicity, as changing the control signals changes the direction of motion, and the robot may have driven out of the field of view. A video of these experiments is attached to the present paper as a supplementary material.

We presented the design of a low-cost mobile robot, as well as an analytical model that closely describes real-world recorded trajectories. To estimate unknown friction coefficients, we designed a hybrid approach that combines the data efficiency of differentiable physics engines with the flexibility of data-driven neural networks. Our proposed method is computationally efficient and more accurate than existing methods. Using our learned model, we also showed the robot drive autonomously along pre-specified paths.

In the future, we would like to improve our analytical model to also account for control signals that are not powerful enough to induce rotation of the wheels. In such cases, the wheel is not completely static and can still turn during robot motion, due to inertia. This causes the robot to change orientations in a manner that cannot be predicted by our current model. Accounting for such control signals would allow for more versatile autonomous control of our robot.