google / brax

Massively parallel rigidbody physics simulation on accelerator hardware.
Apache License 2.0
2.38k stars 257 forks source link

Add gyroscopic term to physics integrator #48

Open o-Oscar opened 3 years ago

o-Oscar commented 3 years ago

Looking at the code in integrator.py, I did not find the right equation of evolution of the angular velocity (ang in the code)

According to Featherstone (and a derivation of mine and my university's textbook), the correct equations of motions should be :

With tau being the external torques, omega the angular velocity of the body, and I the inertia matrix in the local frame.

I guess is does not make much of a difference when most degrees of freedom are constrained, but clean math should be better.

Any plans to include this version of the equations in Brax ? And if no, can you give me some pointers as to how one might go and implement that ?

cdfreeman-google commented 3 years ago

Because we're using "impulsive" torques, these formulas actually end up being the same. But you're totally right that this is the more general expression for the rate of change of the angular velocity. There's a nice little discussion here: https://physics.stackexchange.com/questions/434760/integrating-rigid-body-equations-for-a-game-engine-simulation where someone noticed exactly this issue.

o-Oscar commented 3 years ago

My interpretation of the discussion you mention does not enter in conflict with my previous point.

In the first reply of the discussion, they go through the derivation of the equation that I mentioned, and then discuss its differences with the impulse formulation. The key difference between the two formulation is how the time is evolving.

In the first part of the answer, the time is evolving : the delta t in the following equation is non zero. image between those two different times, the inertia matrix can not be considered constant. Hence the omega x I omega term.

In contrary, in the second part of the answer, they are working with impulses and with instantaneous velocity change.

In the event that the applied torque takes the form of an instantaneous impulse

In this part of the answer, they are working with zero-sized delta t, and going from t0- to t0+. That leads them to say the following :

The space-fixed inertia tensor does not change during this infinitesimal time interval

In brax, what we essentially want to do is :

During the forward physics integration, the correction term I am talking about should be included, because the delta t is non zero. But then, during the constrain solve, with impulse correction, I agree that this term does not appear.

erwincoumans commented 3 years ago

Thanks for the feedback, o-Oscar!

Yes, the gyroscopic term is missing. See also Erin Catto's slides for some background and how to implement it efficiently: https://box2d.org/files/ErinCatto_NumericalMethods_GDC2015.pdf (slide 60 onwards). When implementing the gyroscopic force as an explicit force, it can easily add energy and 'blow up' the simulation, hence an implicit implementation.

See here the implementation of implicit and explicit gyroscopic forces in our Bullet physics engine (C++), it should be straight forward to add to Brax as well (this is the implemention for maximal coordinate rigid bodies): https://github.com/bulletphysics/bullet3/blob/master/src/BulletDynamics/Dynamics/btRigidBody.cpp#L283

Note that for typical gym environments, such as humanoid, ant etc, the velocities and inertias don't cause large gyroscopic forces, so it is not obvious that omitting this term hurts realism significantly. Without the term, you won't get the interesting Dhanibekov effect though: https://www.youtube.com/watch?v=1x5UiwEEvpQ

By the way, I wonder where your minus sign comes from, usually it is a plus. image

o-Oscar commented 3 years ago

Thanks a lot for the response! Ok, so since in the classical RL settings, most bodies have low gyroscopic torques, and due to the added complexity that would arise in order to prevent divergence, It does make sense to not include this term in the simulation. Thanks for clarifying that!

cdfreeman-google commented 3 years ago

Yep, you're totally right o-Oscar. Managed to confuse the heck out of myself, because I know I've measured conservation laws for our codebase, and everything seemed to check out. Turns out, we've only ever simulated primitives that have diagonal inertia tensors, so this term has always vanished.