Triangulations
DelaunayTriangulation.triangulate
— Functiontriangulate(points; segments=nothing, boundary_nodes=nothing, kwargs...) -> Triangulation
Computes the Delaunay triangulation of points
, and then the constrained Delaunay triangulation if any of segments
and boundary_nodes
are not nothing
.
Arguments
points
: The points to triangulate. (This might get mutated for curve-bounded domains.)
Keyword Arguments
For the keyword arguments below, you may like to review the extended help as some of the arguments carry certain warnings.
segments=nothing
: The segments to include in the triangulation. Ifnothing
, then no segments are included.boundary_nodes=nothing
: The boundary nodes to include in the triangulation. Ifnothing
, then no boundary nodes are included, and the convex hull ofpoints
remains as the triangulation. These boundary nodes should match the specification given incheck_args
if a boundary is provided as a set of vertices, meaning the boundary is a piecewise linear curve. To specify a curve-bounded domain, you should follow the same specification, but useAbstractParametricCurve
s to fill out the vector, and any piecewise linear section should still be provided as a sequence of vertices.predicates::AbstractPredicateKernel=AdaptiveKernel()
: Method to use for computing predicates. Can be one ofFastKernel
,ExactKernel
, andAdaptiveKernel
. See the documentation for a further discussion of these methods.weights=ZeroWeight{number_type(points)}()
: The weights to use for the triangulation. By default, the triangulation is unweighted. The weights can also be provided as a vector, with thei
th weight referring to thei
th vertex, or more generally any object that definesget_weight
.IntegerType=Int
: The integer type to use for the triangulation. This is used for representing vertices.EdgeType=isnothing(segments) ? NTuple{2,IntegerType} : (edge_type ∘ typeof)(segments)
: The edge type to use for the triangulation.TriangleType=NTuple{3,IntegerType}
: The triangle type to use for the triangulation.EdgesType=isnothing(segments) ? Set{EdgeType} : typeof(segments)
: The type to use for storing the edges of the triangulation.TrianglesType=Set{TriangleType}
: The type to use for storing the triangles of the triangulation.randomise=true
: Whether to randomise the order in which the points are inserted into the triangulation. This is done usingget_insertion_order
.delete_ghosts=false
: Whether to delete the ghost triangles after the triangulation is computed. This is done usingdelete_ghost_triangles!
.delete_empty_features=true
: Whether to delete empty features after the triangulation is computed. This is done usingclear_empty_features!
.try_last_inserted_point=true
: Whether to try the last inserted point first when inserting points into the triangulation.skip_points=()
: The points to skip when inserting points into the triangulation. Note that, for curve-bounded domains,skip_points
is ignored when usingenrich_boundary!
.num_sample_rule=default_num_samples
: A function mapping a number of pointsn
to a number of samplesm
to use for sampling the initial points during the point location step of the algorithm withinfind_triangle
.rng::Random.AbstractRNG=Random.default_rng()
: The random number generator.insertion_order::Vector=get_insertion_order(points, randomise, skip_points, IntegerType, rng)
: The insertion order to use for inserting points into the triangulation. This is ignored if you are defining a curve-bounded domain.recompute_representative_points=true
: Whether to recompute the representative points after the triangulation is computed. This is done usingcompute_representative_points!
.delete_holes=true
: Whether to delete holes after the triangulation is computed. This is done usingdelete_holes!
.check_arguments=true
: Whether to check the argumentspoints
andboundary_nodes
are valid. This is done usingcheck_args
.polygonise_n=4096
: Number of points to use for polygonising the boundary when considering the polygon hierarchy for a curve-bounded domain usingpolygonise
. Seetriangulate_curve_bounded
.coarse_n=0
: Number of points to use for initialising a curve-bounded domain. Seetriangulate_curve_bounded
. (A value of0
means the number of points is chosen automatically until the diametral circles of all edges are empty.)
Outputs
tri::Triangulation
: The triangulation.
The output from this function is currently not type stable. In particular, the inferred type is only Triangulation
without any other information. If you are depending on the output from triangulate
inside some other function, you should consider putting the output behind a function barrier; information about using function barriers is given here.
Extended help
Here are some warnings to consider for some of the arguments.
points
For curve-bounded domains, points
may get mutated to include the endpoints of the provided curves, and when inserting Steiner points to split segments or refine boundaries.
If your points are defined using non-Float64
coordinates, you may run into precision issues that lead to problems with robustness. The consequences of this could be potentially catastrophic, leading to infinite loops for example. If you do encounter such issues, consider converting your coordinates to Float64
.
segments
When segments are outside of the domain, are if they are not entirely contained with the domain, you may run into issues - especially for curve-bounded domains. It is your responsibility to ensure that the segments are contained within the domain.
The segments
may get mutated in two ways: (1) Segments may get rotated so that (i, j)
becomes (j, i)
. (2) If there are segments that are collinear with other segments, then they may get split into chain of non-overlapping connecting segments (also see below). For curve-bounded domains, segments are also split so that no subsegment's diametral circle contains any other point.
Currently, segments that intersect in their interiors (this excludes segments that only intersect by sharing a vertex) cause problems for triangulating. While there is some support for collinear segments that lie on top of each other (they get split automatically), this is not the case for segments that intersect in their interiors. Moreover, this automatic splitting should not be heavily relied upon, and for curve-bounded domains you should not rely on it at all as it causes problems during the enrichment phase from enrich_boundary!
.
boundary_nodes
While for standard domains with piecewise linear boundaries (or no boundaries) it is fine for points to be outside of the domain (they just get automatically deleted if needed), they may cause problems for curve-bounded domains. Please ensure that all your points are inside the curve-bounded domain if you are providing curves in boundary_nodes
.
For curve-bounded domains, the boundary_nodes
in the resulting Triangulation
will not be aliased with the input boundary nodes.
For curve-bounded domains, note that the triangulation produced from this function is really just an initial coarse discretisation of the true curved boundaries. You will need to refine further, via refine!
, to improve the discretisation, or increase coarse_n
below. See also polygonise
for a more direct approach to discretising a boundary (which might not give as high-quality meshes as you can obtain from refine!
though, note).
DelaunayTriangulation.triangulate_convex
— Functiontriangulate_convex(points, S; delete_ghosts=false, delete_empty_features=true, rng=Random.default_rng(), predicates::AbstractPredicateKernel=AdaptiveKernel(), kwargs...) -> Triangulation
Triangulates the convex polygon S
.
Arguments
points
: The point set corresponding to the vertices inS
.S
: A convex polygon represented as a vector of vertices. The vertices should be given in counter-clockwise order, and must not be circular so thatS[begin] ≠ S[end]
.
Keyword Arguments
delete_ghosts=false
: Iftrue
, the ghost triangles are deleted after triangulation.delete_empty_features=true
: Iftrue
, the empty features are deleted after triangulation.rng=Random.default_rng()
: The random number generator used to shuffle the vertices ofS
before triangulation.predicates::AbstractPredicateKernel=AdaptiveKernel()
: Method to use for computing predicates. Can be one ofFastKernel
,ExactKernel
, andAdaptiveKernel
. See the documentation for a further discussion of these methods.kwargs...
: Additional keyword arguments passed toTriangulation
.
Output
tri::Triangulation
: The triangulated polygon.
DelaunayTriangulation.triangulate_rectangle
— Functiontriangulate_rectangle(a, b, c, d, nx, ny; kwargs...) -> Triangulation
Triangulates the rectangle [a, b] × [c, d]
.
Arguments
a
: The minimumx
-coordinate.b
: The maximumx
-coordinate.c
: The minimumy
-coordinate.d
: The maximumy
-coordinate.nx
: The number of points in thex
-direction.ny
: The number of points in they
-direction.
Keyword Arguments
single_boundary=false
: Iftrue
, then the boundary nodes are stored as a contiguous section. Otherwise, the boundary is split into four sections, in the order bottom, right, top, left.delete_ghosts=false
: Iftrue
, then the ghost triangles are deleted. Otherwise, they are kept.IntegerType::Type{I}=Int
: The type of the vertices.EdgeType::Type{E}=NTuple{2,IntegerType}
: The type of the edges.TriangleType::Type{V}=NTuple{3,IntegerType}
: The type of the triangles.EdgesType::Type{Es}=Set{EdgeType}
: The type of the edges container.TrianglesType::Type{Ts}=Set{TriangleType}
: The type of the triangles container.predicates::AbstractPredicateKernel=AdaptiveKernel()
: Method to use for computing predicates. Can be one ofFastKernel
,ExactKernel
, andAdaptiveKernel
. See the documentation for a further discussion of these methods.
Outputs
tri
: The triangulation of the rectangle.
DelaunayTriangulation.refine!
— Functionrefine!(tri::Triangulation; kwargs...)
Refines the given Triangulation
tri
to meet the given quality constraints.
See the documentation for more information about mesh refinement, e.g. convergence issues and issues with small input-angles.
Arguments
tri::Triangulation
: TheTriangulation
to refine.
Keyword Arguments
min_angle=30.0
: The minimum angle constraint, in degrees.max_angle=180.0
: The maximum angle constraint, in degrees.
Maximum angle constraints are not currently implemented.
min_area=get_area(tri) / 1e9
: The minimum area constraint.max_area=typemax(number_type(tri))
: The maximum area constraint.max_points=max(1_000, num_solid_vertices(tri))^2
: The maximum number of vertices allowed in the triangulation. Note that this refers tonum_solid_vertices
, not the amount returned bynum_points
.seditious_angle=20.0
: The angle at which a triangle is considered seditious, in degrees. Seeis_triangle_seditious
.custom_constraint=(tri, T) -> false
: A custom constraint function that takes aTriangulation
and a triangle, and returnstrue
if the triangle should be refined andfalse
otherwise.use_circumcenter=true
: Whether to insert circumcenters for refining a triangle or generalised Steiner points.
Generalised Steiner points are not yet implemented. Thus, this argument must be true
(and the steiner_scale
keyword below is ignored).
use_lens=true
: Whether to use the diametral lens or the diametral circle for checking encroachment.steiner_scale=0.999
: The perturbation factor to use for generalised Steiner points ifuse_circumcenter=false
. (Not currently used - see above.)rng=Random.default_rng()
: The random number generator to use in case it is needed during point location.concavity_protection=false
: Whether to use concavity protection or not forfind_triangle
. Most likely not needed, but may help in pathological cases.predicates::AbstractPredicateKernel=AdaptiveKernel()
: Method to use for computing predicates. Can be one ofFastKernel
,ExactKernel
, andAdaptiveKernel
. See the documentation for a further discussion of these methods.
Output
The triangulation is refined in-place.
During refinement, points are often deleted, which may often lead to points in get_points(tri)
that do not appear anywhere in the triangulation. (This is why we recommend e.g. each_solid_vertex
over each_point
.) Similarly, since points are deleted, when two triangles have a common circumcenter it might happen (if they are near an input segment) that a point is duplicated inside get_points(tri)
, in case one circumcenter was deleted previously.
DelaunayTriangulation.retriangulate
— Functionretriangulate(tri::Triangulation, points=get_points(tri); kwargs...)
Retriangulates the triangulation tri
using the points points
, returning a new Triangulation
.
Arguments
tri::Triangulation
: The triangulation to retriangulate.points=get_points(tri)
: The points to use for retriangulating the triangulation. By default, this is simplyget_points(tri)
.
Keyword Arguments
skip_points=Set(filter(i -> !has_vertex(tri, i), each_point_index(tri)))
: The points to skip when inserting points into the triangulation.kwargs...
: Extra keyword arguments passed totriangulate
. Other keyword arguments, likesegments
andboundary_nodes
, are automatically passed from the fields oftri
, but may be overridden by passing the corresponding keyword arguments.
DelaunayTriangulation.convert_boundary_points_to_indices
— Functionconvert_boundary_points_to_indices(x, y; existing_points = NTuple{2, Float64}[], check_args=true) -> (Vector, Vector)
convert_boundary_points_to_indices(xy; existing_points = NTuple{2, Float64}[], check_args=true) -> (Vector, Vector)
Converts a boundary represented by (x, y)
or xy
, where the points are combined rather than as separate sets of coordinates, into a set of boundary nodes for use in triangulate
.
Arguments
x
,y
: Thex
andy
-coordinates for the boundary points. The individual vectors must match the specification required for boundaries outlined in the documentation.xy
: As above, except the coordinates are combined rather than given as separate vectors.
Keyword Arguments
existing_points
: The existing points to append the boundary points to. This is useful if you have a pre-existing set of points.check_args
: Whether to check that the arguments match the specification in the documentation.
Outputs
boundary_nodes
: The boundary nodes.points
: The point set, which is the same asexisting_points
but with the boundary points appended to it.
DelaunayTriangulation.check_args
— Functioncheck_args(points, boundary_nodes, hierarchy::PolygonHierarchy, boundary_curves = ()) -> Bool
Check that the arguments points
and boundary_nodes
to triangulate
, and a constructed PolygonHierarchy
given by hierarchy
, are valid. In particular, the function checks:
- The dimension of the points. If the points are not 2D, a warning is thrown.
- The points are all unique. If they are not, a
DuplicatePointsError
is thrown. - There are at least three points. If there are not, an
InsufficientPointsError
is thrown.
If boundary_nodes
are provided, meaning has_boundary_nodes
, then the function also checks:
- If the boundary curves all connect consistently. Here, this means that each section of a boundary curve ends at the start of the next boundary section; for contiguous boundary curves, this means that the start and end boundary nodes are the same.
- If the orientation of the boundary curves are all consistent. This means that the curves are all positively oriented relative to the domain, so that e.g. the exterior boundary curves are all counter-clockwise (relative to just themselves), the next exterior-most curves inside those exteriors are all clockwise (again, relative to just themselves), and so on.
Another requirement for triangulate
is that none of the boundaries intersect in their interior, which also prohibits interior self-intersections. This is NOT checked. Similarly, segments should not intersect in their interior, which is not checked.
DelaunayTriangulation.get_points
— Functionget_points(convex_hull::ConvexHull) -> Points
Returns the underlying point set of convex_hull
.
get_points(tri::Triangulation) -> Points
Return the points of the triangulation. Note that this may include points not yet in tri
.
get_points(boundary_enricher::BoundaryEnricher{P}) -> P
Returns the point set associated with boundary_enricher
.
DelaunayTriangulation.get_triangles
— Functionget_triangles(tri::Triangulation) -> Triangles
Return the triangles of the triangulation. These triangles are all given in counter-clockwise order, and may include ghost triangles.
DelaunayTriangulation.get_boundary_nodes
— Functionget_boundary_nodes(boundary_nodes, mnℓ...)
Given a collection of boundary_nodes
, returns the specified component of the collection. There are several forms for the methods:
get_boundary_nodes(boundary_nodes, m)
: Ifboundary_nodes
has multiple curves, this returns them
th curve. Ifboundary_nodes
has multiple sections, this returns them
th section. Otherwise, this returns them
th boundary node.get_boundary_nodes(boundary_nodes, m, n)
: Ifboundary_nodes
has multiple curves, this returns then
th section of them
th curve. Otherwise, ifboundary_nodes
has multiple sections, this returns then
th boundary node of them
th section.get_boundary_nodes(boundary_nodes, (m, n))
: This is equivalent toget_boundary_nodes(boundary_nodes, m, n)
.get_boundary_nodes(boundary_nodes::A, ::A)
: This just returnsboundary_nodes
.
Examples
julia> using DelaunayTriangulation
julia> get_boundary_nodes([[[1, 2, 3, 4], [4, 5, 1]], [[6, 7, 8, 9], [9, 10, 6]]], 2)
2-element Vector{Vector{Int64}}:
[6, 7, 8, 9]
[9, 10, 6]
julia> get_boundary_nodes([[1, 2, 3, 4], [4, 5, 1]], 1)
4-element Vector{Int64}:
1
2
3
4
julia> get_boundary_nodes([1, 2, 3, 4, 5, 6, 1], 4)
4
julia> get_boundary_nodes([[[1, 2, 3, 4], [4, 5, 1]], [[6, 7, 8, 9], [9, 10, 6]]], 1, 2)
3-element Vector{Int64}:
4
5
1
julia> get_boundary_nodes([[1, 2, 3, 4], [4, 5, 6, 1]], 2, 3)
6
julia> get_boundary_nodes([1, 2, 3, 4, 5, 1], [1, 2, 3, 4, 5, 1])
6-element Vector{Int64}:
1
2
3
4
5
1
DelaunayTriangulation.get_interior_segments
— Functionget_interior_segments(tri::Triangulation) -> Edges
Return the interior segments of the triangulation. These are segments that are forced to be in the triangulation - they are not the same as edges.
DelaunayTriangulation.get_all_segments
— Functionget_all_segments(tri::Triangulation) -> Edges
Return all segments of the triangulation. This includes interior segments and boundary segments. Segments are edges that are forced to be in the triangulation.
DelaunayTriangulation.get_weights
— Functionget_weights(tri::Triangulation) -> Weights
Return the weights of the triangulation. These are the weights of the vertices of the triangulation.
DelaunayTriangulation.get_adjacent
— Functionget_adjacent(adj::Adjacent) -> Dict
Returns the adjacent
map of adj
.
Examples
julia> using DelaunayTriangulation
julia> d = Dict((1, 2) => 3, (2, 3) => 1, (3, 1) => 2);
julia> adj = DelaunayTriangulation.Adjacent(d)
Adjacent{Int64, Tuple{Int64, Int64}}, with map:
Dict{Tuple{Int64, Int64}, Int64} with 3 entries:
(1, 2) => 3
(3, 1) => 2
(2, 3) => 1
julia> get_adjacent(adj)
Dict{Tuple{Int64, Int64}, Int64} with 3 entries:
(1, 2) => 3
(3, 1) => 2
(2, 3) => 1
julia> get_adjacent(adj) == d
true
get_adjacent(adj::Adjacent{I, E}, uv::E) -> Vertex
get_adjacent(adj::Adjacent{I, E}, u, v) -> Vertex
Returns the vertex w
such that (u, v, w)
is a positively oriented triangle in the underlying triangulation, or ∅
if no such triangle exists.
Examples
julia> using DelaunayTriangulation
julia> adj = DelaunayTriangulation.Adjacent(Dict((1, 2) => 3, (2, 3) => 1, (3, 1) => 2, (4, 5) => -1))
Adjacent{Int64, Tuple{Int64, Int64}}, with map:
Dict{Tuple{Int64, Int64}, Int64} with 4 entries:
(4, 5) => -1
(1, 2) => 3
(3, 1) => 2
(2, 3) => 1
julia> get_adjacent(adj, 4, 5)
-1
julia> get_adjacent(adj, (3, 1))
2
julia> get_adjacent(adj, (1, 2))
3
julia> get_adjacent(adj, 17, 5)
0
julia> get_adjacent(adj, (1, 6))
0
get_adjacent(tri::Triangulation) -> Adjacent
Returns the adjacency map of the triangulation. This is a map from each edge (u, v)
to a vertex w
such that (u, v, w)
is a positively oriented triangle in tri
.
See also Adjacent
.
get_adjacent(tri::Triangulation, uv) -> Vertex
get_adjacent(tri::Triangulation, u, v) -> Vertex
Returns the vertex w
such that (u, v, w)
is a positively oriented triangle in the underlying triangulation, or ∅
if no such triangle exists.
get_adjacent(vorn::VoronoiTessellation) -> Adjacent{Index,Edge}
Gets the adjacency information of the Voronoi tessellation vorn
as an Adjacent
object. This object maps oriented edges to the polygons that they belong to.
get_adjacent(vor::VoronoiTessellation, ij) -> Index
get_adjacent(vor::VoronoiTessellation, i, j) -> Index
Gets the polygon index associated with the oriented edge ij
in the Voronoi tessellation vor
.
DelaunayTriangulation.get_adjacent2vertex
— Functionget_adjacent2vertex(adj2v::Adjacent2Vertex) -> Dict
Returns the adjacent2vertex
map of adj2v
.
Examples
julia> using DelaunayTriangulation
julia> e1 = Set(((1, 2), (5, 3), (7, 8)));
julia> e2 = Set(((2, 3), (13, 5), (-1, 7)));
julia> d = Dict(9 => e1, 6 => e2);
julia> adj2v = DelaunayTriangulation.Adjacent2Vertex(d)
Adjacent2Vertex{Int64, Set{Tuple{Int64, Int64}}} with map:
Dict{Int64, Set{Tuple{Int64, Int64}}} with 2 entries:
6 => Set([(13, 5), (-1, 7), (2, 3)])
9 => Set([(1, 2), (7, 8), (5, 3)])
julia> get_adjacent2vertex(adj2v)
Dict{Int64, Set{Tuple{Int64, Int64}}} with 2 entries:
6 => Set([(13, 5), (-1, 7), (2, 3)])
9 => Set([(1, 2), (7, 8), (5, 3)])
julia> get_adjacent2vertex(adj2v) == d
true
get_adjacent2vertex(adj2v::Adjacent2Vertex, w) -> Edges
Returns the set of edges E
such that (u, v, w)
is a positively oriented triangle in the underlying triangulation for each (u, v) ∈ E
.
Examples
julia> using DelaunayTriangulation
julia> adj2v = DelaunayTriangulation.Adjacent2Vertex(Dict(1 => Set(((2, 3), (5, 7), (8, 9))), 5 => Set(((1, 2), (7, 9), (8, 3)))))
Adjacent2Vertex{Int64, Set{Tuple{Int64, Int64}}} with map:
Dict{Int64, Set{Tuple{Int64, Int64}}} with 2 entries:
5 => Set([(1, 2), (8, 3), (7, 9)])
1 => Set([(8, 9), (5, 7), (2, 3)])
julia> get_adjacent2vertex(adj2v, 1)
Set{Tuple{Int64, Int64}} with 3 elements:
(8, 9)
(5, 7)
(2, 3)
julia> get_adjacent2vertex(adj2v, 5)
Set{Tuple{Int64, Int64}} with 3 elements:
(1, 2)
(8, 3)
(7, 9)
get_adjacent2vertex(tri::Triangulation) -> Adjacent2Vertex
Returns the Adjacent2Vertex
map of the triangulation tri
. This is a map from a vertex w
to a set of all edges (u, v)
such that (u, v, w)
is a positively oriented triangle in tri
.
get_adjacent2vertex(tri::Triangulation, w) -> Edges
Returns the set of all edges (u, v)
in tri
such that (u, v, w)
is a positively oriented triangle in tri
.
DelaunayTriangulation.get_graph
— Functionget_graph(tri::Triangulation) -> Graph
Returns the Graph
of the triangulation tri
. This is an undirected graph.
DelaunayTriangulation.get_neighbours
— Functionget_neighbours(graph::Graph) -> Dict{Vertex, Set{Vertex}}
Returns the neighbours
map of graph
.
get_neighbours(G::Graph, u) -> Set{Vertex}
Returns the set of neighbours of u
in G
.
get_neighbours(tri::Triangulation) -> Dict{Vertex, Set{Vertex}}
Returns the neighbours
map of tri
. Note that, if has_ghost_triangles(tri)
, then some of the neighbours and vertices will be ghost vertices.
get_neighbours(tri::Triangulation, u) -> Set{Vertex}
Returns the set of neighbours of u
in tri
. Note that, if has_ghost_triangles(tri)
, then some of the neighbours and vertices will be ghost vertices.
DelaunayTriangulation.get_boundary_curves
— Functionget_boundary_curves(tri::Triangulation) -> NTuple{N, Function}
Returns the functions defining the boundaries of tri
. If !is_curve_bounded(tri)
, then this returns an empty Tuple
. Otherwise, this returns a Tuple
of functions, one for each section of the boundary, where the i
th element of the Tuple
corresponds to the i
th section of the boundary, which corresponds to the ghost vertex -i
. For curves that are defined by boundary nodes rather than by a function, the function is PiecewiseLinear
. For the other functions, these are all defined by t -> NTuple{2, Number}
, where t ∈ [0, 1]
and the NTuple{2, Number}
is the coordinate on the curve at that t
.
get_boundary_curves(boundary_enricher::BoundaryEnricher{P,B,C}) -> C
Returns the boundary curves associated with boundary_enricher
.
DelaunayTriangulation.get_boundary_edge_map
— Functionget_boundary_edge_map(tri::Triangulation) -> Dict
Returns the boundary edge map of the triangulation tri
. This is a Dict
that maps a boundary edge (u, v)
to its position in get_boundary_nodes(tri)
. In particular, the returned value is a Tuple
(position, index)
so that boundary_nodes = get_boundary_nodes(tri, position)
are the boundary nodes associated with the section that (u, v)
resides on, and u = get_boundary_nodes(boundary_nodes, index)
and v = get_boundary_nodes(boundary_nodes, index + 1)
.
See also construct_boundary_edge_map
.
get_boundary_edge_map(boundary_enricher::BoundaryEnricher{P,B,C,I,BM}) -> BM
Returns the boundary edge map associated with boundary_enricher
.
get_boundary_edge_map(boundary_enricher::BoundaryEnricher, i, j)
Returns the value from the key (i, j)
in the boundary edge map of boundary_enricher
. The returned value is a Tuple
(position, index)
so that boundary_nodes = get_boundary_nodes(get_boundary_nodes(boundary_enricher), position)
are the boundary nodes associated with the section that (i, j)
resides on, and i = get_boundary_nodes(boundary_nodes, index)
and j = get_boundary_nodes(boundary_nodes, index + 1)
.
get_boundary_edge_map(tri::Triangulation, ij)
get_boundary_edge_map(tri::Triangulation, i, j)
Returns the value from the key (i, j)
in the boundary edge map of tri
. The returned value is a Tuple
(position, index)
so that boundary_nodes = get_boundary_nodes(tri, position)
are the boundary nodes associated with the section that (i, j)
resides on, and i = get_boundary_nodes(boundary_nodes, index)
and j = get_boundary_nodes(boundary_nodes, index + 1)
.
DelaunayTriangulation.get_ghost_vertex_map
— Functionget_ghost_vertex_map(tri::Triangulation) -> Dict
Returns the ghost vertex map of the triangulation tri
. This is a Dict
that maps ghost vertices to their associated section in boundary_nodes
. There are three cases; below, I
is integer_type(tri)
:
has_multiple_curves(tri)
Returns dict::Dict{I, NTuple{2, I}}
, mapping ghost vertices i
to Tuple
s (m, n)
so that get_boundary_nodes(tri, m, n)
are the boundary nodes associated with i
, i.e. the n
th section of the m
th curve is associated with the ghost vertex i
.
has_multiple_sections(tri)
Returns dict::Dict{I, I}
, mapping ghost vertices i
to n
so that get_boundary_nodes(tri, n)
are the boundary nodes associated with i
, i.e. the n
th section of the boundary is associated with the ghost vertex i
.
otherwise
Returns dict::Dict{I, A}
, mapping the ghost vertex i
to get_boundary_nodes(tri)
, where A = typeof(get_boundary_nodes(tri))
.
See also construct_ghost_vertex_map
.
DelaunayTriangulation.get_ghost_vertex_ranges
— Functionget_ghost_vertex_ranges(tri::Triangulation) -> Dict
Returns the ghost vertex ranges map of the triangulation tri
. This is a Dict
that maps ghost vertices i
to the range of all other ghost vertices associated with the curve that i
is associated with.
See also construct_ghost_vertex_ranges
.
DelaunayTriangulation.get_convex_hull
— Functionget_convex_hull(tri::Triangulation) -> ConvexHull
Returns the convex hull of the points in tri
. This is given as a ConvexHull
object, where the vertices are sorted counter-clockwise and defined so that the first and last vertices are equal.
DelaunayTriangulation.get_representative_point_list
— Functionget_representative_point_list(tri::Triangulation) -> Dict
Returns the Dict
of RepresentativeCoordinates
of tri
, mapping curve indices i
to the representative point for that curve. These representative points are how we interpret ghost triangles relative to that curve.
DelaunayTriangulation.num_curves
— Functionnum_curves(tri::Triangulation) -> Integer
Returns the number of curves in tri
.
num_curves(boundary_nodes) -> Integer
Get the number of curves in boundary_nodes
.
Examples
julia> using DelaunayTriangulation
julia> DelaunayTriangulation.num_curves([1, 2, 3, 1])
1
julia> DelaunayTriangulation.num_curves([[1, 2, 3], [3, 4, 1]])
1
julia> DelaunayTriangulation.num_curves([[[1, 2, 3], [3, 4, 1]], [[5, 6, 7, 8, 5]]])
2
DelaunayTriangulation.num_sections
— Functionnum_sections(tri::Triangulation) -> Integer
Assuming tri
only has one curve, returns the number of sections in tri
.
num_sections(boundary_nodes) -> Integer
Assuming boundary_nodes
has only one curve, get the number of sections in boundary_nodes
.
Examples
julia> using DelaunayTriangulation
julia> DelaunayTriangulation.num_sections([1, 2, 3, 4, 5, 1])
1
julia> DelaunayTriangulation.num_sections([[1, 2, 3, 4], [4, 5, 1]])
2
julia> DelaunayTriangulation.num_sections([[1, 2, 3], [3, 4, 5, 6, 7, 8], [8, 9], [9, 1]])
4
DelaunayTriangulation.get_right_boundary_node
— Functionget_right_boundary_node(tri::Triangulation, k, ghost_vertex) -> Vertex
Returns the boundary node to the right of the vertex k
in tri
.
Arguments
tri::Triangulation
: TheTriangulation
.k
: The boundary vertex.ghost_vertex
: The ghost vertex associated with the boundary section thatk
is on.
Outputs
r
: The vertex right ofk
on the boundary.
DelaunayTriangulation.get_left_boundary_node
— Functionget_left_boundary_node(tri::Triangulation, k, ghost_vertex) -> Vertex
Returns the boundary node to the left of the vertex k
in tri
.
Arguments
tri::Triangulation
: TheTriangulation
.k
: The boundary vertex.ghost_vertex
: The ghost vertex associated with the boundary section thatk
is on.
Outputs
ℓ
: The vertex left ofk
on the boundary.
DelaunayTriangulation.get_all_boundary_nodes
— Functionget_all_boundary_nodes(tri::Triangulation) -> Set{Vertex}
Returns the set of all boundary vertices in tri
, in no specific order.
DelaunayTriangulation.contains_boundary_edge
— Functioncontains_boundary_edge(tri::Triangulation, ij) -> Bool
contains_boundary_edge(tri::Triangulation, i, j) -> Bool
Returns true
if the boundary edge (i, j)
is in tri
, and false
otherwise. Orientation matters here.
DelaunayTriangulation.has_vertex
— Functionhas_vertex(G::Graph, u) -> Bool
Returns true
if u
is a vertex in G
, and false
otherwise.
has_vertex(tri::Triangulation, u) -> Bool
Returns true
if u
is a vertex in tri
, and false
otherwise.
DelaunayTriangulation.delete_ghost_vertices_from_graph!
— Functiondelete_ghost_vertices!(G::Graph)
Deletes all ghost vertices from G
.
delete_ghost_vertices_from_graph!(tri::Triangulation)
Deletes all ghost vertices from the graph of tri
.
DelaunayTriangulation.has_ghost_vertices
— Functionhas_ghost_vertices(G::Graph) -> Bool
Returns true
if G
has ghost vertices, and false
otherwise.
has_ghost_vertices(tri::Triangulation) -> Bool
Returns true
if tri
has ghost vertices, and false
otherwise.
DelaunayTriangulation.get_convex_hull_vertices
— Functionget_convex_hull_vertices(tri::Triangulation) -> Vector{Vertex}
Returns the vertices on the convex hull of tri
, in counter-clockwise order.
See also ConvexHull
.
DelaunayTriangulation.is_exterior_ghost_vertex
— Functionis_exterior_ghost_vertex(tri::Triangulation, i) -> Bool
Returns true
if the ghost vertex i
in tri
is an exterior ghost vertex, and false
otherwise.
See also is_ghost_vertex
and is_interior_ghost_vertex
.
Extended help
An exterior ghost vertex is a ghost vertex corresponding to a curve or section appearing on the exterior boundary.
DelaunayTriangulation.is_interior_ghost_vertex
— Functionis_interior_ghost_vertex(tri::Triangulation, i) -> Bool
Returns true
if the ghost vertex i
in tri
is an interior ghost vertex, and false
otherwise.
See also is_ghost_vertex
and is_exterior_ghost_vertex
.
Extended help
An interior ghost vertex is a ghost vertex corresponding to a curve or section appearing on the interior boundary.
DelaunayTriangulation.get_curve_index
— Functionget_curve_index(dict, ghost_vertex) -> Int
get_curve_index(ghost_vertex) -> Int
Given a Dict
from construct_ghost_vertex_map
and a ghost_vertex
, returns the index of the curve corresponding to that ghost vertex. The second method maps ghost_vertex
to 1
if it is an Integer
or a Vector
, and ghost_vertex[1]
if it is a Tuple
.
Examples
julia> using DelaunayTriangulation
julia> DelaunayTriangulation.get_curve_index(-1)
1
julia> DelaunayTriangulation.get_curve_index((5, 3))
5
julia> gv_map = DelaunayTriangulation.construct_ghost_vertex_map([[[1, 5, 17, 18, 1]], [[23, 29, 31, 33], [33, 107, 101], [101, 99, 85, 23]]])
Dict{Int64, Tuple{Int64, Int64}} with 4 entries:
-1 => (1, 1)
-3 => (2, 2)
-2 => (2, 1)
-4 => (2, 3)
julia> DelaunayTriangulation.get_curve_index(gv_map, -1)
1
julia> DelaunayTriangulation.get_curve_index(gv_map, -2)
2
julia> DelaunayTriangulation.get_curve_index(gv_map, -3)
2
julia> DelaunayTriangulation.get_curve_index(gv_map, -4)
2
DelaunayTriangulation.get_section_index
— Functionget_section_index(dict, ghost_vertex) -> Int
get_section_index(ghost_vertex) -> Int
Given a Dict
from construct_ghost_vertex_map
and a ghost_vertex
, returns the index of the section corresponding to that ghost vertex. The second method maps ghost_vertex
to itself if it is an Integer
, 1
if it is a Vector
, and ghost_vertex[2]
if it is a Tuple
.
Examples
julia> using DelaunayTriangulation
julia> DelaunayTriangulation.get_section_index((2, 3)) # 3rd section of the 2nd curve
3
julia> DelaunayTriangulation.get_section_index(4)
4
julia> DelaunayTriangulation.get_section_index([1, 2, 3, 4, 5, 1])
1
julia> gv_map = DelaunayTriangulation.construct_ghost_vertex_map([[[1, 5, 17, 18, 1]], [[23, 29, 31, 33], [33, 107, 101], [101, 99, 85, 23]]])
Dict{Int64, Tuple{Int64, Int64}} with 4 entries:
-1 => (1, 1)
-3 => (2, 2)
-2 => (2, 1)
-4 => (2, 3)
julia> DelaunayTriangulation.get_section_index(gv_map, -1)
1
julia> DelaunayTriangulation.get_section_index(gv_map, -2)
1
julia> DelaunayTriangulation.get_section_index(gv_map, -3)
2
julia> DelaunayTriangulation.get_section_index(gv_map, -4)
3
DelaunayTriangulation.map_ghost_vertex
— Functionmap_ghost_vertex(tri::Triangulation, ℓ) -> Vertex
Given a ghost vertex ℓ
in tri
, returns the corresponding section in the boundary_nodes
of tri
. See also get_ghost_vertex_map
.
DelaunayTriangulation.get_ghost_vertex_range
— Functionget_ghost_vertex_range(tri::Triangulation, ℓ) -> UnitRange
Given a ghost vertex ℓ
of tri
, returns the range of all ghost vertices corresponding to the same curve or section as ℓ
does.
DelaunayTriangulation.all_ghost_vertices
— Functionall_ghost_vertices(tri::Triangulation) -> KeySet
Returns the set of all ghost vertices in tri
.
DelaunayTriangulation.num_points
— Functionnum_points(points) -> Integer
Returns the number of points in points
.
Examples
julia> using DelaunayTriangulation
julia> points = [(1.0, 1.0), (2.3, 1.5), (0.0, -5.0)];
julia> DelaunayTriangulation.num_points(points)
3
julia> points = [1.0 5.5 10.0 -5.0; 5.0 2.0 0.0 0.0];
julia> DelaunayTriangulation.num_points(points)
4
DelaunayTriangulation.get_representative_point_coordinates
— Functionget_representative_point_coordinates(tri::Triangulation, curve_index) -> NTuple{2, Number}
Returns the coordinates of the representative point of the curve_index
th curve in tri
.
DelaunayTriangulation.compute_representative_points!
— Functioncompute_representative_points!(tri::Triangulation; use_convex_hull=!has_boundary_nodes(tri), precision=one(number_type(tri)))
Computes a new set of representative points for tri
.
Arguments
tri::Triangulation
: TheTriangulation
for which to compute the representative points.
Keyword Arguments
use_convex_hull=!has_boundary_nodes(tri)
: Iftrue
, then the representative points are computed using the convex hull of the triangulation. Otherwise, the representative points are computed using the boundary nodes of the triangulation.precision=one(number_type(tri))
: The precision to use when computing the representative points viapole_of_inaccessibility
.
Output
There are no outputs as tri
is updated in-place, but for each curve the representative point is computed using pole_of_inaccessibility
.
While get_exterior_curve_indices(tri)
does store the curves corresponding to exterior curves, this function still treats the first curve as the most important exterior curve, computing the representative point so that it is in no holes. In particular, other exterior curves might have representative points that are in a hole of one of their interior holes. This isn't much of a problem, indeed it wouldn't be a significant problem even if we had the representative point in a hole of the first curve, but it is something to be aware of.
DelaunayTriangulation.add_weight!
— Functionadd_weight!(weights, w)
Pushes the weight w
into weights
. The default definition for this is push!(weights, w)
.
add_weight!(tri::Triangulation, w)
Pushes the weight w
into the weights of tri
.
DelaunayTriangulation.get_weight
— Functionget_weight(weights, i) -> Number
Gets the i
th weight from weights
. The default definition for this is weights[i]
, but this can be extended - e.g., ZeroWeight
uses get_weight(weights, i) = 0.0
.
If i
is not an integer, then the default definition is is_point3(i) ? i[3] : zero(number_type(weights))
.
get_weight(tri::Triangulation, i) -> Number
Gets the i
th weight from tri
.
DelaunayTriangulation.ZeroWeight
— TypeZeroWeight{T}()
Struct used for indicating that a triangulation has zero weights. The weights have type T
, which is Float64
by default.
DelaunayTriangulation.get_area
— Functionget_area(r::BoundingBox) -> Float64
Returns the area of r
, i.e. hspan(r) * vspan(r)
.
get_area(stats::TriangulationStatistics)
Returns the area field from the TriangulationStatistics
stats
.
get_area(stats::TriangulationStatistics, T)
Returns the area field from the individual triangle statistics for the triangle T
in the TriangulationStatistics
stats
.
get_area(vor::VoronoiTessellation, i) -> Number
Gets the area of the i
th Voronoi polygon.
get_area(tri::Triangulation) -> Number
Returns the area of tri
.
DelaunayTriangulation.iterated_neighbourhood
— Functioniterated_neighbourhood(tri::Triangulation, i, d) -> Set{Vertex}
Returns the set of vertices of tri
in the iterated neighbourhood of the vertex i
of depth d
.
Extended help
The $d$-times iterated neighbourhood is defined by
\[N_i^d = \bigcup_{j \in N_i^{d-1}} N_j \setminus \{i\},\]
where $N_i^1 = N_i$ is the set of neighbours of $i$.
DelaunayTriangulation.dist
— Functiondist(p, q) -> Number
Assuming p
and q
are two-dimensional, computes the Euclidean distance between p
and q
.
dist(tri::Triangulation, p) -> Number
Given a point p
, returns the distance from p
to the triangulation, using the conventions from distance_to_polygon
:
δ > 0
: If the returned distance is positive, thenp
is inside the triangulation.δ < 0
: If the returned distance is negative, thenp
is outside the triangulation.δ = 0
: If the returned distance is zero, thenp
is on the boundary of the triangulation.
Where we say distance, we are referring to the distance from p
to the boundary of the triangulation.
DelaunayTriangulation.get_insertion_order
— Functionget_insertion_order(points, randomise, skip_points, ::Type{I}, rng) where {I} -> Vector{I}
get_insertion_order(tri::Triangulation, randomise, skip_points, rng) -> Vector{I}
Gets the insertion order for points into a triangulation.
Arguments
points
: The points to insert.randomise
: Iftrue
, then the insertion order is randomised. Otherwise, the insertion order is the same as the order of the points.skip_points
: The points to skip.I::Type{I}
: The type of the vertices.rng::Random.AbstractRNG
: The random number generator to use.
Output
order
: The order to insert the points in.
This order
might be mutated (by circshift!
) in get_initial_triangle
.
DelaunayTriangulation.get_vertices
— Methodget_vertices(tri::Triangulation) -> Set{Vertex}
Returns the set of all vertices in tri
. Note that, if has_ghost_triangles(tri)
, then some of these vertices will be ghost vertices.
See also each_vertex
, each_solid_vertex
, and each_ghost_vertex
.
DelaunayTriangulation.get_edges
— Methodget_edges(tri::Triangulation) -> Set{NTuple{2,Vertex}}
Returns the set of all edges in tri
. Orientation is ignored, so that only one of (i, j)
and (j, i)
will appear in the result. Note that, if has_ghost_triangles(tri)
, then some of these edges will be ghost edges.
See also each_edge
, each_solid_edge
, and each_ghost_edge
.
DelaunayTriangulation.num_neighbours
— Functionnum_neighbours(G::Graph, u) -> Integer
Returns the number of neighbours of u
in G
.
num_neighbours(tri::Triangulation, u) -> Integer
Returns the number of neighbours of u
in tri
. Note that, if has_ghost_triangles(tri)
, then some of the neighbours counted might be ghost vertices if u
is a boundary vertex.
DelaunayTriangulation.validate_triangulation
— Functionvalidate_triangulation(tri::Triangulation; print_result=true, predicates::AbstractPredicateKernel=ExactKernel()) -> Bool
Tests if tri
is a valid Triangulation
. Returns true
if so, and false
otherwise. If print_result=true
and tri
is not a valid triangulation, all the issues with tri
will be printed.
Use the predicates
keyword argument to control the method used for computing predicates. Can be one of FastKernel
, ExactKernel
, and AdaptiveKernel
. See the documentation for a further discussion of these methods.