Mesh Generation

Triangulation

Triangulation is 2D is often solved using Shewchuk’s triangle library. It is both robust and flexible. We provide a pythonic wrapper over Shewchuk’s triangle that exposes most of its powers.

class pymesh.triangle

Wrapper around Shewchuk’s triangle.

points

numpy.ndarray – 3D or 2D points to be triangulated. If points are embedded in 3D, they must be coplanar.

segments

numpy.ndarray – n by 2 matrix of indices into points. together points and segments defines a Planar Straight Line Graph (PSLG) that triangles accepts.

triangles

numpy.ndarray – m by 3 matrix of indices into points. When segments is empty and triangles is non-empty, use triangle to refine the existing triangulation.

holes

numpy.ndarray – h by dim matrix of points representing hole points. Alternatively, one can set auto_hole_detection to True to infer holes from the input PSLG’s orientation.

min_angle

float – Lower bound on angle in degrees. Default is 20 degress. Setting min_angle > 20.7 will loose the theoretical guarentee of termination, although it often works fine in practice. However, setting min_angle > 34 tends to cause triangle to not terminate in practice.

max_area

float – Max triangle area. Default is unbounded.

max_areas

numpy.ndarray – Max area scalar field. It should have the same length as triangles. Not used by default.

keep_convex_hull

boolean – Whether to keep all triangles inside of the convex hull. Default is false.

conforming_delaunay

boolean – Whether to enforce conforming Delaunay triangulation. Default is false (i.e. use constrained Delaunay triangulation).

exact_arithmetic

boolean – Whether to use exact predicates. Defeault is true.

split_boundary

boolean – Whether to allow boundary to be split. Default is false.

max_num_steiner_points

int – The maximum number of Steiner points. Default is -1 (i.e. unbounded).

verbosity

int – How much info should triangle output?

  1. no output
  2. normal level of output
  3. verbose output
  4. vertex-by-vertex details
  5. you must be debugging the triangle code
algorithm

str – The Delaunay triangulation algorithm to use. Choices are:

  • DIVIDE_AND_CONQUER: Default. Implementation of [1].
  • SWEEPLINE: Fortune’s sweep line algorithm [2].
  • INCREMENTAL: Also from [1].
auto_hole_detection

boolean – Whether to detect holes based on the orientation of PSLG using winding number. Default is False.

vertices

numpy.ndarray – Vertices of the output triangulation.

faces

numpy.ndarray – Faces of the output triangulation.

mesh

Mesh – Output mesh.

voronoi_vertices

numpy.ndarray – Vertices of the output Voronoi diagram. Only generated when no input segments and triangles are provided.

voronoi_edges

numpy.ndarray – Voronoi edges. Negative index indicates infinity.

Example

>>> vertices = np.array([
...     [0.0, 0.0],
...     [1.0, 0.0],
...     [1.0, 1.0],
...     [0.0, 1.0],
...     ]);
>>> tri = pymesh.triangle();
>>> tri.points = vertices;
>>> tri.max_area = 0.05;
>>> tri.split_boundary = False;
>>> tri.verbosity = 0;
>>> tri.run(); # Execute triangle.
>>> mesh = tri.mesh; # output triangulation.

References:

[1](1, 2) Leonidas J. Guibas and Jorge Stolfi, Primitives for the Manipulation of General Subdivisions and the Computation of Voronoi Diagrams, ACM Transactions on Graphics 4(2):74-123, April 1985.
[2]Steven Fortune, A Sweepline Algorithm for Voronoi Diagrams, Algorithmica 2(2):153-174, 1987.

Tetrahedralization

In contrast with 2D, tetrahedralization in 3D is a much hard problem. Many algorithms tries to tackle this problem from different angles. No single algorithm or package standouts as the best. We therefore offer a number of different engines for our users.

pymesh.tetrahedralize(mesh, cell_size, radius_edge_ratio=2.0, facet_distance=-1.0, feature_angle=120, engine='auto', with_timing=False)

Create a tetrahedral mesh from input triangle mesh.

Parameters:
  • mesh (Mesh) – Input triangular mesh.
  • cell_size (float) – Max radius of the circumscribed sphere of the output tet.
  • radius_edge_ratio (float) – Max radius of the circumscribed sphere to the shortest edge length of each tet.
  • facet_distance (float) – Upper bound on the distance from the circumcenter of a facet to the center of its “Delaunay ball”, where a Delaunay ball is defined as the smallest circumscribed sphere with center on the surface of the domain.
  • feature_angle (float) – Angle threshold (in degrees) for feature extraction.
  • engine (string) –

    The tetrahedralization engine to use. Valid options are:

    • auto: default to tetgen
    • cgal: CGAL 3D mesh generation, using Polyhedron domain with auto feature extraction.
    • cgal_no_features: CGAL 3D mesh generation, using Polyhedron domain without feature extraction.
    • cgal_implicit: CGAL 3D mesh generation, using implicit domain with winding number as oracle.
    • tetgen: TetGen from Hang Si.
    • quartet: Quartet from Robert Bridson and Crawford Doran
    • delpsc: DelPSC from Tamal K Dey , Joshua A. Levine, Andrew Slatton
    • vegafem: Tet mesher provided by VegaFEM library.
    • mmg: Implicit domain meshing from MMG3D.
    • tetwild: TetWild engine based on our Siggraph paper.
  • with_timing (boolean) – whether to output timing info.
Returns:

Tetrahedral mesh (and running time if with_timing is True).

In addition to pymesh.tetraheralize(), we also provide a more complete wrapper around Si’s awesome TetGen package.

class pymesh.tetgen

Wrapper around Si’s TetGen.

All attributes, except vertices, faces, voxels and mesh, are either input geometry or configuration parameters.

points

numpy.ndarray – n by 3 list of points to be tetrahedralized.

triangles

numpy.ndarray – m by 3 matrix of indices into points. Together, points and triangles defined PLC.

tetrhaedra

numpy.ndarray – t by 4 matrix of indices into points. Used for refinement.

point_markers

numpy.ndarray – List of integer point markers of size n. Point marker cannot be 0.

point_weights

numpy.ndarray – List of point weights. Used for weight Delaunay tetrahedralization.

triangle_marker

numpy.ndarray – List of integer triangle markers of size t.

split_boundary

bool – whether to split input boundary. Default is true.

max_radius_edge_ratio

float – Default is 2.0.

min_dihedral_angle

float – Default is 0.0.

coarsening

bool – Coarsening the input tet mesh. Default is false.

max_tet_volume

float – Default is unbounded.

optimization_level

int – Ranges from 0 to 10. Default is 2.

max_num_steiner_points

int – Default is unbounded.

coplanar_tolerance

float – Used for determine when 4 points are coplanar. Default is 1e-8.

exact_arithmetic

bool – Whether to use exact predicates. Default is true.

merge_coplanar

bool – Whether to merge coplanar faces and nearby vertices. Default is true.

weighted_delaunay

bool – Compute weighted Delaunay tetrahedralization instead of conforming Delaunay. Default is false. This option requires point_weights.

keep_convex_hull

bool – Keep all tets within the convex hull. Default is false.

verbosity

int – Verbosity level. Ranges from 0 to 4:

  1. no output
  2. normal level of output
  3. verbose output
  4. more details
  5. you must be debugging the tetgen code
vertices

numpy.ndarray – Vertices of the output tet mesh.

faces

numpy.ndarray – Faces of the output tet mesh.

voxels

numpy.ndarray – Voxels of the output tet mesh.

mesh

Mesh – Output tet mesh.

Example

>>> input_mesh = pymesh.generate_icosphere(1.0, [0.0, 0.0, 0.0]);
>>> tetgen = pymesh.tetgen();
>>> tetgen.points = input_mesh.vertices; # Input points.
>>> tetgen.triangles = input_mesh.faces; # Input triangles
>>> tetgen.max_tet_volume = 0.01;
>>> tetgen.verbosity = 0;
>>> tetgen.run(); # Execute tetgen
>>> mesh = tetgen.mesh; # Extract output tetrahedral mesh.