new convex hull implementation

[GoogleCodeExporter]

I implemented the convex hull algorithm of Preparata and Hong. This has a time 
complexity of 
O(n log n), which is much better than what is currently implemented (although I 
do not know the 
exact time complexity of the current code, I got speed-ups by a factor of about 
500 for 6000 
input vertices).

As the algorithm is sensitive to rounding errors, all computations are done 
with integer math. 
The input vertices are discretized to a grid with 10216 positions in each 
direction (that is the 
maximum value with which no overflow occurs). There is also the option to 
shrink the convex 
hull to compensate for the collision margin.

The time complexity could be improved by Chan's algorithm to O(n log h), where 
h is the number 
of output vertices in the convex hull (see 
http://en.wikipedia.org/wiki/Convex_hull_algorithms).

The algorithm was tested in-depth by applying it to more than a million meshes.

I use the algorithm only in my own custum convex hull shape implementation, so 
by applying 
the patch the algorithm is not yet used by Bullet. But integrating it should be 
easy, and it really 
pays off!

btConvexHullShape could be improved by using not only the points of the convex 
hull, but also 
their topology. Then in localGetSupportingVertexWithoutMargin, one does not 
have to go 
through all vertices, but one starts at some vertex and moves along the edges 
in the direction 
which maximizes the dot-product with the given vector. In my custom convex 
shape, I use this 
code:

        const Vector* p = w + edges->getSourceVertex();
        btScalar maxDot = p->x * vx + p->y * vy + p->z * vz;
        setVector(out, *p);
        const btConvexHull::Edge* e = edges;
        do
        {
            p = w + e->getTargetVertex();
            btScalar dot = p->x * vx + p->y * vy + p->z * vz;
            if (dot > maxDot)
            {
                maxDot = dot;
                setVector(out, *p);
                edges = e->getReverseEdge();
                e = edges;
            }
            e = e->getNextEdgeOfVertex();
        } while (e != edges);

Vector is comparable to btVector3, w is a pointer to the Vector-array 
containing the output 
vertices, out is the result (the supporting vertex).

Original issue reported on code.google.com by ol...@arcor.de on 12 Sep 2009 at 6:57

Attachments:

• convexhull_patch.txt

[GoogleCodeExporter]

I forgot to mention that the algorithm handles all degenerate cases including 
the flat 2D case. It also computes 
the faces of the convex hull, which may be arbitrary n-gons in degenerate cases 
(not just triangles).

Original comment by ol...@arcor.de on 13 Sep 2009 at 9:40

[GoogleCodeExporter]

Thanks for the contribution. We are busy on some other projects and 
presentation, but 
will get to it before Bullet 2.76 release.

Original comment by erwin.coumans on 17 Sep 2009 at 9:06

[GoogleCodeExporter]

Original comment by erwin.coumans on 17 Sep 2009 at 9:06

[GoogleCodeExporter]

Original comment by erwin.coumans on 22 Sep 2009 at 11:46

• Added labels: Type-Enhancement
• Removed labels: Type-Defect

[GoogleCodeExporter]

Oh, I just detected a bug in "int 
btConvexHullInternal::Rational64::compare(const Rational64& b) const". Line 
926 has to be
                     : "=&b"(result), [tmp] "=&r"(tmp), "=a"(dummy)
instead of
                     : "=b"(result), [tmp] "=&r"(tmp), "=a"(dummy)

Well, is this a contradiction to my statement that I have tested the code 
in-depth? Yes and no, I simply made 
some few-lines-optimization of the asm statement AFTER testing, this introduced 
the bug. 

Original comment by ol...@arcor.de on 29 Sep 2009 at 1:50

[GoogleCodeExporter]

the current Bullet code should perform pretty well for btConvexHullShape with 
less than 100 vertices. We 
don't recommend over 100 vertices for game purposes, and for less than 100 
vertices, your implementation 
isn't likely make any difference.

6000 input vertices ;-) ?

Can you implement it as a new btDenseConvexHullShape (or something like that), 
so we can provide it as a 
special option?

Original comment by erwin.coumans on 3 Nov 2009 at 7:11

[GoogleCodeExporter]

OK, for a small number of vertices, it will be faster to just go through the 
list of vertices in 
btConvexHullShape. But the contribution is also about a new algorithm to 
compute the convex hull. I had 
some problems with the old one for cases which were nearly degenerated (and 
they had less than 100 
vertices). The Preparata&Hong-algorithm which I implemented is not only faster, 
but also more robust.

The algorithm is already in use in CINEMA 4D MoDynamics, and there the users 
often work with meshes of 
6000 and more vertices. Is there a problem with convex hulls of about 1000 
vertices except speed? I mean, in 
Bullet a built-in shape like a cylinder is handled by the same collision 
algorithm, and as a convex hull a 
cylinder has an infinite number of supporting vertices.

Original comment by ol...@arcor.de on 4 Nov 2009 at 10:04

[GoogleCodeExporter]

Can you please create a new patch using a new shape type (for example 
btDenseConvexHullShape), rather than modifying the existing btConvexHullShape?

By the way, Bullet doesn't build a convex hull for a btCylinder, it uses a 
function to 
compute the supporting vertex. Apart from performance, there should be no 
problem using 
many vertices.

Original comment by erwin.coumans on 13 Dec 2009 at 3:48

[GoogleCodeExporter]

I'll do so once I can find some time.

Yes, of course Bullet doesn't use a convex hull for a cylinder. I only 
mentioned the cylinder to point out that 
there is no problem besides performance with many vertices. And the performance 
is improved by my 
implementation. I'll make a comparison to get some idea when it pays off to 
move along the edges towards the 
maximum dot product.

Original comment by ol...@arcor.de on 14 Dec 2009 at 7:35

[GoogleCodeExporter]

Hi,

>>the contribution is also about a new algorithm to compute the convex hull

I'm curious where exactly the convex hull is computed? Bullet doesn't compute 
the convex hull for 
btConvexHullShape (it simply uses all points, and interior points just don't 
contribute to the convex hull, only 
cost performance. There should be no problem/robustness issues there.

I hope you can find some time in 2010 for a btDenseConvexHullShape.
Thanks!
Erwin

Original comment by erwin.coumans on 19 Jan 2010 at 8:24

[GoogleCodeExporter]

>> I'm curious where exactly the convex hull is computed?

In LinearMath/btConvexHull, which seems to be a third-party contribution ported 
to Bullet. I implemented a new 
algorithm to compute the convex hull of a point cloud which is much faster than 
the current one and which also 
handles all degeneracies. Maybe that's not so interesting for the use of Bullet 
as a game engine because then the 
convex hull is probably be computed in advance by some authoring tool. But for 
all those who use Bullet in an 
application and have to compute the convex hull in that same application (like 
we have to do for MoDynamics in 
CINEMA 4D), a good convex hull algorithm is really helpful.

Original comment by ol...@arcor.de on 22 Jan 2010 at 8:04

[GoogleCodeExporter]

Very nice code.

I can tell that people loading "external" (from project) or "online" objects in 
simulation would be happy to use that.

Even If I do serialize and save convex  hull computation results for each new 
object, first time is long, so if computing times can be lowered, it's great.

Some more happy people "fast, robust, scales well"
http://digestingduck.blogspot.com/2010/12/computational-geometry-sucks.html

Original comment by tuan.kur...@gmail.com on 10 Dec 2010 at 4:57

[GoogleCodeExporter]

Finally added it in LinearMath/btConvexHullComputer.*
http://code.google.com/p/bullet/source/detail?r=2353

The #define USE_X86_64_ASM should be added in the build system (or manually 
enabled) to avoid compilation issues on exotic platforms.

Thanks a lot and sorry for the delay.

Original comment by erwin.coumans on 22 Mar 2011 at 12:54

• Changed state: Fixed

[GoogleCodeExporter]

This sounds great.  Can the ConvexDecomposition code be adapted to take 
advantage of this?  It could really use a performance boost. (even just the 
Demo's shrink part for margins would benefit)

Original comment by markri...@hotmail.com on 26 Mar 2011 at 9:59

[GoogleCodeExporter]

The comment for getNextEdgeOfFace() states it will return a "clockwise list of 
all edges of a face", however in my tests it sometimes returns them clockwise, 
sometimes counter-clockwise, inconsistently depending on the input points for 
different convex outputs (although all faces of the same convex mesh will be 
consistently clockwise or counter-clockwise).

So either that is a bug, or the comment should be changed.  Thanks.

Other than that it's working very well so far, thank you very much for this, it 
will be very useful.

Original comment by markri...@hotmail.com on 28 Mar 2011 at 2:55

[GoogleCodeExporter]

Clockwise / counter-clockwise reproduction case:

btConvexHullComputer* convexhc = new btConvexHullComputer();
btAlignedObjectArray<btVector3> vertices;
vertices.resize(10);

vertices[0] = btVector3(0.92, 1.85, -0.59);
vertices[1] = btVector3(0.92, 1.85, 0.24);
vertices[2] = btVector3(0.92, 0.45, -1.76);
vertices[3] = btVector3(0.92, -0.15, -1.76);
vertices[4] = btVector3(0.69, 0.28, -1.76);
vertices[5] = btVector3(-0.44, 1.85, 0.24);
vertices[6] = btVector3(0.69, -0.15, -1.76);
vertices[7] = btVector3(-0.85, 1.55, 0.24);
vertices[8] = btVector3(0.92, -0.15, 0.24);
vertices[9] = btVector3(-0.85, -0.15, 0.24);

// Face-edge output is clockwise
convexhc->compute(&vertices[0].getX(), sizeof(btVector3), 
vertices.size(),0.0,0);

vertices[0] = btVector3(1.85, 0.56, 0.94);
vertices[1] = btVector3(1.85, -0.84, -0.24);
vertices[2] = btVector3(1.85, 0.56, -0.24);
vertices[3] = btVector3(1.61, -1.01, -0.24);
vertices[4] = btVector3(0.48, 0.56, 1.76);
vertices[5] = btVector3(-0.15, 0.56, 1.76);
vertices[6] = btVector3(-0.15, 0.26, 1.76);
vertices[7] = btVector3(0.07, 0.27, 1.76);
vertices[8] = btVector3(-0.15, -1.12, -0.24);
vertices[9] = btVector3(-0.15, 0.56, -0.24);

// Face-edge output is counter-clockwise
convexhc->compute(&vertices[0].getX(), sizeof(btVector3), 
vertices.size(),0.0,0);

Original comment by markri...@hotmail.com on 28 Mar 2011 at 3:33

[GoogleCodeExporter]

Actually I got some meshes where the winding of the hull is inconsistent (both 
clockwise and counter clockwise).

So I post-processed the faces and detect the winding and adjust it. We could 
add this into the btConvexHullComputer as well. See my snippet here:
http://tinyurl.com/3fys5nn

Original comment by erwin.coumans on 7 Apr 2011 at 5:41

[GoogleCodeExporter]

[deleted comment]

[GoogleCodeExporter]

Oh yes, there's really a bug with the order, thanks for finding! The order is 
consistent in the internal algorithm, but the algorithm permutes the XYZ-axes 
(internal Y-axis has highest extent, Z-axis has lowest extent). So depending on 
the permutation, the btConvexHullComputer gets a wrong orientation.

I attached a patch which should solve the issue. Now face-edges are listed 
counter-clockwise, vertex-edges clockwise. Please verify that it's working 
consistently now.

Even with the old version the orientation should have been consistent within a 
single hull. Are you sure, Erwin, that you had an inconsistent orientation 
within a single hull?

Original comment by ol...@arcor.de on 12 Apr 2011 at 9:01

Attachments:

• patch.txt

[GoogleCodeExporter]

Great!  Appears to be working very well so far (double precision) for me with a 
few thousand random point objects.  Nice.  Thanks.

Original comment by markri...@hotmail.com on 12 Apr 2011 at 1:30

[GoogleCodeExporter]

Hey Ole,

It was a bug on my side, things were consistent one way or the other.

I found another issue though: sometimes planar triangles are not merged into a 
single polygon. For example, if you change the vertex cound of 
Demos\InternalEdgeDemo\Taru.mdl\Taru.mdl in 18, like:

#define TaruVtxCount 18

and compile Demos/InternalEdgeDemo, you can see that one of the flat sides of 
the convex hull is tesselated. Is there a threshold to control this?

Original comment by erwin.coumans on 13 Apr 2011 at 12:10

[GoogleCodeExporter]

The point coordinates are converted to internal integer coordinates with an 
approx. range of -5100 ... 5100, and the convex hull computation is done for 
the integer coordinates. Probably the triangles aren't exactly co-planar any 
longer in the integer coordinates.

There's no threshold to make nearly co-planar faces co-planar.

Original comment by ol...@arcor.de on 13 Apr 2011 at 6:49

[GoogleCodeExporter]

I had a look at the values and the convex hull computation is extremely 
sensitive,
changing 2.47058e-007f into 0.f made it work (see vertices below)

What would be the best way to modify the convex hull computation to deal with 
this numerical issues? Or should we pre-process the input vertices and weld 
x/y/z coordinates that are very similar (within a small tolerance)?

static float TaruVtx[] = {
1.08664f,-1.99237f,0.0f,
0.768369f,-1.99237f,-0.768369f,
1.28852f,0.f,-1.28852f,
1.82224f,0.f,0.0f,
0.0f,-1.99237f,-1.08664f,
0.0f,0.0f,-1.82224f,
0.0f,-1.99237f,-1.08664f,
-0.768369f,-1.99237f,-0.768369f,
-1.28852f,0.f,-1.28852f,
0.0f,0.0f,-1.82224f,
-1.08664f,-1.99237f,0.f,
-1.82224f,0.f,0.f,
-0.768369f,-1.99237f,0.76837f,
-1.28852f,0.f,1.28852f,
0.f,-1.99237f,1.08664f,
0.f,0.f,1.82224f,
0.768369f,-1.99237f,0.768369f,
1.28852f,0.f,1.28852f, //1.28852f,2.47058e-007f,1.28852f,
};

Original comment by erwin.coumans on 13 Apr 2011 at 3:13

[GoogleCodeExporter]

The algorithm itself is designed to work with exact arithmetics. There are no 
epsilons, and I think this is why it works so reliably.

So if you want to solve your issue with co-planar faces, you either have to 
make sure that the input data is such that the faces are exactly co-planar 
after conversion to the internal integer coordinates (which seems to be 
difficult as you don't know the faces in advance if you start with a point 
cloud -- that's the purpose of the convex hull algorithm), or you have to merge 
nearly co-planar faces of the convex hull result. That seems to be easier.

But is it really a problem if the convex hull result contains some nearly 
co-planar faces? The collision algorithm uses only the vertices of the convex 
hull, so for this purpose only nearly co-planar vertices (which may also occur) 
can lead to reduced performance.

Original comment by ol...@arcor.de on 14 Apr 2011 at 7:21

[GoogleCodeExporter]

The new polyhedral contact clipping uses the two nearest faces (one on each 
object) for clipping. For example: when a tesselated cube rests on a larger 
tesselated cube (centered) it is best to generate 4 contact points. If 
co-planar triangles are not detected/merged, the clipping code will not produce 
those 4 points in one go.

We could merge nearly co-planar faces in the convex hull result, but perhaps 
there is another way?

I can imagine that introducing interval arithmetic in the convex hull 
computation (when comparing normals for example) could fix this issue as well. 
I'm just not familiar enough with the implementation (yet).

Original comment by erwin.coumans on 14 Apr 2011 at 4:59

[GoogleCodeExporter]

Hi Ole, I applied your latest patch here:

http://code.google.com/p/bullet/source/detail?r=2389

Thanks again!

Original comment by erwin.coumans on 15 Apr 2011 at 6:38

[GoogleCodeExporter]

I implemented the removal of co-planar faces in the output of the 3d convex 
hull, using the 2d Graham scan convex hull. This improves stability and 
performance for convex polyhedra.

http://code.google.com/p/bullet/source/detail?r=2397

Original comment by erwin.coumans on 20 May 2011 at 12:31

[GoogleCodeExporter]

[deleted comment]

[GoogleCodeExporter]

The algorithm gives wrong results when I test it with the following vertices:

btConvexHullComputer* convexhc = new btConvexHullComputer();
btAlignedObjectArray<btVector3> vertices;
vertices.resize(8);

vertices[0] = btVector3(2.0f, 3.535533905932738f, 3.535533905932737f);
vertices[1] = btVector3(4.0f, 2.0f, 0.0f);
vertices[2] = btVector3(0.0f, 2.0f, 0.0f);
vertices[3] = btVector3(1.0f, 0.0f, 0.0f);
vertices[4] = btVector3(4.0f, 1.414213562373095f, 1.414213562373095f);
vertices[5] = btVector3(0.0f, 1.414213562373095f, 1.414213562373095f);
vertices[6] = btVector3(3.0f, 0.0f, 0.0f);
vertices[7] = btVector3(2.0f, 5.0f, 0.0f);

Original comment by mojtahed...@gmail.com on 2 Sep 2014 at 5:26

[GoogleCodeExporter]

@mojtahedzadeh Thanks for the report. We moved to github, so this issue tracker 
is not active any more. Moved the issue here: 
https://github.com/bulletphysics/bullet3/issues/262

How did you run into this problem? Was it specifically generated to test the 
convex hull algorithm/implementation?

Original comment by erwin.coumans on 8 Oct 2014 at 8:30

[GoogleCodeExporter]

@erwin.coumans: Thanks for moving my comment into the active webpage of the 
project.

The vertices are designed to test 3D convex hull algorithms.

Original comment by mojtahed...@gmail.com on 5 Nov 2014 at 2:50

[GoogleCodeExporter]

mojtahedzadeh: I reproduced the issue and it seems to give the right results. 
What makes you think it goes wrong?
Could you continue the discussion in the github tracker at 
https://github.com/bulletphysics/bullet3/issues/262?

Original comment by erwin.coumans on 18 Nov 2014 at 11:03