Delaunay triangularization of Polyhedron (Python)

Question

I'm trying to get the Delaunay Triangulation of a polyhedron in python so that I can calculate the centroid. I see that there is a Delaunay function in scipy.spatial and that it works in n-dimensions. The trouble is that the documentation shows 2D use and gives me no indication of what to do with higher dimensions. Being able to decompose this object into an array would probably solve this issue for me, but I don't know how to do that.

The problem I'm running into is that I do not know how to verify that this is working correctly as it is outputting an object. I can find nothing on Google about how to graph a polyhedron or how to use this object that scipy is spitting back.

If I do

import numpy as np
from scipy.spatial import Delaunay

points = np.array([[0,0,0],[1,0,0],[1,1,0],[1,0,1],[1,1,1],[0,1,0],[0,1,1],[0,0,1]])
Delaunay(points)

I really would just like to be able to get back the coordinates of these tetrahedrons so that I can calculate the centroids of the polyhedrons. It would also be really nice if I were able to graph the tesselated polyhedron too. I saw in MATLAB that I can do this with a fuction called trimesn, and I found one from matplotlib but it seems to be really different and its documentation is not great.

from matplotlib.collections import TriMesh TriMesh.__doc__ 

u'\n      Class for the efficient drawing of a triangular mesh using\n   
Gouraud shading.\n\n    A triangular mesh is a
:class:`~matplotlib.tri.Triangulation`\n    object.\n    '

You do not need the Delaunay triangulation to compute the centroid. You could compute the appropriate weighted sum of tetrahedra centroids. — Joseph O'Rourke, Jan 16 '16 at 16:01
How do I get the tetrahedrons? In this particular case it is simple, but for an arbitrary hexahedron? — Steven Walton, Jan 18 '16 at 15:56

score 4 · Accepted Answer · answered Jan 15 '16 at 23:23

What tess = Delaunay(pts) returns is an object of the Delanauy class. You can check the tetrahedrons as tess.simplices. It has different attributes and methods. In 2D, for example, it can plot you triangulation, convex hull and Voronoi tesselation.

Regarding the visualization of the final collection of tetrahedrons I didn't find a straightforward way of doing it, but I managed to get a working script. Check the code below.

from __future__ import division
import numpy as np
from scipy.spatial import Delaunay
import matplotlib.pyplot as plt
from mpl_toolkits.mplot3d import Axes3D
from mpl_toolkits.mplot3d.art3d import Poly3DCollection, Line3DCollection
from itertools import combinations


def plot_tetra(tetra, pts, color="green", alpha=0.1, lc="k", lw=1):
    combs = combinations(tetra, 3)
    for comb in combs:
        X = pts[comb, 0]
        Y = pts[comb, 1]
        Z = pts[comb, 2]
        verts = [zip(X, Y, Z)]
        triangle = Poly3DCollection(verts, facecolors=color, alpha=0.1)
        lines = Line3DCollection(verts, colors=lc, linewidths=lw)
        ax.add_collection3d(triangle)
        ax.add_collection3d(lines)

pts = np.array([
            [0,0,0],
            [1,0,0],
            [1,1,0],
            [1,0,1],
            [1,1,1],
            [0,1,0],
            [0,1,1],
            [0,0,1]])
tess = Delaunay(pts)

fig = plt.figure()
ax = fig.add_subplot(111, projection='3d')
for k, tetra in enumerate(tess.simplices):
    color = plt.cm.Accent(k/(tess.nsimplex - 1))
    plot_tetra(tetra, pts, color=color, alpha=0.1, lw=0.5, lc="k")
ax.scatter(pts[:, 0], pts[:, 1], pts[:, 2], c='k')
plt.savefig("Delaunay.png", dpi=600)
plt.show()

The resulting image is

Beat me to it! But much better than what I would've come up with:) — Andras Deak -- Слава Україні, Jan 15 '16 at 23:25

score 2 · Answer 2 · answered Jan 18 '16 at 17:27

You do not need the Delaunay triangulation to compute the centroid of a polyhedron. The centroid is a weighted sum of tetrahedra centroids, where the weight is the volume of each tetrahedron.

You do not need to partition the polyhedron into tetrahedra. First, triangulate the faces of the polyhedron, i.e., quads get partitioned into two coplanar triangles, etc. Next, pick an arbitrary point p in space, say, the origin. Now, for each triangle face (a,b,c), compute the signed volume of the tetrahedron (p,a,b,c). This works if all triangles are oriented counterclockwise. The signed volume takes care of everything via cancellation. Use the signed volume as the weight to multiply the tetrahedra centroids.

Tetrahedra signed volume is explained in Chapter 1 of my book, "Computational Geometry in C."

Delaunay triangularization of Polyhedron (Python)

2 Answers2