Planet SoC

The Monad.Reader Issue 12 is now out. It contains my article on physics engine implementation as well as articles from other SoCers, Max Bolingbroke and Neil Mitchell, and wonderful poetry by Wouter Swierstra, the editor.
Thanks to Thomas Schilling and Roman Leshchinskiy for help with the article.

As I wrote before, I added complete support for polyhedra and fixed collision handler.
Also, I made a little optimization to collision detector: bounding spheres. The idea is that testing spheres for collision is faster than testing polyhedra. Also, if we consider sphere with center in body’s center of mass, it’s rotationally invariant.
Then I spent some time trying to compile head GHC and DPH. My (not so pleasant) experience is described here.
When I succeeded, I wrote simple benchmark to evaluate gains from using bounding volumes. Result was disappointing — simple 2-step simulation with 10 bodies took 9 seconds and almost 1Gb of memory to complete. After consulting with my mentors we concluded that vectoriser is not ready yet (delay is due to issues with GHC build system) and that I should postpone benchmarks for a while.

Thanks to my fellow physicist, Oleg Matveychuk, I’ve done with the collision
handler bug I said in the last status report. That was indeed a mistake in the
paper (Collision
Detection and Response for Computer Animation).
One of the collision equation in that paper looks like

I’m back from IMC at Bulgaria with second prize, and ready to work!

While I was away I thought about generic convex polyhedra and actually
started implementing them. I hope to finish that in a few days. This
not only makes engine practically more useful but also makes some
computations more cheap and robust.

Also I discovered some problem with my last collision handler. It
conforms to one paper, and now I doubt that paper is correct. I’ll
consult local physicists in the nearest time.

Finally, I read several papers on broad-phase collision detection. The
one I liked most is Broad-Phase Collision Detection Using Semi-Adjusting
BSP-trees. So I’ll start implementing it as soon as possible.

Last week I’ve implemented full handling of rigid bodies collision. When two
bodies collide, not just they velocities, but also their angular velocities are
modified, as it is in real world.
Now, supposing engine is correct (it needs more testing, though), need to
make it fast (after all, this is what Data Parallel Haskell is about).
One of the obvious optimization is to perform broad-phase collision
detection, that is, to filter out pairs of bodies which certainly cannot collide
and run collision detection algorithm (which is then called narrow-phase) on the
pair of bodies left.
So, I’ll need to study the subject and choose algorithm which is best
suitable for nested data parallelism. And there will be good chance for
that — I’m leaving to participate in IMC.

With help and patches from kpierre and excid, Hpysics’ visulization code now performs double buffering, which makes it much more pleasant for eyes. Also this allowed me to record a demo.

• Ogg Theora version [2.1M]
• YouTube version (but their player doesn’t go well on this video for some reason)

Some time ago I mentioned a problem I’ve run into with QuickCheck: the following code, intended for generation of non-zero random vectors, does not terminate:

nonzero :: Gen Vec
nonzero = do
  v <- arbitrary
  if norm v > 0 then return v else nonzero

Now it’s time to give the answer to this puzzle, which was found by Dmitry Astapov (ADEpt). Test data generators have an implicit size parameter, which is passed down in monadic computation. So if we got zero size, we will always get zero vectors — “arbitrary” vectors with components from interval [0;0]. So, solution may look like

nonzero = do
  v <- arbitrary
  if norm v > 0 then return v else resize 1 nonzero

Well, it was too optimistic to think that V-Clip is working properly. Number of bugs, errors and just stupid mistakes was sufficient to keep me busy till now. (Although they were not only related to my V-Clip implementation, but to simulation algorithm too, and also I found — and fixed in my implementation — an error in V-Clip paper itself — I’ll probably write about this later.) But I have good news :)

{- Here a demo video supposed to be, but unfortunately I can’t get it of acceptable quality, so you’re encouraged to build and run vis.hs. -}

As you can see, V-Clip isn’t only consuming CPU time, but really detects collisions. Although this demo might seem physically realistic, it isn’t — only linear, but not angular velocities are updated during collisions. This is what I’m going to do this week.

Finally, I’ve got V-Clip working. At least now (after one bugfix) it always terminates — quite encouraging for such algorithm!

*Properties> quickCheck prop_VClipTerminates 
OK, passed 100 tests.
32% Done, ("Edge","Edge").
14% Done, ("Edge","Vertex").
11% Done, ("Vertex","Edge").
10% Done, ("Face","Vertex").
8% Penetration, ("Face","Vertex").
8% Done, ("Vertex","Face").
7% Done, ("Vertex","Vertex").
5% Penetration, ("Vertex","Face").
3% Penetration, ("Face","Edge").
2% Penetration, ("Edge","Face").

The next step in testing will be visualizing cubes collisions — and I’m still moving towards rotation!

This week I was implementing Mirtich’s V-Clip algorithm for collision detection.
Also I found several bugs in my functions, which took some time to fix. To check the code I used QuickCheck properties. I expect to finish V-Clip during the first half of the following week, if everthing goes fine.
I also found time to cabalize and haddockize Hpysics. Note that you need haddock 2 to build docs, because earlier versions don’t like parallel arrays syntax.
Also, there was an interesting challenge with QuickCheck. To test some things I need to generate non-zero vectors. I tried the following:

nonzero :: Gen Vec
nonzero = do { v <- arbitrary; if norm v > 0 then v else nonzero }

(i.e. if we got zero vector, just generate another one until we get what we want), but for some reason this code doesn’t terminate. To understand why, I probably should learn more about how QuickCheck’s rng works, but I won’t complain if someone explains this effect :)