Improve simulation stability when object are far from world origin

[LeshaInc]

During integration and constraint solving, position gets incremented by very small amounts, especially at high substep count. When the objects are far from origin, small increments to position are rounded to zero. To mitigate that, I've added AccumulatedTranslation component, which gets incremented instead of Position. Then, at the end of each substep, a system updates the actual position.

This change seems to be the magic bullet for "far from origin" kinds of problems. But it's not perfect yet.

Here's the cubes example with origin set to (10000, 10000). When using f32, the next float after 10000 is 10000.0009766, which gives precision of roughly 1e-3 .

Here's the before:

before.mp4

Here's the after:

after.mp4

As you can see, it starts very promising, but then explodes. I don't know what causes that, but If I change substeps from 12 to 6 (not part of PR), it doesn't explode:

after_6_substeps.mp4

Also, I've changed contacts to be queried in local space, since that's what parry does internally, only then to transform it back to world space (and back to local in PenetrationConstraint). This improves precision.

Additionally, I've found a few places where values of different orders of magnitude were subtracted (for example when computing relative motion of contact points): I've reordered operands to try and keep relative magnitudes close.

Future work

In the future, perhaps it would be beneficial to update actual position not on every substep, but on every step instead. This way, even on very large substep counts (where dt is very tiny) everything should work fine.

[LeshaInc]

After merging #83 cubes example with origin set to (10000, 10000) still explodes, but after 10 seconds instead of 4.

explosion.mp4

I've tried increasing AABBS like so:

mins -= half_extents * 0.1;
maxs += half_extents * 0.1;

I've waited for a few minutes and so no explosion, and so I think the problem is caused by missed collision pair after some amount of substeps, since the broad phase is run only once.

Do you have any clues regarding this explosive behavior?

[Jondolf]

I'm still getting explosions with AABBs expanded further like that, although I'm not getting the issue when the cubes are initially stacked (we might have different OSs though, I'm on Windows).

In your case I'm guessing all the cubes were marked as sleeping, so they weren't simulated and couldn't explode even if the issue is still there. If I keep moving the cubes frequently, I still get explosions:

2023-07-18.12-50-58.mp4

I don't know what's causing the explosive behaviour, but I think a couple things are implemented slightly incorrectly.

1. I'm pretty sure constraints and solve_vel should use local_r1/local_r2 rotated to match the current pose instead of world_r1/world_r2. Using the original contact points is wrong, since the body is often rotated during the constraint solver. When I previously tried this, I was getting mixed results, so I'm not entirely sure.
2. Dynamic friction is wrong, since it causes heavier bodies to come to a stop earlier. See Dynamic friction fix #52.
3. Static friction also feels slightly wrong, but I'm not entirely sure without more testing. I have found several slightly different implementations and I'm not sure which one is correct.
4. Probably unrelated, but currently the narrow phase is inside the constraint solver instead of a separate step. Once a collision is detected, the constraint is solved immediately, instead of computing all collisions beforehand. The current implementation seems to be more stable, but it would be nice to have a separate narrow phase for parallelism and modularity.

These are probably not directly related, since disabling friction didn't seem to help when I tried a while ago, and the other changes cause mixed results. Using f64 seems to help a lot, so I suppose it could be another precision issue somewhere.

[Jondolf]

As for this PR, it's definitely a massive improvement when using high masses, small dt or a world origin that is far away.

Overall, I l like the AccumulatedTranslation approach, although it would be nice to have more docs regarding what it's used for and why it exists, at least for other contributors.

It might also make sense to change it to pub(crate) or to remove it from the prelude, since I consider it to be an implementation detail of the solver, and users will almost never need to use it directly.

[LeshaInc]

It might also make sense to change it to pub(crate) or to remove it from the prelude, since I consider it to be an implementation detail of the solver, and users will almost never need to use it directly.

It should be public for people who want to implement their own constraints.

[Jondolf]

Yeah fair, it should be public as long as it's clear what it's used for and why.

[LeshaInc]

I tried removing AccumulatedTranslation from prelude, but it turned out to be very cumbersome due to glob imports everywhere.

[Jondolf]

It would need to be outside of components, maybe in solver or something. But it's probably fine in the prelude for now. The docs are good so I think I can merge this.

I would maybe change the cubes example back to what it was, since the origin is unnecessary for the actual example. Could eventually add a separate example with a very far away origin and/or high mass ratios though. Also FYI, you can convert a Vector to f32 using .as_f32().