New plane solver#905
Conversation
This uses a quasirandom sphere cover generator to sample input velocities, because the number of iterations of the new implementation is dependent on the direction of the input velocity. The existing one can early-out for some input directions (cases where no projection is required) as well.
|
Ok, I've managed to get a decent benchmark set up for this. Performance of both methods vary based on where the velocity is relative to the normal cone (the old one can early-out if all constraints are fine, the new one takes varying numbers of steps depending on whether the initial velocity is allowed, or the end result is constrained by one plane or two). Because of that I came up with a thing for sampling an endless sequence of points on the sphere evenly (based on a quasirandom sequence) to use a different velocity direction for each iteration. Uniform sampling around the whole sphere of directions might not be the most realistic distribution of directions for real-world use, but it's better than nothing. Scaling for the new method looks better than the old one on a linear scale:
On log-log it appears that the old one is actually closer to quadratic than cubic (the brighter lines without markers are power laws, n for red and n^2 for green), and the new one is roughly linear (in case anyone reading this is unfamiliar with log-log plots, power laws with exponent α on linear plots translate to lines with gradient α on log-log plots:
|
|
Explainer up at https://benspiers.co.uk/Games/Velocity-Projection in case it's useful to anyone |
…eren't accurate Gotta account for the cross product being unnormalised when computing distances, or they're meaningless
Guaranteeing that the triple product x0 . (n1 x n2) is nonnegative for a SimplicialCone::Wedge(n1,n2) allows us to skip some sign checks when handling it later. The alpha == 0.0 case when handling Rays gets picked up on the next iteration because the resulting search vector is zero
Not sure if this is too much tbh, but can figure that out in review
This should maybe be equal to normals.len(), or I think we could fairly safely just cap it at 4?
This one has perfect agreement \o/
…ion functions from move_and_slide Not sure if this is the best structure long-term but it gets me out of having to write module docs for now
This means doing a bunch of feature-gating and using Scalar/Vector/Dir and .adjust_precision() around the place
|
Fixed a maths error, and hopefully got CI passing again? Tried to run the steps locally but the tests fill up the memory and then the storage of the laptop I have with me while I'm home for the holidays, then fail because the disk is full 😅 I've added a "test" (it can't fail) that finds the worst agreement between the old method and the new for each of the sets of normals used in the benchmark, and prints what those worst cases are. It (after some searching and head-scratching) pointed me to a problem with my maths, which I've now fixed. The "test" now reports good agreement between the two methods. I also did that optimisation around the ordering in Still needs a pass for docs and cleanup before ready for merging, but it's getting towards a good place now I think. Not sure if we'll want to keep the benchmark around - I see that avian has its own bench system but not sure if it this kind of benchmark would fit into it. It might not be of much continued value now I've demonstrated the performance and scaling anyway. The test isn't any use in an automated suite so it'll either need redoing to make it fallible or removing. Checking that the solutions it finds is optimal might be tricky (it'd basically involve running the old brute-force method to check against, and keeping two implementations in the code seems like it might be weird?) but checking that they're at least feasible ( |
- if one of the constraints is not satisfied (to within DOT_EPSILON) - if the quantity we're supposed to be minimising (|v_out - v_in|) is worse (by more than DOT_EPSILON) for the new method than bruteforce Not sure if the DOT_EPSILON is the right tolerance to use for the second criterion, but it seems to work
…her than me to understand what's going on here
Make the version in velocity_projection pub so it can be used in the benchmark as well, document it because it's pub now, and add the "test" feature to the [[bench]] so that the struct is included when compiling the benchmark
Jondolf
left a comment
There was a problem hiding this comment.
Huge thanks for this work! The perf win is small in absolute terms, but a win is a win 😄 and the algorithm itself is super cool!
I did a quick pass to clean up some of the docs and naming. This looks good to me now, so I'll merge tomorrow unless something pops up
- Rename shadows_enabled to shadow_maps_enabled in examples - Replace EarthlikeAtmosphere resource with Atmosphere::earthlike() method - Rebase on upstream main (new plane solver avianphysics#905)
- Rename shadows_enabled to shadow_maps_enabled in examples - Replace EarthlikeAtmosphere resource with Atmosphere::earthlike() method - Rebase on upstream main (new plane solver avianphysics#905)
- Rename shadows_enabled to shadow_maps_enabled in examples - Replace EarthlikeAtmosphere resource with Atmosphere::earthlike() method - Rebase on upstream main (new plane solver avianphysics#905)
Opening as a draft because it still needs a docs pass, benchmarking to ensure it's actually faster/scales better than the existing algorithm and removing the old code if so.
Objective
The current plane solver is O(n^3) in the number of collision planes. It's not likely to break the bank but it'd be preferable to bring this down to O(n^2) or lower
Solution
Testing
I'm intending to benchmark the new method against the existing one for some arrangements of maybe 1 to 10 normals? to try and establish scalings empiricallyI think there are a few more micro-optimisations I could squeeze out, like imposing a specific ordering of the normals inSimplicialCone::Wedge(n1, n2)relative to the new search direction, which would save a few comparisons and dot products in the relevant branch ofSimplicialCone::project_point.SimplicialConeenum exists to keep the algorithm allocation-free (as opposed to keeping the current simplex points in a vec. Alternatively the first couple of iterations of the loop could be manually unrolled, which would get rid of the match and might speed things up a bit