GJK penetration issue

[Azrael3000]

Hi,

I have just experimented a bit with your library and found a corner case in which gjk breaks. You can find the test case here: https://github.com/Azrael3000/ccd-gjk-issue

The case is simply two cubes (sidelength 2) on top of each other (one with pos 0,0,0 the other with 0,0,1.9).

Output is: GJK: intersect: 0 depth: 0.000000e+00 GJK: dir: -nan -nan -nan GJK: pos: 0.000000e+00 0.000000e+00 9.500000e-01 MPR: intersect: 0 depth: 1.000000e-01 MPR: dir: 0.000000e+00 0.000000e+00 1.000000e+00 MPR: pos: 0.000000e+00 0.000000e+00 9.500000e-01

As can be seen MPR works without a problem but GJK yields the wrong depth and dir. Note that ccd was compiled with double precision accuracy. Have a look if you can reproduce it, otherwise I might try and find some time to actually run it with a fpe trap.

Cheers, Azrael3000

[danfis]

Hi, yes, I can reproduce it and it looks like a bug. What I have tried so far is changing the "first_dir" direction from (1,0,0) which is default to (1,1,0) and it suddenly gives the correct answer. I'm not sure how much time I can spend on this, but I will try to look into it. If you find out where could be the problem let me know.

Cheers, Dan

[Azrael3000]

Hi Dan,

I finally had some time to look into this. Unfortunately I have no solution yet, but some further information:

~~The following happens if USE_DOUBLE is OFF: The segfault ultimately comes from the pointToSegement distance algorithm (https://github.com/danfis/libccd/blob/master/src/vec3.c#L103). The issue is that the points x0 and b are zero. These two points can be traced back to https://github.com/danfis/libccd/blob/master/src/ccd.c#L793 . These points are computed just above. Note that in this case MPR fails as well (by giving completely wrong results)!~~

If USE_DOUBLE is ON: The segfault occurs at https://github.com/danfis/libccd/blob/master/src/vec3.c#L314 which is called from https://github.com/danfis/libccd/blob/master/src/ccd.c#L193

One big difference between double on/off is that the simplex is of size 4 if double is used and 2 with float.

It seems I would need to understand more of the algorithm to actually debug this properly. I do not know whether I will find time for this in the near future. As we use the library in our simulation software it will remain on my TODO list, but the effort for fixing seems quite high unfortunately.

Please ignore all USE_DOUBLE OFF references, that was my mistake.

[Azrael3000]

The simplexToPolytope4 function checks if the origin is on any face of the 4-simplex (i.e. tetrahedron). Now in the case that is investigated, this happens to be the case. However, this distance https://github.com/danfis/libccd/blob/master/src/ccd.c#L574 is -1.776357e-15 which is greater than DBL_EPSILON (at -2e-16 or so). Because the origin is slightly out of the tetrahedron this causes the FPE.

Actually if DBL_EPSILON is increased then no 4-simplex is created as the algorithm already gets a "touch" in the doSimplex3 routine where the ccdVec3PointTriDist2 function is called as well.

So the culprit clearly is this point to triangle distance function. This function only uses basic arithmetic but - and / are not well behaved in a floating point context. So what should be done? The simplest solution is to modify the ccdIsZero function that follows these function calls. If they operate on x*(return value) with 10 < x < 100 the issue would vanish. Obviously x could be determined by properly analysing the error propagation inside the ccdVec3PointTriDist2 function. But I don't see a big issue in conservatively setting it to 100. Opinions are welcome @danfis