Delaunay Tetrahedralization and its dual Voronoi Diagrams

Contents

1 Introduction

2 Related work

3 Mesh Sampling
3.1 Ray Intersection
3.2 Signed Distance Field

4 Constructing DT
4.1 Initialisation
4.2 Predicates
4.2.1 Robustness
4.3 Point Location
4.3.1 Walking
4.4 Incremental flip based Algorithm
4.4.1 Two Dimension
4.4.2 Three dimension
4.4.3 Degeneracies
4.4.4 Time Complexity
4.5 Duality between DT and VD

5 Implementation
5.1 Algorithm : InsertOnePoint(T,p)
5.2 Algorithm : FLIP(T,Ta)
5.3 Robustness
5.4 Results

6 Conclusion
6.1 Limitations
6.2 Future goals

Abstract

The Delaunay tetrahedralization is one of the most popular and common method used for solving problems related to meshes. It is either used for generating a mesh or for breaking it up, as Voronoi diagrams, dual of the DT is a commonly used process for that. The main task of this project is to implement a robust Delaunay Tetrahedralization structure, with a set of points generated from sampling a given 3D Mesh.

Points within the volume of the mesh can be obtained by several methods. We present two such methods and discuss the result obtained. These points serve as vertices for the tetrahedrons that are a part of the combinatorial structure DT.

3D Delaunay Tetrahedralization(DT) is not as optimal as 2D Delaunay triangulations. Implementing them gives rise to several degeneracies, which are quite difficult to handle. In this project, we have implemented a simple Incremental Insertion Algorithm based on the paper presented by Ledoux(2007), inorder to construct the DT structure. Correctness of the structure is given utmost importance rather than its speed.

Acknowledgements

I would like to thank my programme co-ordinator Jon Macey for his assistance. I would also like to thank Dr. Xiaosong Yang for rendering his valuable time to assist me with my project. I would also like to thank Mathieu Sanchez for his guidance in this project and for allowing me to use his SDF library in my project.

Chapter 1

Introduction

Meshes composed of triangles or tetrahedra are used in various applications like terrain databases, Geographical Information Systems, but most demandingly, in Mesh Generation and obtaining partial differential solution using Finite Element Method. There have been several algorithms that explains how to generate the triangles or tetrahedra given a set of points. Although they are used for various applications, the central focus of this thesis is computing Tetrahedrons within a 3D Mesh and thereafter obtain its dual, the Voronoi diagram, using Delaunay Tetrahedralization for the purpose of Fracturing a 3D Mesh.

Fracturing of objects(Destruction) is becoming a key aspect in any major visual effects pipeline. A typical destruction pipeline’s basic step involves preparing the 3d object for fracture. Although, there are several ways to prepare the 3D object, the most common and widely used is shattering using Voronoi diagrams(VD). The VD can be computed directly given a set of points using several algorithms. Pipelines that already exist are either based on Rigid Body Simulations or on a new solution, Finite Element Analysis(FEA) . The majority of the algorithms used in these pipelines do not compute the VD directly as additional computations on intermediate voronoi vertices are performed which greatly reduces the speed of the algorithm. Instead the Voronoi Diagram needed is extracted from its dual, the Delaunay triangulation directly without the need for additional computations and hence has the advantage of speeding up algorithms.

The process of breaking a 3d Mesh into triangles is called Triangulation. This can be achieved by a method called Delaunay Triangulation(DT). DT for a set of points can be obtained by checking the empty circle criterion for any triangle in DT(P). The triangle is said to be Delaunay only if the circumcircle of every triangle is empty, i.e has no other points or vertices. This can be extended to three and higher dimension, when circumsphere is considered. In three dimensions the primitives are no longer triangles, but tetrahedra, hence the process is termed Tetrahedralization. A circumsphere check is considered here to check if the Tetrahedron is Delaunay.

This thesis focuses specifically on the creation of Tetrahedrons using Delaunay Tetrahedralization and thus extracting its dual Voronoi Diagram.

In Chapter 2, we present several work done related to Delaunay Tetrahedralization. We further extend the discussion to its applications.

In Chapter 3, we present the main idea behind generation of points inside the Mesh using two different methods Ray Casting and SDF.

In Chapter 4, the construction of Delaunay Tetrahedrons and the various algorithms used for generating it are discussed. The main algorithm used in this project is explained in the succeeding sub-chapters. The degeneracies that will occur for the algorithm chosen are also discussed.

In Chapter 5, the implementation details of the algorithm are explained. The problems occured while implementing them and the limitations of the algorithm are presented.

In Chapter 6, we evaluate the work done. The drawbacks and advantages of the algorithm implemented are portrayed. Furthermore, its use in recent scenario and future work are briefly discussed.

Chapter 2

Related work

Delaunay triangulation is one of the most popular and most often used methods in problems mainly related to computational geometry. It was discovered in 1934 by the French mathematician Boris Nikolaevich Delone or Delaunay. Many of the Delaunay properties were intensively studied only in 2D for many years. It was extensively researched in the engineering community since the mid - 1980s. For instance, Frey in 1987 presented a solution for eliminating badly shaped triangles from their triangulation by inserting a new vertex at their circumcenters . Similarly Weatherill (1992) gives an alternate solution by inserting vertices at their centroids. These ideas were very helpful when problems related to Mesh generation started evolving in the early 1990s.

Delaunay Triangulations have been constructed in two and three dimensions using several interesting algorithms : sweep- plane, divide-and-conquer and incremental insertion. All these algorithms yields an optimal solution in two dimensions. In 1985, Guibas and Stolfi (1985) proposed a Divide and Conquer algorithm which achieves the optimal bound of O(n log n). Later, Steve Fortune (1987) proposed a Sweepline algorithm for Voronoi diagrams which also achieves this bound. However, when it comes to three-dimensions, things get a bit more complicated. Cignoni et.al(1998) manages to compute DeWall algorithm for constructing the DT in any dimensions, based on the divide-and- conquer algorithm (Fortune 1987). But, Shewchuck(1997) suggests using plan-sweep paradigm for the construction of constrained DT as it is sub-optimal when compared to the former. Unfortunately, these algorithms requires significant programming complexity. Hence, incremental insertion is preferred over the two as it is quite simple to implement and less-error prone.

Incremental insertion can be implemented using two algorithms : Bowyer-Watson Algorithm (1981) and flip-based Algorithm. Flip based algorithm is chosen over Bowyer-Watson as the latter is more error-prone because of the creation of holes, which occur when deleting a conflicting Tetrahedra that does not satisfy the Delaunay criterion. Lawson (1977) developed the first flip-based algorithm in two dimensions. While using the same concept of Lawson for higher dimensions, Joe (1989)proves that the algorithm fails when the union of the two tetrahedra is concave. However, Joe (1991) solves this problem and proves that flip based algorithm works for higher dimensions as well.

In recent years, several Mesh Generators namely Tetgen and Netgen have emerged. Tetgen is a quality tetrahedral Mesh Generator developed by HangSi (2004). It is designed in C++ and is used for constructing Delaunay tetrahedralization, voronoi diagram and convex hull for three-dimensional point sets. Tetgen initially was implemented using a randomized incremental flip-based algorithm of Edelsbrunner and Shah(1996). The algorithm is fast and memory efficient. The latest version of Tetgen as updated on 2009, uses a new implementation of the Bowyer-Watson Algorithm for Delaunay Tetrahedralization. Hang Si, proves that it is faster than incremental flip based algorithm.

Similarly, Netgen is an automatic mesh generation tool in two and three dimension developed by Joachim Schoberl(1994). It uses Constructive Solid Geometry(CSG) as its domain input format for generating tetrahedral meshes in 3D.

Computation Geometry Algorithms Library(CGAL) is an open source library offers data structures and algorithms for creating Triangulations in 2D and 3D. CGAL also uses incremental insertion algorithm, but uses Bowyer-Watson instead of flipping. For locating a point in the structure, it employs Delaunay hierarchy (Devillers 2002).

Chapter 3

Mesh Sampling

A Delaunay tetrahedralization, takes as input a set of points to construct the DT structure. These points are obtained by sampling a given 3D Mesh.

There are several algorithms to obtain points on the surface of a mesh. The easiest being, generating random points along the triangles of a polygonal mesh using Barycentric co-ordinates. But, for certain applications like Fracturing, points are needed in positions along the volume of the source mesh. Inorder to place samples, not just on the surface but also in the volume enclosed by a mesh, we present two methods.

3.1 Ray Intersection

Ray intersection is a very straightforward way of obtaining points inside the volume of a mesh. To begin with, the bounding box of the given mesh is calculated. And then, several random points are generated on the surface of this bounding box. The method we used is the same as obtaining sample points on a triangle using Barycentric co-ordinates as depicted in Figure 3.11(Jose 2011). The bounding box is triangulated and random points are then obtained on the triangles of the bounding box.

Abbildung in dieser Leseprobe nicht enthalten

Figure 3.11 Random points on a triangle

We start with two variables u and v with random values between 0 and 1 such that u + v < = 1 . The resulting point is obtained by the weighted sum of the triangle vertices, such as,

Si = u*A + v*B + w*C where,

w = 1 — (u + v)

After generating random points on the bounding box, we shoot rays iteratively joining two random points on opposite faces of the bounding box. Each ray will intersect the geometry in a number of segments indicated by red dashed lines in Figure 3.12. Every segment is defined by an entry and an exit point. The pros of this approach is that it fits the mesh perfectly and precisely obtain the intersection points of the ray on the mesh . The cons of using this approach is that, it may be a bit slow and the distribution is far from ideal for few samples-although it tends to converge decently when the number of requested samples increases. This can be managed by increasing the density of rays instead of the number of samples generated per segment of the ray. The number of samples on a segment depends on a ‘linear density factor’. It is calculated as length of segment * density, where density depends on the number of samples requested and the volume of the bounding box (Jose 2011).

Abbildung in dieser Leseprobe nicht enthalten

Figure 3.12 Ray intersection

3.2 Signed Distance Field

A distance field function provides the shortest distance to an object from an arbitrary point in space. The sign of the function determines if the point is inside or outside the object. This technique can be used to generate random points within the mesh. As done in the previous method, random points are generated within the bounding box using a simple random function.

Distance field function computes the shortest distance of every point to the mesh, by comparing the distance from the point to each triangle. The sign of the function is computed by a scanning method. A grid is scanned starting from a corner that is definitely outside the object and if a cell is crossed by the surface, the sign is changed. However, this solution works best only for a discrete function. For a continuous function, the sign can be computed using the sign of the dot product between the surface normal at the closest point and the vector from P to the closest point.

The downside in using SDF is that, computing the signed distance field on larger meshes can become time consuming. Computing the distance to every single polygon of the mesh is largely inefficient and does not take advantage of the spatial coherence. Optimization can be done by using spatial partitioning, hierarchy trees and square distances.

Chapter 4

Constructing DT

There are several paradigms in computational geometry for computing the DT in two and three dimensions. However, constructing a DT in three dimensions is quite complicated. In this thesis we have used Incremental insertion algorithm, as this has the complexity O(n 2) which is worst case- optimal. With this algorithm, every point is inserted into the DT one at a time. The tetrahedra containing the newly added point is partitioned and new tetrahedrons are created thus updating the DT structure.

By inserting each point one at a time, we observe that only the Delaunay tetrahedra local to the point is altered and not the entire structure. Whereas in the other two paradigm, divide-and-conquer and sweep plane, the entire structure is built in a single operation. Inserting a new vertex would result in computing the whole structure again. Hence, incremental insertion is mandatory when a dynamic spatial model is built . The insertion of a single point can be achieved with two incremental insertion algorithms

1) the Bowyer-Watson algorithm and
2) flip-based algorithm.

The idea behind Bowyer-Watson(1981) is very simple. Tetrahedrons that does not satisfy the Delaunay criterion when a point is inserted, are deleted. However, this method is prone to errors as it forms a ‘gap’ or ‘hole’ when the tetrahedrons are deleted and affects the overall topological structure as well. Hence, flip based algorithm is preferred over the former. The rest of this chapter we discuss the Incremental Flip-based algorithm.

4.1 Initialisation

Let S be the point set. We first start with a single boundary tetrahedron, constructed in such a way that it encompasses the whole set of points, S. The initial tetrahedron constructed should be several times larger than the range of S as depicted in Figure 4.11 which illustrates a two-dimensional example. This serves as the initial structure for the construction of DT(S), then the points are inserted one by one and the structure is updated. This approach always ensures that the point to be inserted is always within an existing tetrahedron, which avoids unnecessary operations of dealing with points outside the

[...]