Voronoi Tessellations

DelaunayTriangulation.voronoiFunction
voronoi(tri::Triangulation; clip=false, smooth=false, kwargs...) -> VoronoiTessellation

Computes the Voronoi tessellation dual to a triangulation. If the triangulation is weighted, this computes the power diagram.

Arguments

  • tri: The triangulation.

Keyword Arguments

  • clip=false: If true, then the Voronoi tessellation is clipped to the convex hull of the triangulation. Otherwise, the Voronoi tessellation is unbounded.
  • clip_polygon=nothing: If clip=true, then this is the polygon to clip the Voronoi tessellation to. If nothing, the convex hull of the triangulation is used. The polygon should be defined as a Tuple of the form (points, boundary_nodes) where the boundary_nodes are vertices mapping to coordinates in points, adhering to the usual conventions for defining boundaries.
Convex

The clipping polygon must be convex. If it is not, unexpected results may occur.

  • smooth=false: If true, then the Voronoi tessellation is smoothed into a centroidal tessellation. Otherwise, the Voronoi tessellation is not smoothed. Must have clip=true if smooth=true.
  • predicates::AbstractPredicateKernel=AdaptiveKernel(): Method to use for computing predicates. Can be one of FastKernel, ExactKernel, and AdaptiveKernel. See the documentation for a further discussion of these methods.
  • kwargs...: Additional keyword arguments passed to centroidal_smooth if smooth=true.

Output

source
DelaunayTriangulation.centroidal_smoothFunction
centroidal_smooth(vorn::VoronoiTessellation; maxiters=1000, tol=default_displacement_tolerance(vorn), rng=Random.default_rng(), predicates::AbstractPredicateKernel=AdaptiveKernel(), kwargs...) -> VoronoiTessellation

Smooths vorn into a centroidal tessellation so that the new tessellation is of a set of generators whose associated Voronoi polygon is that polygon's centroid.

Arguments

Keyword Arguments

  • maxiters=1000: The maximum number of iterations.
  • clip_polygon=nothing: If clip=true, then this is the polygon to clip the Voronoi tessellation to. If nothing, the convex hull of the triangulation is used. The polygon should be defined as a Tuple of the form (points, boundary_nodes) where the boundary_nodes are vertices mapping to coordinates in points, adhering to the usual conventions for defining boundaries. Must be a convex polygon.
  • tol=default_displacement_tolerance(vorn): The displacement tolerance. See default_displacement_tolerance for the default.
  • rng=Random.default_rng(): The random number generator.
  • predicates::AbstractPredicateKernel=AdaptiveKernel(): Method to use for computing predicates. Can be one of FastKernel, ExactKernel, and AdaptiveKernel. See the documentation for a further discussion of these methods.
  • kwargs...: Extra keyword arguments passed to retriangulate.

Outputs

Extended help

The algorithm is simple. We iteratively smooth the generators, moving them to the centroid of their associated Voronoi polygon for the current tessellation, continuing until the maximum distance moved of any generator is less than tol. Boundary generators are not moved.

source
DelaunayTriangulation.get_generatorsFunction
get_generators(vorn::VoronoiTessellation) -> Dict{Vertex,Point}

Gets the generators of the Voronoi tessellation vorn as a Dict, mapping vertices to their coordinates. These coordinates are given as Tuples of the form (x, y).

source
DelaunayTriangulation.get_polygon_pointsFunction
get_polygon_points(vorn::VoronoiTessellation) -> Vector{Point}

Gets the polygon points of the Voronoi tessellation vorn. These are the vertices of the Voronoi polygons, and are given as Tuples of the form (x, y).

source
DelaunayTriangulation.get_polygonsFunction
get_polygons(vorn::VoronoiTessellation) -> Dict{Index,Vector{Vertex}}

Gets the polygons of the Voronoi tessellation vorn as a Dict, mapping polygon indices to their vertices, where the vertices refer to points in get_polygon_points(vorn). The vertices are given in counter-clockwise order and the first and last vertices are equal.

source
DelaunayTriangulation.get_circumcenter_to_triangleFunction
get_circumcenter_to_triangle(vorn::VoronoiTessellation) -> Dict{Index,Triangle}

Gets the circumcenters of the Voronoi tessellation vorn as a Dict, mapping circumcenter indices to their corresponding triangles. The triangles are sorted so that the minimum vertex is last.

source
get_circumcenter_to_triangle(vor::VoronoiTessellation, i) -> Triangle

Gets the triangle associated with the ith circumcenter. The triangle is sorted so that the minimum vertex is last.

source
DelaunayTriangulation.get_triangle_to_circumcenterFunction
get_triangle_to_circumcenter(vorn::VoronoiTessellation) -> Dict{Triangle,Index}

Gets the triangles of the Voronoi tessellation vorn as a Dict, mapping triangle indices to their corresponding circumcenters. The circumcenters are given as their vertices, referring to points in get_polygon_points(vorn).

source
get_triangle_to_circumcenter(vor::VoronoiTessellation, T) -> Index

Gets the circumcenter index associated with the triangle T. The triangle should be sorted so that the minimum vertex is last.

source
DelaunayTriangulation.get_cocircular_circumcentersFunction
get_cocircular_circumcenters(vorn::VoronoiTessellation) -> Set

Gets the cocircular circumcenters of the Voronoi tessellation vorn as a Set of circumcenter indices. These are circumcenters that come from triangles that are cocircular with another adjoining triangle.

source
DelaunayTriangulation.get_generatorFunction
get_generator(vor::VoronoiTessellation, i) -> NTuple{2, Number}
get_generator(vor::VoronoiTessellation, i...) -> NTuple{length(i), NTuple{2, Number}}

Gets the coordinates for the generators i..., returned as Tuples of the form (x, y) for each generator.

source
DelaunayTriangulation.get_polygon_pointFunction
get_polygon_point(vor::VoronoiTessellation, i) -> NTuple{2, Number}
get_polygon_point(vor::VoronoiTessellation, i...) -> NTuple{length(i), NTuple{2, Number}}

Gets the coordinates corresponding to the vertices i... of the polygons, returned as Tuples of the form (x, y) for each vertex.

source
DelaunayTriangulation.get_polygonFunction
get_polygon(vor::VoronoiTessellation, i) -> Vector{Vertex}

Gets the vector of vertices corresponding to the ith polygon, given in counter-clockwise order and with the first and last vertices equal. To obtain the coordinates, see get_polygon_point.

source
DelaunayTriangulation.num_polygon_verticesFunction
num_polygon_vertices(vor::VoronoiTessellation) -> Integer

Returns the number of polygon vertices in the Voronoi tessellation vor. This might include duplicate vertices if get_polygon_points(vor) has duplicates.

source
DelaunayTriangulation.get_centroidFunction
get_centroid(stats::TriangulationStatistics, T)

Returns the centroid field from the individual triangle statistics for the triangle T in the TriangulationStatistics stats.

source
get_centroid(vor::VoronoiTessellation, i) -> NTuple{2, Number}

Gets the centroid of the ith Voronoi polygon, given as a Tuple of the form (x, y).

source
DelaunayTriangulation.polygon_boundsFunction
polygon_bounds(vorn::VoronoiTessellation, unbounded_extension_factor=0.0; include_polygon_vertices=true) -> (Number, Number, Number, Number)

Gets the bounding box for the Voronoi tessellation vorn.

Arguments

  • vorn::VoronoiTessellation: The Voronoi tessellation.
  • unbounded_extension_factor=0.0: The factor by which to extend the bounding box for unbounded polygons.

Keyword Arguments

  • include_polygon_vertices=true: If true, then the bounding box will also include the polygon vertices. Otherwise, only the generators are included.

Output

  • xmin: Given by xmin′ - unbounded_extension_factor * (xmin′ - xmin′), where xmin′ is the original minimum x-coordinate of the computed bounding box and similarly for xmax′.
  • xmax: Given by xmax′ + unbounded_extension_factor * (xmax′ - xmax′), where xmax′ is the original maximum x-coordinate of the computed bounding box and similarly for xmin′.
  • ymin: Given by ymin′ - unbounded_extension_factor * (ymin′ - ymin′), where ymin′ is the original minimum y-coordinate of the computed bounding box and similarly for ymax′.
  • ymax: Given by ymax′ + unbounded_extension_factor * (ymax′ - ymax′), where ymax′ is the original maximum y-coordinate of the computed bounding box and similarly for ymin′.
source
polygon_bounds(points, boundary_nodes, check_all_curves=Val(false)) -> (Number, Number, Number, Number)

Computes the bounding box of the polygon defined by (points, boundary_nodes). The boundary_nodes must match the specification in the documentation and in check_args. If check_all_curves is true, then the bounding box of the union of all curves of the polygon is computed instead of just the first curve.

source
DelaunayTriangulation.get_polygon_coordinatesFunction
get_polygon_coordinates(vorn::VoronoiTessellation, i, bounding_box=nothing; predicates::AbstractPredicateKernel=AdaptiveKernel()) -> Vector{NTuple{2,Number}}

Returns the coordinates of the polygon with index i in vorn. If bounding_box is provided, then the polygon is clipped to the bounding box. If the polygon is unbounded, then bounding_box must be provided.

See also get_unbounded_polygon_coordinates and get_bounded_polygon_coordinates.

Arguments

  • vorn: The VoronoiTessellation.
  • i: The index of the polygon.
  • bounding_box=nothing: The bounding box to clip the polygon to. If nothing, then the polygon is not clipped. If the polygon is unbounded, then bounding_box must be provided.

Keyword Arguments

  • predicates::AbstractPredicateKernel=AdaptiveKernel(): Method to use for computing predicates. Can be one of FastKernel, ExactKernel, and AdaptiveKernel. See the documentation for a further discussion of these methods.

Outputs

  • coords: The coordinates of the polygon. This is a circular vector.
source
DelaunayTriangulation.polygon_featuresFunction
polygon_features(vor::VoronoiTessellation, i) -> (Number, NTuple{2, Number})

Gets the area and centroid of the polygon with index i in vor.

source
polygon_features(points, boundary_nodes) -> (Number, NTuple{2, Number})

Computes the signed area and centroid of the polygon defined by (points, boundary_nodes). The boundary_nodes must match the specification in the documentation and in check_args.

source