bepu / bepuphysics2

Pure C# 3D real time physics simulation library, now with a higher version number.
Apache License 2.0
2.25k stars 261 forks source link

Cylinder Convex Hull Tester #263

Closed Wishyou-D closed 1 year ago

Wishyou-D commented 1 year ago

Hi,

I'm reading the CylinderConvexHullTester and I'm quite lost in terms of trying to understand the code. Could you give an explanation of how the Test function works? I would really appreciate it. Thanks in advance!

RossNordby commented 1 year ago

You picked a fun one; that's one of the most complicated pairs in the library!

First, note that these pair testers operate on bundles of pairs, not one pair at a time. The bundles contain Vector<float>.Count lanes, and each lane works independently. You can kind of think of it like running a compute shader, except the lanes of the warp are explicit.

Most pairs have two main stages:

  1. Compute the axis of minimum penetration and the depth along that axis,
  2. Compute contacts between features aligned with that axis.

In the CylinderConvexHullTester, the DepthRefiner does #1. DepthRefiner is an implementation of tootbird search: https://www.bepuentertainment.com/blog/2022/3/10/seeking-the-tootbird

The appropriate cylinder feature (either the side or one of the two circular caps) is then clipped against the most aligned convex hull face. There's some annoyance involved with generating contacts when a cylinder cap needs contacts inside of a convex hull face, but otherwise it's pretty similar to polygon clipping.

Once the contact candidates have been collected for both inner contacts and edge intersections, they're reduced into a manifold of no more than 4 contacts.

Wishyou-D commented 1 year ago

How does the clipping work in detail against the side or cap respectively? What is the numerator and denominator for in the case when the cylinder feature is the side? Does this case also generate 4 contacts?

RossNordby commented 1 year ago

Clipping is effectively a simplified special case of Sutherland-Hodgman. For the cap intersection, the convex face is tested against the cap circle; for the side, the line segment on the side of the cylinder is tested against the convex face.

The numerator/denominator tracking is just there to defer divisions. It avoids needing to compute the true t value for candidates. As far as optimizations go, it's highly questionable and if you were reimplementing something like it, I'd skip it to start with- it's more error prone, more complicated, and the performance benefit is minimal.

Side contacts only generate up to 2 contacts.

Wishyou-D commented 1 year ago

Can you explain "for the side, the line segment on the side of the cylinder is tested against the convex face"? How can you clip a face to a line segment?

RossNordby commented 1 year ago

It's more intuitive if you flip it around- clip the line segment against the convex face.

Wishyou-D commented 1 year ago

When would the tester go into this part? Can you give an example?


if (tMin < tMax && tMin > 0 && candidateCount < maximumCandidateCount)
{
    //Create min contact.
    var newContactIndex = candidateCount++;
    ref var candidate = ref candidates[newContactIndex];
    Unsafe.As<float, Vector2>(ref candidate.X) = hullEdgeOffset * tMin + previousVertex;
    candidate.FeatureId = baseFeatureId + startId;
}`
RossNordby commented 1 year ago

That's a part of the hull face vs cylinder cap test. Each edge of the convex hull face (projected onto the cylinder's cap) is tested against the cylinder cap's bounding circle. Each hull face edge can contribute up to 2 contacts but may sometimes generate only 0 or 1.

Any contacts associated with a face edge are created at the endpoints of a line segment bounded by the original face edge and the circle. That span is represented by an interval from tMin to tMax.

If tMin is 0, that means the previous edge must have extended all the way to the vertex. Since the previous tMax case already created a contact, there's no value in adding another contact at the current edge's tMin = 0; the two contacts would be in exactly the same spot. That's why there's a tMin > 0 condition.

So the tester would only enter that part when the earliest intersection between the edge segment and the circle is somewhere in the middle of the edge.

Wishyou-D commented 1 year ago

Thanks! That was super helpful. What about this part? Could you give an example?


 if (denominatorSquared < max * edgePlaneNormalLengthSquared)
  {
      //As the angle between the axis and edge plane approaches zero, the axis should unrestrict.
      //angle between capsule axis and edge plane normal = asin(dot(edgePlaneNormal / ||edgePlaneNormal||, capsuleAxis))
      //sin(angle)^2 * ||edgePlaneNormal||^2 = dot(edgePlaneNormal, capsuleAxis)^2
      var restrictWeight = (denominatorSquared / edgePlaneNormalLengthSquared - min) * inverseSpan;
      if (restrictWeight < 0)
          restrictWeight = 0;
      else if (restrictWeight > 1)
          restrictWeight = 1;
      var unrestrictedNumerator = a.HalfLength[slotIndex] * denominator;
      if (denominator < 0)
          unrestrictedNumerator = -unrestrictedNumerator;
      numerator = restrictWeight * numerator + (1 - restrictWeight) * unrestrictedNumerator;
  }
RossNordby commented 1 year ago

This is a weird one that's a bit hard to explain. If I'm remembering correctly, it's covering a corner case where the cylinder's side edge is perpendicular to the currently-being-tested convex hull face's edge's plane.

Attempting to expand that sentence into something comprehensible:

  1. Let HullFace be face of the hull that was selected by the normal of minimum penetration.
  2. Let HullFace.Edges[i] be the i-th edge in the selected face.
  3. Let CylinderEdge be the line segment on the side of the cylinder which is furthest along the normal of minimum penetration.
  4. Let GetEdgePlane(i) be a function that takes an edge index (into the HullFace.Edges array) and returns the plane offset and normal. The plane's normal is created from the edge offset (from its vertex A to vertex B) crossed with the normal of minimum penetration. You can think of this as extruding the hull face into a polytope along the normal of minimum penetration, and the edge planes are the planes on the sides of that polytope.

To clip CylinderEdge against HullFace, CylinderEdge is tested against every GetEdgePlane(i) in sequence. Each edge face may bound the interval of intersection.

When the CylinderEdge is near-perfectly flat on a given edge plane, the interval of intersection can become hard to compute numerically. Rather than try to deal with extremely tiny denominators that could explode into NaNs or infinities, such edges can effectively be skipped by "unrestricting" the bounds implied by that edge. Instead of trusting the near-infinity values, it instead just treats that edge plane as not contributing to the intersection interval bounds.

But there's one more concern: discontinuous changes during contact generation are a big source of jitter. Even rare corner cases can show up frequently if your pile of stuff is big enough; if you completely change contact locations instantly at some fixed threshold value, you'll probably get contacts popping into and out of existence from frame to frame.

To help avoid that, the transition between restricted and unrestricted happens gradually over a nonzero width range. An estimate of the angle between the CylinderEdge and the plane surface is used to compute restrictWeight, so the bounds implied by that edge vary smoothly rather than discontinuously.