Coplanarity is teH eVil!!!11!1!!

I sort of implied this in passing in my triangle-triangle tests post, but let’s state it officially: coplanarity (or collinearity) in your inputs will almost always screw you over in one way or another! Here’s another example of what can go wrong: this time performing continuous collision detection between a moving point and a moving triangle. (In case this looks familiar it’s probably because you’ve seen me mention it before in other forums, but hey, self-plagiarization is OK as long as you’re not trying to publish!)

Let’s start by looking at the methods proposed for this problem. One of the first papers to address the continuous collision detection problem of moving point vs. moving triangle was [Moore88]. Their solution was to solve the parametric vector equation

p+tv = p0 + tv0 + u((p1p0)+t(v1v0)) + v((p2p0)+t(v2v0))

where p is the point with velocity v over the time step, the pi and vi terms are the triangle vertices and their velocities, and t, u, and v are the parametric variables. This resolves to a 5-th order polynomial in t. Moore and Wilhelms solved the polynomial using a binary search technique (to guarantee convergence). The full test only needs to be performed if the signed distance to the triangle plane of the point in its start and end position have different signs.

A 5-th degree polynomial isn’t exactly appetizing, so in [Provot97], Xavier Provot suggested a different method to solve the same problem. Provot notes that a necessary (but not sufficient) condition for the point to intersect the triangle is for it to lie in the plane of the triangle. We can test coplanarity using the trusty scalar triple product. (Hmm… someone, remind me to blog about the scalar triple product.)

Given a point p with velocity v and a triangle abc with vertex velocities va, vb, and vc, we can define the parametric vectors

ap(t) = (p+tv) – (a+tva) = (pa)+t(vva),
ab(t) = (ba)+t(vbva), and
ac(t) = (ca)+t(vcva).

p will lie in the plane of abc when

Dot(ap(t), Cross(ab(t), ac(t))) = 0,

that is, when the vector ap lies in the plane of ab and ac.

In other words, what Provot suggests is that we solve the above equation for t, and for each t-value found, we bring the point and triangle to their positions for that t-value and perform a containment test to see if the point lies in the triangle. The smallest t for which the point lies inside the triangle is returned as the time of intersection. This approach of Provot’s has been been referenced/used in a lot of places, e.g. in [Bridson02], [Redon04], and [Hutter07].

So, well, this is great and all, except for one niggly little thing… it doesn’t work!

The problem is that when the point and the triangle move together in the same plane you cannot solve for a unique t value for which the point lies in the plane of the triangle, and the method breaks down. (This is also the case for the Moore-Wilhelms method, btw.) In this context it’s somewhat ironic that [Bridson02] has ‘robust’ in the title and chides Provot for “not properly account[ing] for rounding errors”, but suffers from unaccounted-for robustness issues in its own use of Provot’s method in exactly the case I outlined. Oops! To be robust you’d have to special-case this situation, which turns your reasonably elegant code pretty ugly because you now have to solve the problem in 2D after detecting this degenerate case. Plus you need tolerances for deciding when you drop to 2D, etc. It’s u-g-l-y.

Unfortunately, I’m not aware of a better solution. You could probably use interval arithmetic and solve the problem that way, perhaps with interval-Newton, but that doesn’t exactly categorize as a better solution in my book (no pun intended). Well, it would be completely robust I guess, so better in that sense at least.

BTW, note that an approach very similar to Provot’s was described in [Schömer95] (also without noting the coplanarity problem), so perhaps the method should not be attributed to Provot, as is typically the case, but to Schömer.

OK, so lots of verbiage from Christer again, but what’s the takeway? Well, this: if you don’t have special-case code for coplanarity and collinearity in your geometric calculations, odds are you have (possibly severe) bugs in your code!

References

  • [Bridson02] Bridson, Robert. Ronald Fedkiw. John Anderson. “Robust treatment of collisions, contact and friction for cloth animation,” Proceedings of SIGGRAPH’02, San Antonio, Texas, pp. 594-603, 2002. http://graphics.stanford.edu/papers/cloth-sig02/
  • [Hutter07] Hutter, Marco. Arnulph Fuhrmann. “Optimized Continuous Collision Detection for Deformable Triangle Meshes,” Journal of WSCG, Vol.15, No.1-3, pp. 25-32, February 2007. http://wscg.zcu.cz/wscg2007/Papers_2007/journal/D11-full.pdf
  • [Moore88] Moore, Matthew. Jane Wilhelms. “Collision Detection and Response for Computer Animation,” Computer Graphics (Proceedings of SIGGRAPH ’88), vol. 22, pp. 289-298, 1988.
  • [Provot97]Provot, Xavier. “Collision and self-collision handling in cloth model dedicated to design garments,” Computer Animation and Simulation ’97, pp. 177-190, 1997
  • [Redon04] Redon, Stephane. Young Kim. Ming Lin. Dinesh Manocha. “Fast Continuous Collision Detection for Articulated Models,” ACM Symposium on Solid Modeling and Applications, June 2004. http://gamma.cs.unc.edu/Articulate/
  • [Schömer95] Schömer, Elmar. Christian Thiel. “Efficient collision detection for moving polyhedra,” Proceedings of the Eleventh Annual Symposium on Computational Geometry, pp. 51-60, 1995. http://citeseer.ist.psu.edu/article/omer95efficient.html

9 thoughts on “Coplanarity is teH eVil!!!11!1!!”

  1. Hi Christer,

    Conservative Advancement doesn’t suffer from planar problems, for continuous collision detection.

    Also, Stephane Redon deals with the planar case in his paper “An Algebraic Solution to the Problem of Collision Detection for Rigid Polyhedral Objects”. Scroll down at the end of http://i3d.inrialpes.fr/~redon, and search for ‘planar’ in the paper. I implemented this, and it seems to work fine.

    Hope this helps,
    Erwin
    http://bulletphysics.com

  2. Took me a while to reply to this, but here goes…

    Erwin, I’m not sure what you mean when you claim that “Conservative Advancement doesn’t suffer from planar problems, for continuous collision detection” and how it relates to the problem I talked about.

    Conservative advancement to me refers to iterative methods whereby you compute a tentative, conservative time of impact (TOI), advance the simulation to this time, determine if a termination condition has been met, and if not, you repeat for another iteration.

    The issue I outlined in the post is that as long as you treat the coplanar situation as a 3D problem, the computed TOI is always zero! When the TOI is zero, conservative advancement breaks down because a step of zero doesn’t advance you any.

    Here, to compute a TOI that isn’t zero, you have to detect that what you have is a degenerate situation and switch from treating the problem as a 3D problem to instead treating it as a 2D problem (or do something else that is equivalent to this).

  3. No, Conservative advancement includes all the iterations to calculate the actual time of impact.
    Only the internal steps are conservative, the endresult is the actual TOI. See the FAST or CATCH paper by Young Kim for more info about this process.

    http://graphics.ewha.ac.kr/fast/
    http://graphics.ewha.ac.kr/CATCH/

    The TOI using CA is not zero for this planar case, even in 3D.

    On the algebraic version, Stephane Redon detects this 2D planar case as one of the special cases.

  4. Thanks Erwin. Your point is most definitely a valid one, but before agreeing I first have to disagree some on terminology!

    I think it is pretty safe to say that the term “conservative advancement” (CA) describes a generic approach, and not a particlar approach. This is supported by Mirtich’s “Timewarp Rigid Body Simulation” paper, which (in Section 2.2) refers to CA as “an alternative to RD based on the idea of never integrating over a discontinuity.” This is in contrast to what he calls retroactive detection: “Retroactive detection (RD) is the most common approach to handling discontinuities. The simulator takes small steps forward and checks for discontinuities after each step.” Clearly both are generic concepts that do not prescribe how the actual steps are performed (other than being conservative, for CA).

    In this generic sense, Provot’s method is a CA method. Of course, it’s now clear to me that you’re talking specifically about the “Mirtich-style” CA of e.g. Mirtich or of Zhang et al and not about the Provot method at all (which is how I first interpreted what you said).

    So, terminology aside, you’re absolutely right. If we switch to a method that does conservative advancement proportional to the minimum distance between the point and the triangle we avoid the problem of Provot’s method (but at the cost of computing the closest points for each iteration).

    In addition to Mirtich’s and Zhang et al’s methods that you point at, the interval halving method I describe in Section 5.5.1 of RTCD can also solve the moving point vs. triangle problem thusly.

    BTW, I don’t have a reference in my book for the interval halving method, but after writing it I learned it is pretty much identical to one by Kenji Shimada (US Patent 4,888,707, Dec. 19, 1989), which I believe is also described in “Collision detection among many objects in simulation” by Shimada, Okano, and Kawabe. The Shimada method is the earliest published CA method I’m aware of (although I wouldn’t be surprised if there’s some early (non-)linear programming paper that dealt with all this).

    Clive Johnson of EA Canada pointed me to a recent paper by Schwarzer, Saha, and Latombe on adaptive bisection, which is also effectively a reinvention of Shimada’s method.

  5. Some reflections:

    Conservative Advancement (“Mirtich-style”) can be used to find the TOI between edge-edge or point-triangle in the case where these primitives are “padded” with a “collision radius” (i.e. capped cylinder vs. capped cylinder and sphere vs. “capped triangle slab”). That is, the TOI where the surfaces of these padded primitives touch, instead of where the infinitely thin center skeletons (triangle, line segment, point) coincide or cross each other. I am not (as yet) aware of how to find that TOI with a “Provot-style” method.

    In your post on triangle-triangle tests, you mention that many papers “typically just return 0” for degenerate cases. Now, I just found this kind of “sloppiness” in a piece of code in your book … ;-) Namely the ClosestPtSegmentSegment() function in chapter 5.1.9 “Closest points of two line segments”. If the two line segments are found to be parallel, you pick an arbitrary value (here 0) for the s parameter on segment S1.

    Though not causing robustness problems in the algorithm itself, choosing 0 for s here has implication if I use this calculation in, for example, a capped cylinder collision routine.

    If my two capped cylinders are exactly parallel in direction and partly overlapping when projected onto this common direction, then the pair of closest points is reported at an arbitrary end of one of the cylinders. We can in fact do better in this special case. I’d like the closest points to lie at half the overlapping interval when projected onto the common direction, like this (fits into your code) :

    if (denom != 0.0f) {
    s = …
    } else {
    // Segments are parallel, find midpoint of overlapping segment parts.
    // (Not optimized since this is a special case…)

    float s1 = Clamp((b*0.5f – c) / a, 0.0f, 1.0f);
    float t1 = Clamp((b*0.5f + f) / e, 0.0f, 1.0f);

    float s2 = Clamp((b*t1 – c) / a, 0.0f, 1.0f);
    float t2 = Clamp((b*s1 + f) / e, 0.0f, 1.0f);

    s = 0.5f * (s1 + s2);
    t = 0.5f * (t1 + t2);

    return Dot(c1 – c2, c1 – c2);
    }

    This way a collision response could be correctly applied to all four segment endpoints according to the s and t values returned from this routine (otherwise the weighting of the response would be heavily tilted towards the point with parameter (s = 0) on segment S1).

    Best regards,
    Göran W.
    http://www.surgical-science.com/

  6. And, yes, I can say “Oops!?”. The code snippet I posted above was way wrong…

    Here is a corrected version of the code for the parallel case:

    // Segments are parallel, find midpoint of overlapping segment parts.

    // First, compute the points on segment S2 closest to the endpoints of segment S1.
    // Then project these points back onto segment S1.

    float s1 = 0.0f;
    float t1;
    t1 = Clamp((b*s1 + f) / e, 0.0f, 1.0f);
    s1 = Clamp((b*t1 – c) / a, 0.0f, 1.0f);

    float s2 = 1.0f;
    float t2;
    t2 = Clamp((b*s2 + f) / e, 0.0f, 1.0f);
    s2 = Clamp((b*t2 – c) / a, 0.0f, 1.0f);

    // Finally, compute the midpoints of these intervals on each segment.
    s = 0.5f * (s1 + s2);
    t = 0.5f * (t1 + t2);

    return Dot(c1 – c2, c1 – c2);

    Hope I got it right this time! ;-)

    /Göran W.

  7. Indeed, the terminology conservative advancement is just refering to the general approach, of taking a conservative large step, using the projected linear and angular motion.

    >>(but at the cost of computing the closest points for each iteration).
    Note that Gino van den Bergen came up with a way to merge the two iteration loops of closest points (gjk) and conservative advancement.

    Unless I miss something (I probably do, as I only skimmed through it), the paper you link to http://robotics.stanford.edu/~latombe/papers/adaptive-bisection/final.pdf uses bisection, dividing the segment in the middle, see Figure 2, instead of using the projection.
    Isn’t that an important difference?

    Erwin

Leave a Reply