Comparison to SD/FAST
Data from: MuJoCo: A physics engine for model-based control
E. Todorov, T. Erez, and Y. Tassa (2012). International Conference on Intelligent Robots and Systems
SD/FAST is a cross-compiler that generates model-specific C code for simulating multi-joint dynamics. Even though it has not been updated in nearly 20 years, it has remained the gold standard in many ways. This is partly because model-specific C code is hard to compete with, and partly because at the time people used to calculate the exact number of floating point operations and optimize numerical performance aggressively. This is why when we began working on MuJoCo we did not aim to outperform SD/FAST, but rather to come close to it in terms of multi-joint dynamics while adding contact dynamics and other important features.
We compared MuJoCo and SD/FAST on four model systems: isolated bodies, a chain, a hand, and a humanoid model. To our pleasant surprise, the outcome was a tie: MuJoCo was faster on some tests while SD/FAST was faster on others. The table below shows the number of evaluations per second, rounded to 1000, on one core of a now-outdated Intel processor (average over X9650, X5570, i7 940, X5860).
|Model (DOF)||SD/FAST eval/s||MuJoCo eval/s|
Each evaluation was performed at a randomly-generated state (same for both engines), and included computation of the joint-space inertia matrix and the vector of Coriolis, centripetal and gravitational forces. Once these quantities are available, physics simulation comes down to solving a linear system. Note that we did not use the O(N) solver in SD/FAST, because that solver does not compute the inertia matrix which is essential for simulating contact dynamics. Instead we used the solver which implements similar algorithms to those used in MuJoCo.
Comparison to Bullet, Havoc, ODE and PhysX
Data from: Simulation tools for model-based robotics: Comparison of Bullet, Havok, MuJoCo, ODE and PhysX
T. Erez, Y. Tassa and E. Todorov (2015). International Conference on Robotics and Automation
Bullet, Havoc, ODE and PhysX are among the most popular engines used today for gaming, animation and multi-body dynamics simulation. They come from the tradition of the now-retired MathEngine, and represent the system in Cartesian coordinates (with 6 DOFs per body) while joints are represented as equality constraints imposed numerically. This is natural for systems involving many floating bodies. But for robotic and biomechanical systems where the DOF-to-body ratio is typically close to 1, a joint coordinate representation (as in MuJoCo and SD/FAST) is more efficient computationally, and also more accurate - because the joints cannot dislocate and thus the engine does not have to fix constraint violations at each step.
Our empirical results are consistent with the above reasoning. For multi-joint systems relevant to robotics and biomechanics MuJoCo is the fastest, while for systems composed of many floating bodies MuJoCo is the slowest and ODE is the fastest. Simulation speed however is not the most relevant measure of performance. Speed can always be increased by using a larger timestep, but this comes at the expense of accuracy and stability. So the more relevant measure is the speed-accuracy trade-off of each engine. We define accuracy by running the simulation with a very small timestep (0.016 ms) where all engines should be very accurate, and then measuring the deviations from the engine-specific reference as the timestep increases. We find that even for systems that are better suited for gaming engines, MuJoCo is the fastest at a fixed accuracy level and the most accurate at a fixed speed level. This is because it can use much larger timesteps without loosing accuracy.
Note that the above measure of accuracy is with respect to an engine-specific reference and not an analytical benchmark. This is because we want to have a measure that can be applied to complex systems where analytical results are not available. As a result, an engine can simulate "physics" that deviate from reality without being penalized. For example, both Havoc and PhysX ignore Coriolis forces. But since they do that for all timesteps, our accuracy measure (which really reflects self-consistency) does not penalize them.
The figure below shows four test systems we simulated in each engine (left column), raw speed (middle column), and speed-accuracy trade-off (right column). MuJoCo has two flavors depending on the integrator being used: semi-implicit Euler vs. 4th-order Runge Kutta. For one of the test systems (planar chain) we were also able to use the new version of Bullet (labeled BulletMB) which operates in joint coordinates; for the other systems it seemed to need more work before it can be used reliably. PhysX also has a version aimed at simulating multi-joint systems, but it does not use joint coordinates; instead it relies on the underlying Cartesian representation and has similar behavior as the standard version.
Some of the speed-accuracy curves are incomplete because when a given engine diverges at a given timestep we cannot obtain a measurement. The hand grasping an object was the most challenging test system, and is also very relevant to robotics. To measure the maximum timestep at which each engine could simulate this system without dropping the object, we started at timestep 1/64 ms and increased it by factor of 2 until the object was dropped. The table below shows the last successful timestep for each engine. Note the order-of-magnitude difference between MuJoCo and the best-performing gaming engine in this test (PhysX).
Scaling to many bodies
The first three tests above involved systems that are challenging to simulate because of large mass ratios, elaborate kinematic structures and contact configurations that propagate forces in complex ways. The numerical challenges are different for systems with large numbers of floating bodies. They do not require the accuracy needed for robotic systems and do not go unstable as easily, but in order to manage the large number of bodies it is essential to use algorithms with linear scaling. Such algorithms are relatively straightforward to implement in Cartesian coordinates, and considerably more difficult to implement in joint coordinates. The latter requires exploiting the sparsity of the joint-space inertia matrix and the contact Jacobian, and implementing the necessary indexing machinery.
In addition to the default dense solver used in the above tests, MuJoCo has a sparse solver designed for such situations. The figure below illustrates a test where we drop N capsules, with N between 50 and 250, and measure the CPU time for running a 10 sec simulation at 10 ms timestep. The sparse solver has close-to-linear scaling and is therefore suitable for such systems. The dense solver on the other hand has between quadratic and cubic scaling, and should not be used in this case. Note the different scale on the vertical axis in the two plots. The curve labeled "no contact" corresponds to a simulation where all contacts are disabled and the objects simply fall through the ground. These results are obtain on one core of an Intel i7-3930K processor.