Opinion: Beauty In Triangles

[In this reprinted #altdevblogaday-opinion piece, Geomerics' head of technology Sam Martin explains Delaunay triangulations, Voronoi diagrams, and more in exploring the beauty in triangles.] I went through a period where I was continuously dreaming about triangles. Whenever I have a hard problem I find myself mulling it over at night, but around 2003 I was working on navigation meshes and seeing triangles in my sleep every night. I have since broken away from computational geometry and love my new work in lighting, but I would like to reminisce about Delaunay triangulation for old times' sake. Delaunay triangulations are more of a natural discovery than a mathematical construction. They are a 'rich' structure, and giving a single definition seems to under sell them, but we must start somewhere:

Such a triangulation always exists but is only unique if no more than 3 vertices exist on each circumcircle. For example, if 4 points form a square then all 4 points are on the same circumcircle and you can triangulate the square 2 ways. By our definition, both options would be valid Delaunay triangulations. The same principle applies to the vertices of a pentagon, hexagon and so on – they are "co-circular" points. In practice, triangulation algorithms simply impose some form of tie-breaker to ensure they output triangles. This arbitrary decision is a bit messy, so I prefer to think of there being a more natural structure you might call Delaunay Polygonisation which would leave co-circular regions as polygons. Delaunay triangulation with circumcircles Delaunay triangulations have several nice properties, the main one being that they produce the most maximally convex triangles. Avoiding thin long triangles is a universally good thing, although it should be noted that the Delaunay property doesn't prevent thin triangles but it will avoid them where possible. From a practical viewpoint, Delaunay triangulations are a good basis for interpolation and helpful in most mesh structures while being fairly inexpensive to compute. But in my eyes, their beauty comes from the surprising relations that connection them to other structures. A Voronoi diagram is the closest structure and the dual of a Delaunay triangulation. By "dual" we mean that for every Delaunay polygon we can find a corresponding vertex in our Voronoi diagram, and likewise for every Voronoi cell we can find a vertex in our Delaunay triangulation. This connection is so close you can think of them as different representations of the same information. A Voronoi diagram encodes nearest neighbor information explicitly and is composed of convex "cells" rather than triangles. (Cells are like polygons but can be open). The duality between Voronoi and Delaunay allows us to get much of the useful information we may like in a Voronoi diagram, but while working with a triangulation and avoiding the fiddly line-line intersection accuracy issues that can arise in constructing Voronoi diagrams in practice. A Delaunay triangulation and its Voronoi diagram (in red) Voronoi diagrams have a somewhat simpler and more natural definition than Delaunay triangulation. The cells in a Voronoi diagram describe the region of the plane that is closest to a given point. The edges are therefore the boundaries in the plane where we are equidistant from two points. Vertices are locations where we are equidistant from more than 2 points. If we have a Voronoi diagram and a point of interest, by locating the cell that contains the point we can immediately state our closest point and our closest neighbours. Before describing an algorithm to compute these diagrams, I'd like to show you how to draw them by hand. Grab an envelope and try placing some random dots on a page, then draw a convex hull around them. Another remarkable property of Delaunay triangulations is that their boundary is always the convex hull of the set of points. So if you compute the Delaunay triangulation you get the convex hull 'for free'. Random points with convex boundary You might now be able to see triples of vertices, which if you drew a circle through them would not enclose any other points. If you find one, draw in the triangle. There is always only 1 circle you can draw through 3 points but it can be hard to draw free hand for large circles or where more than 3 points are close to being co-circular. However you can improve your accuracy by first finding the centre of the circle. If you take two of the points and draw a separating line half way between them your circumcircle must have a centre on this line. With a bit of imagination and/or scribbling you should be able to see that for every point on this separating line you can draw a circle centred on the line that passes through your two points. The third point then constrains the choice down to a single circle. Points with a circumcircle and separating line If this approach becomes too inaccurate, you can build on the separating lines to construct the Voronoi diagram locally, and use the duality to find the Delaunay triangulation. As you might expect, the separating lines are the basis for edges in the Voronoi diagram. If you draw these relatively accurately, you can inspect where they intersect. At intersections, the nearest neighbors change, and you can manually erase segments that are no longer relevant, giving you the Voronoi diagram. The duality tells us that for every edge in this Voronoi diagram we must have an edge in the Delaunay triangulation – even if the edge in the Voronoi diagram is not directly between the two points that created it, which is sometimes the case. So for every edge you find in your Voronoi diagram, join the two points that created it. Delaunay triangulation with some Voronoi edges When you have co-circular points you will see three or more separating lines intersecting at the same location. If you only joined points based on edges in the Voronoi diagram this would leave you with polygons around these intersection points. So vertices with degrees greater than 3 in the Voronoi diagram correspond to polygons in the Delaunay triangulation. By far the most popular algorithm for computing Delaunay triangulations is Steve Fortune's sweep line algorithm, which is O(n log n) in the number of points. This running time comes from the fact it must sort the points. The algorithm is conceptually elegant and worthy of its own article, although sadly robust implementations are necessarily full of obscure conditions for edge cases and far from easy to read. I actually find the disparity between the elegance of the theory of the sweep algorithm and its real world incarnation interesting. It appears we do not yet have a way to express the edge cases without enumerating them by hand. I would hope that some more general abstraction exists that could generate the implementation. Conveniently, there is a simpler O(n^2) algorithm that exploits another fascinating property of these structures: The Delaunay triangulation of a set of 2D points can be found by computing the 3D convex hull of the same points with a generated z coordinate. Intuitively, what links these structures is the notion of 'nearest'. The shrink-wrapping effect produced by a convex hull turns out to be intimately related to the nearest neighbor property of the Voronoi diagram and similarly, the maximally convex property of the Delaunay triangulation. It is this tight inter-relationship and natural construction that makes structures them so interesting. Given our points in the plane, we shall add a new z coordinate based on their squared distance from the origin. The location of the origin itself is unimportant, but it is common to locate it in the centre of the set of points where we can most clearly see that the points now form a bowl. Perhaps surprisingly, all the downward faces in the convex hull of this new set of points is a triangle (or polygon) in the Delaunay triangulation. So once we have the convex hull, we can simply read off the Delaunay triangulation. Relation of Delaunay triangulations to convex hulls You may now be able to get some intuition as to where the circumcircle comes from. The circle is in some sense the region in which any further point would affect the shrink wrapping of the convex hull. The relation of convex hulls to Delaunay triangulations also offers some explanation for why the boundary of a Delaunay triangulation must always be a convex hull – the boundary of any projection of the 3D convex hull to a 2D plane will always remain convex. This relation between the Delaunay triangulation and its convex hull in one dimension higher is more general than it may first appear. Generalizations of Delaunay and Voronoi structures exist in higher dimensions, and in all cases you can compute them via a convex hull in the dimension above. As a means to a practical implementation the "curse of dimensionality" ensures this insight is not particularly practical, but it is very elegant. Much like Fortune's sweep line algorithm, the elegance of most convex hull implementations suffer from a smorgasbord of edge cases. But if we are to forego any sense of practically, a very small and elegant algorithm does exist, and if you like you can generalize it to higher dimensions. It's downfall is that it's O(n^4). The algorithm is simple: Enumerate every possible triangle we can make from the given set of points; for each triangle construct its 3D plane using the extra z coordinates, and test this plane against all other 3D points. If the plane is in front of all other points it must be on the boundary of the convex hull, and if it faces downwards this triangle is in the Delaunay triangulation. So simple and elegant, but definitely not that practical. This algorithm can be found in "Computational Geometry in C" by Joseph O'Rourke, which I would immediately recommend as the best text book on the subject. Just beware that when he largely brushes the problem of line-line intersection accuracy and stability under the carpet he is side stepping one of the most difficult problems in geometry. [This piece was reprinted from #AltDevBlogADay, a shared blog initiative started by @mike_acton devoted to giving game developers of all disciplines a place to motivate each other to write regularly about their personal game development passions.]