So far, I’ve mostly been working on rotational mechanics – multibody models where each mass only has one rotational degree of freedom and two state variables (angle and angular velocity). Every now and then, I use some commercially available tools to simulate more complex 2d and 3d models – often being surprised by how slow simulation is. At the same time, I’ve been wondering a lot:
- how I’d actually model these problems by hand – especially when additional boundary conditions and/or contacts are involved.
- why some video games seem to be able to include complex dynamics with detailed contact and collisions involved, simulating in real time (while commercial multibody simulation often takes hours to calculate a few seconds of dynamic behavior).
An important answer to these questions is the topic of impulse-based dynamics. I stumbled upon it in some amazing literature from Prof. Jan Bender (see http://www.interactive-graphics.de/) and so far I’ve been experimenting a bit with the method myself. This post today features a small and practical introduction to modelling boundary conditions in impulse-based mechanics.
The concept of impulse-based modelling differs compared to the force-based approach I am used to. See the following graphic for a simple comparison:
So let’s look at a simple example. Consider the model of a mathematical pendulum: The equations for a simple pendulum can be derived very easily (see the wikipedia article). Things however get a little bit more complicated when we consider a combination of pendulums (starting with a double pendulum). Let’s get started with impulse-based dynamics by modelling an n-pendulum as a set of masses with 2 degrees of freedom and additional constraints as pendulum levers:
As I mentioned at the beginning, I so far was basically clueless on how to model the problem given above (unless using reduced coordinates – which works fine here but gets extremely difficult with more complex models).
The impulse-based approach calculates each timestep with the following rules:
- calculate impulse for external forces
- calculate new velocity + position and check if all boundary conditions are fulfilled. If not: Calculate compensatory impulse (iterate)
- update velocity and position
- do something with your model (if required – for example draw or save current state)
This is not very different from a simple euler integration method. The concept of fulfilling boundary conditions by iteration however allows us the apply the levers for our model as boundary conditions in the iteration section.
My code for simulating a multi-pendulum looks like this (This is from GNU Octave, but it should serve as convenient pseudo-code for any programming language of your choice):
I’ve hidden some of the boundary condition checking just to show the overall scheme – which is actually very, very simple (the actual simulation happens in ~30 lines of code if comments and whitespace are omitted). The collision and boundary conditions checking works like this:
- calculate new positions (from impulse)
- check if any boundary collitions are violated
- if yes: calculate a compensatory impulse and reiterate
In actual code:
The calculation of the compensatory impulse is where the real magic happens – compared to a force-based approach, where the calculation of compensatory forces for boundary conditions can be difficult and complex, the iteration method with the impulse is very simple and straight-forward. Let’s look at the results (here: plots from GNU Octave, converted into an animated gif):
Let’s do a short summary:
The approach works surprisingly well – the pendulum moves more or less in the way I’d expect it to do and it provides a structure flexible and simple enough to work with (for example add more elements). It simulates conveniently fast (considering the fact that the code is not speed-optimized at all). So basically, impulse-based dynamics are an amazing answer to the questions I started with!
The downside with impulse-based dynamics seems to be a lack of decent accuracy control. Multi-pendulums are known for exhibiting chaotic behavior, but even considering this, different simulation parameters show a high influence on results and the numerical damping in my model (so far, I haven’t added any damping in the model itself). See the following plot for a comparison of model parameters:
The first setting that I’ve initially chosen seems to provide the expected behavior (black curves in the plot – these parameters are also used for the animation above). In the animation above, the number of iterations required for the last step are also included and it appears that changing the parameters towards a higher amount of iterations (basically the result of either a smaller tolerance or a larger stepsize) seems to add damping to the model.
I’ll definitely continue to play with impulse-based dynamics. I’m especially looking forward to modelling friction with this approach – here my main driving force is building large virtual towers of boxes and have them collapse.
The code I’ve been using for this example can be found here. It’s supposed to work with both GNU Octave and Matlab.