Writing a physics engine (Part 1: Introduction)

Several months ago, I was interested in training RL agents in physics-based environments (e.g., robotics). Of course, that meant that I needed an a physics simulator. Not having much of a background in this space, I was immediately overwhelmed with just how many physics engines and solvers¹ there are out there. There’s a whole handful of different engines that are used in graphics, robotics, biomechanics. To list a few, I’ve looked at:

• MuJoCo. Widely popular for robotics and graphics.
• XPBD. Popular for graphics.
• AVBD. State of the art method for rigid body dynamics in graphics.
• OpenSim/Simbody. Widely popular in biomechanics for studying human movement.
• And some of the solvers implemented in Solver2D <- A really great resource!

Honestly, I didn’t really understand the difference between all of these solvers. It wasn’t really clear to me what each solver had to offer over the others… so I decided to implement all of them. After doing so, it was a lot easier to understand what all these solvers are doing. I honestly think that implementing my own engine was the fastest way for me to understand how they work. It took some time to just get started and become familiar with the notation but I think it’s less daunting then it may seem.

Now, this will not be a blog post comparing each engine, but instead a series of blog posts going through the process of how to implement a physics engine. Specifically, we’ll go over how to implement a velocity-based solver (like MuJoCo) that will be able to simulate things like:

If you want to follow along, I won’t assume that you will have any background in rigid body physics or physics-based animation, but having a couple things in your back pocket will be helpful:

• A basic understanding of linear algebra (e.g., linear operators, spatial transforms) and optimization (nonlinear optimization, iterative methods).
• Some physics background: you understand \(F=ma\). If you’re comfortable with \(\frac{d}{dt} (\frac{\partial \mathcal L}{\partial \dot{q}}) - \frac{\partial \mathcal L}{\partial q} = 0\) that’s cool too.
• If you want to follow along and implement a physics engine, I assume you already have contact detection written.

I’ll plan to structure the series of posts in something like the following order (I haven’t written anything yet so maybe this will also hold me accountable):

1. Introduction (you’re reading this now).
2. Rigid body dynamics (no articulated joints)
   1. Rotations and orientations in 3D
   2. Mass matrix and contact Jacobian
   3. Collision resolution
3. Articulated rigid body dynamics
   1. Generalized coordinates and velocities
   2. Redefining the mass matrix and contact Jacobian
4. MuJoCo’s fast constraint solver
5. Integration methods
   1. Implicit damping
   2. Adaptive time-stepping

Whether you’re interested in building simulators for robotics, implementing a physics engine into your video game, or just trying to better understand how these physics solvers work under the hood, hopefully these blog posts will be helpful!