Tetrahedralization

In the three-dimensional case the simplex is a tetrahedron, so the analogue of triangulation is tetrahedralization: the volumetric domain must be partitioned into tetrahedra so that they do not overlap, cover the domain completely and meet pairwise only along shared faces, edges or vertices. As in the two-dimensional case, we require the Delaunay condition — now for the circumscribed sphere of a tetrahedron: it must contain no other mesh vertices. This condition helps avoid some badly elongated elements, but by itself it does not guarantee a good shape for every tetrahedron; this is why we also prescribe quality criteria and run mesh optimization.

The minimal three-dimensional simplex is built on four points:

Why the three-dimensional case is harder

The classical divide-and-conquer scheme, described for triangulation in the appendix, becomes much more involved in three dimensions. When stitching two parts in the plane, the gap between them is a pair of edge chains, and filling it is easy. In space the gap is formed by a whole set of faces, and it has to be closed with tetrahedra correctly. Moreover, the local Delaunay restructuring in the plane is a simple diagonal swap of a quadrilateral (the 2–2 flip), while in space the elementary restructurings (the 2–3 and 3–2 flips) change the very number of tetrahedra and require careful admissibility checks. Finally, unlike the plane, not every three-dimensional domain can be partitioned into tetrahedra without adding new points — the classical example is the Schönhardt polyhedron.

For these reasons we do not implement the low-level tetrahedralization by hand but delegate it to the same CGAL library (see below), focusing on the problem statement: describing the domain, its boundary and the mesh quality requirements.

The domain as a piecewise linear complex

The three-dimensional computational domain is specified by its boundary — a closed, oriented, triangulated surface. Such a description is called a piecewise linear complex (PLC): a set of vertices, edges and planar faces consistently covering the body's surface. It is the direct three-dimensional analogue of the boundary contours of the two-dimensional triangulation: there the boundary consisted of closed polylines, here — of closed triangular meshes.

For three characteristic bodies the surface is built explicitly:

• cube — eight vertices and twelve triangular faces (two per each of the six sides);
• ball — an icosphere: an icosahedron whose faces are repeatedly subdivided into four and projected onto the sphere;
• hollow cylinder — the outer and inner lateral surfaces plus two annular end caps; the result is a closed surface of genus 1 enclosing a body with a through hole.

Constraints: edge and face recovery

The faces and edges of the boundary surface play the role of constraints — just like the non-crossable edges in the two-dimensional constrained Delaunay triangulation. They must be present in the final mesh: tetrahedra must not cross the body's boundary through their interiors, and sharp edges (the cube's edges, the rims of the cylinder caps) must be preserved as mesh edges rather than be “shaved off”. Sharp edges are detected automatically by the dihedral angle between adjacent faces and are recovered separately, after which the surface remains embedded in the volumetric mesh.

Quality criteria

Unlike the planar case, where the density was set merely by the number and placement of points, the volumetric mesh is driven by sizing fields and shape criteria:

• a bound on the tetrahedron size controls the mesh density: in CGAL it is set by the cell_size parameter, an upper bound on the circumradius of a mesh tetrahedron. In the demonstration utility the target scale is estimated as $\sqrt[3]{V/n}$, where $V$ is the body volume and $n$ is the prescribed number of points;
• a bound on the ratio of the circumscribed sphere radius to the shortest edge (the radius–edge ratio) prevents many elongated “needle” tetrahedra; the flat “slivers” that this criterion does not catch are removed by the final optimization;
• boundary-approximation criteria (the angle, size and distance error of boundary facets) control the accuracy of the surface mesh, while for the ball and the cylinder the accuracy also depends on how well the initial polyhedral surface approximates the intended curved surface.

The algorithm inserts interior points (Delaunay insertion) and refines the mesh until all criteria are met. A final optimization step then improves element shapes through vertex perturbation and sliver exudation, that is, vertex reweighting used to remove nearly flat tetrahedra.

Implementation

The low-level construction is delegated to the Mesh_3 module of the CGAL library — the same one used for the two-dimensional triangulation. The boundary surface of each body is represented using Mesh_polyhedron_3<K>::type, the domain by Polyhedral_mesh_domain_with_features_3 (with automatic sharp-edge detection), and the tetrahedralization itself is performed by calling make_mesh_3 with the criteria described above (Mesh_criteria_3). The exact-predicates kernel Exact_predicates_inexact_constructions_kernel is used, so the in-sphere and orientation checks are evaluated robustly. The output is a set of tetrahedral simplices of the form (4.3)$\displaystyle T$ $\displaystyle = (p_1, p_2, p_3, p_4).$, from which the local matrices are later built and the global finite-element system is assembled.

Examples

All meshes below are real outputs of the tetrahedralize utility of the computational core; the text gives the parameters of the full tetrahedral mesh. The interactive models show its surface shell: the surface can be rotated with the mouse or a finger, and the lighting-based shading emphasizes boundary facets of the tetrahedra.

The figure “Cube” shows the tetrahedralization of a unit cube with $n$ $= 600$: 674 vertices, 2809 tetrahedra. The cube faces are flat, so all the mesh “curvature” goes inside, while the surface shows the triangular faces of the near-boundary tetrahedra; the cube edges are preserved as sharp features.

The figure “Ball” demonstrates surface approximation: the icosphere boundary of a ball ($R$ $= 1$, $n$ $= 800$; 737 vertices, 3417 tetrahedra) approximates the sphere with the prescribed accuracy, with no sharp edges.

Finally, the figure “Hollow Cylinder” shows a hollow cylinder ($R$ $= 1$, $r$ $= 0.4$, height $1.5$) with $n$ $= 1000$: 1118 vertices, 4588 tetrahedra. The through hole and the cap rims are preserved by the constraints: the sharp edges are detected by the dihedral angle and recovered in the mesh.