Calculation of centroid and polyhedron volume when specifying vertices

Question

Given the location of the vertices of a convex polyhedron (3D), I need to calculate the centroid and polyhedron volume. The following code is available on the Mathworks website .

function C = centroid(P) k=convhulln(P); if length(unique(k(:)))<size(P,1) error('Polyhedron is not convex.'); end T = delaunayn(P); n = size(T,1); W = zeros(n,1); C=0; for m = 1:n sp = P(T(m,:),:); [null,W(m)]=convhulln(sp); C = C + W(m) * mean(sp); end C=C./sum(W); return end 

The code is elegant but very slow. I need to calculate the volume and centroid of thousands of polyhedra hundreds of times. Using this code in its current state is not possible. Does anyone know a better approach or can make this code faster? There are some minor changes that I can think of, for example, replacing mean expression for the mean.

Answer 1

There is a much simpler way to calculate volume with minimal effort. The first fragrance uses 3 local topological information sets of a polyhedron, a tangent unit vector of edges, unit vectors of the normal plane on this tangent, and a unit vector of the face itself (which are very easy to extract from the vertex). See "Polyhedron volume" for details.

The second fragrance uses facial regions, normal vectors and facial barycenters to calculate polyhedron volume in accordance with this Wikipedia article . Both algorithms are quite simple and very easy to implement, and through a simple summation structure you can easily vectorize as well. I believe that both approaches will be much faster than performing a full tessellation of the polyhedron.

Then the centroid of the polyhedron can be calculated by applying the divergence theorem, which translates the integration over the full volume of the polyhedron into the integration over the surface of the polyhedron. A detailed description can be found in Calculating the volume and center of gravity of a polyhedron in 3d . I did not check whether tessellation of the polyhedron into triangles is really required, or whether it is possible to work with more complex polygonal surfaces of the polyhedron, but in any case, the surface tessellation area is much simpler than volumetric tessellation. In general, such a combined approach should be much faster than the volumetric approach.

Answer 2

Reflecting on your only option, if quickhull is not good enough, cudahull if you want accurate solutions. Although, even then you will only get a maximum maximum of 40 times (it seems).

I assume that you have the convex hulls you make have at least 10 vertices (if it is much smaller, you cannot do this). If you do not mind "close enough" decisions. You can create a quickhull version that limits the number of vertices per polygon. The number of vertices that you limit to calculation also allows you to calculate the maximum error if necessary.

The fact is, since the number of vertices on a convex hull approaches infinity, you get a sphere. This means that due to the fast hull, each additional vertex added to the convex hull has less effect * than the previous ones.

* Depending on how quickhull is encoded, this can only be true in the general sense. To put this into practice, you will need to modify the fast recursion algorithm, therefore, while the "next vertex" is always calculated (except that after adding the last vertex or no points for this section) the vertices are actually added to the convex hull in an order that maximizes the increase the volume of polyhedra (probably in order of the most distant to the least distant). You will incur some execution cost to track the order to add the top, but as long as the ratio of pending convex points of the hull to pending points is large enough, it's worth it. As for the error, the best option is probably to stop the algorithm either when the actual convex hull is reached, or when the maximum volume is increased to some part of the current total volume. If performance is more important, then simply limit the number of convex points of the hull to one polygon.

You can also look at various approximate convex shell algorithms, but the method described above should work well for approximating the volume / centroid with the possibility of determining errors.

Answer 3

How much you can speed up your code depends on how you want to calculate your center of gravity. See this centroid calculation answer for your options. It turns out that if you need a centroid of a solid polyhedron, you are mostly out of luck. If, however, only the vertices of the polyhedron have weight, you can simply write

 [k,volume] = convhulln(P); centroid = mean(P(k,:));