So I have been lately somewhat involved in modeling Coulombic friction in my (perhaps game-) physics simulations project, as you can see from my post on writing test code and the first O’Caml rant.

Just today I spent some time in augmenting the test-case for the static friction case, and it occurred to me, that you, my loyal reader and supporter, might be interested in what is involved with the math and code when implementing these kind of models, so here we go again:

First of all, let us set the background: We are interested in this case friction between two bodies connected with some kind of constraint – let this constraint be realized by the equation: The function could be almost anything involving the coordinates – for example in my test-case has , which just tells that whatever the rest of the dynamics are, the point-body described by the coordinates should always lie on the negative sine-curve.

Then basic Lagrangian mechanics tells us that there exists a quantity , which ensures that when one applies the force(s) to the body (bodies), then the constraint can be satisfied. This force is the contact or normal force of the constraint and it has the special property that it **doesn’t perform any physical work**, because it is always orthogonal to the movement-vector of the body.

Then the typical, and quite accurate model, for friction force caused by the constraint is that the movement of the body (/ the relative movement of the bodies) is slowed down by a force that is linearly proportional in magnitude to the normal force of the constraint . Then of course the magnitude of the constraint force is and the proportionality constant of the friction force is slightly higher in case the body is not (bodies are not) moving.

Until now, I have mostly concentrating on the dynamical friction, where the friction force is completely opposite to the velocity vector of the bodies :

where is the coefficient of the friction.

Now we get to the actual predicament: How to do this in practice? The force is proportional to another force, but alas, it is the absolute value that it is proportional to! We would get perfect model for linear dependence – meaning something like , but we do have linear dependency on the velocity – *almost*: it is linearly dependent on the direction of the velocity.

In order to study this, I devised three models how to numerically solve them:

1) Make the linear dependency on the velocity linear and use other values from previous known state – this results effectively in:

, where is the previous velocity, is the previous Lagrange multiplier and are the coordinates from previous state. The nice thing about this model is that it is really simple and computationally inexpensive to implement, but it has the weird property that the velocity used comes from two different points in time, but according to my tests the model is quite ok (except for really smally velocities).

2) Make the force linearly dependent on the Lagrange multiplier ,

Here of course is always true, except when the constraint force changes directions between two time-steps, but in these situations normally anyways the bodies can drift apart (in case we have a normal contact-constraint ) and the model works nicely as it is completely up-to-date with the constraint force, but alas it is somewhat more complex to compute.

The third option, that I haven’t tried out yet, is to just throw away all quantities of current step and just throw the entire force on the right-hand-side of the dynamics equation, making it basically a constant force over each time-step and this just might be completely sufficient for my purposes but all this still requires some work and testing, as the dynamics model I’m using is somewhat novel. Hopefully I can share more about that later on. For now, I bid you good night.

Btw, I’m planning on modeling the rest-friction case as just conditional constraint in the tangent-space of the original constraint. More on that to come. π