Triangulation Operations
DelaunayTriangulation.add_boundary_information!
— Functionadd_boundary_information!(tri::Triangulation)
Updates tri
so that the ghost triangle information defined by the boundary nodes in tri
is added to the triangulation.
DelaunayTriangulation.add_ghost_triangles!
— Functionadd_ghost_triangles!(tri::Triangulation)
Adds all the ghost triangles to tri
.
DelaunayTriangulation.add_point!
— Functionadd_point!(c::RepresentativeCoordinates, p)
Treating c
as an arithmetic average, updates the coordinates of c
to include p
.
add_point!(tri::Triangulation, new_point; kwargs...) -> Triangle
add_point!(tri::Triangulation, x, y; kwargs...) -> Triangle
Arguments
tri::Triangulation
: TheTriangulation
.new_point
: The point to be added to the triangulation. The second method uses(x, y)
to represent the new point instead. Ifnew_point
is an integer, then the point added isget_point(tri, new_point)
.
Keyword Arguments
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.point_indices=each_solid_vertex(tri)
: The indices of the points to be used in thefind_triangle
algorithm for selecting the initial point.m=default_num_samples(length(point_indices))
: The number of samples (without replacement) to be used in thefind_triangle
algorithm for selecting the initial point.try_points=()
: Additional points to try for selecting the initial point, in addition to them
sampled.rng::Random.AbstractRNG=Random.default_rng()
: The random number generator to be used infind_triangle
.initial_search_point=integer_type(tri)(select_initial_point(tri, new_point; point_indices, m, try_points, rng))
: The initial point to be used infind_triangle
.update_representative_point=false
: Whether to update the representative point of the triangulation after adding the new point.store_event_history=Val(false)
: Whether to store the event history of the triangulation from adding the new point.event_history=nothing
: The event history of the triangulation from adding the new point. Only updated ifstore_event_history
is true, in which case it needs to be anInsertionEventHistory
object.concavity_protection=false
: Whether to use concavity protection for findingV
below. Seeconcavity_protection_check
. This is only needed if your triangulation is not convex.V=find_triangle(tri, get_point(tri, new_point); m=nothing, point_indices=nothing, try_points=nothing, k=initial_search_point, concavity_protection, rng)
: The positively oriented triangle containing the point being added.
In cases where your triangulation is not convex and !concavity_protection
, this V
may not be correct, and you may encounter errors - errors either during add_point!
or separately when you try to use the triangulation. In such cases, you should set concavity_protection=true
to ensure that V
is correct.
peek=Val(false)
: Whether the point should actually be added into the triangulation, or just 'peeked' at so that the events that would occur from its addition can be added intoevent_history
.
Outputs
The triangulation is updated in-place, but we do return
V
: The triangle containing the point being added.
In cases where (x, y)
is outside of the triangulation, it will be added successfully but note that the convex_hull
field of tri
will no longer be accurate. You can use convex_hull!
to fix it.
add_point!(tri::Triangulation, x, y, w; kwargs...)
Adds the point (x, y)
into tri
with weight w
. This function requires that add_weight!
is defined on the weights stored in tri
. The kwargs
match those from add_point!(tri::Triangulation, ::Any)
.
DelaunayTriangulation.add_segment!
— Functionadd_segment!(tri::Triangulation, segment; predicates::AbstractPredicateKernel=AdaptiveKernel(), rng::Random.AbstractRNG=Random.default_rng())
add_segment!(tri::Triangulation, i, j; predicates::AbstractPredicateKernel=AdaptiveKernel(), rng::Random.AbstractRNG=Random.default_rng())
Adds segment = (i, j)
to tri
.
Arguments
tri::Triangulation
: TheTriangulation
.segment
: The segment to add. The second method uses(i, j)
to represent the segment instead.
Keyword Arguments
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.rng::AbstractRNG=Random.default_rng()
: The RNG object.
Outputs
There is no output, but tri
will be updated so that it now contains segment
.
DelaunayTriangulation.add_triangle!
— Functionadd_triangle!(T, V...)
add_triangle!(T, i, j, k)
Add the triangles V...
or V = (i, j, k)
to the collection of triangles T
.
Examples
julia> using DelaunayTriangulation
julia> T = Set(((1, 2, 3), (4, 5, 6)))
Set{Tuple{Int64, Int64, Int64}} with 2 elements:
(4, 5, 6)
(1, 2, 3)
julia> add_triangle!(T, (7, 8, 9));
julia> add_triangle!(T, (10, 11, 12), (13, 14, 15));
julia> add_triangle!(T, 16, 17, 18);
julia> T
Set{Tuple{Int64, Int64, Int64}} with 6 elements:
(7, 8, 9)
(10, 11, 12)
(4, 5, 6)
(13, 14, 15)
(16, 17, 18)
(1, 2, 3)
DelaunayTriangulation.clear_empty_features!
— Functionclear_empty_features!(tri::Triangulation)
Clears all empty features from the triangulation tri
.
DelaunayTriangulation.delete_ghost_triangles!
— Functiondelete_ghost_triangles!(tri::Triangulation)
Deletes all the ghost triangles from tri
.
Ghost vertices are still used in the keys
of the Adjacent2Vertex
of tri
, and are still present in the Graph
. If you want to delete the ghost vertex keys
from the Adjacent2Vertex
, you need to use delete_adjacent2vertex!
. For deleting the ghost vertices from the Graph
, you need delete_ghost_vertices_from_graph!
. Additionally, edges in Adjacent
can still map to ghost vertices. If you also want to delete those, you need to filter through the values
of the Adjacent
map that are ghost vertices, and use delete_adjacent!
.
DelaunayTriangulation.delete_holes!
— Functiondelete_holes!(tri::Triangulation)
Deletes all the exterior faces to the boundary nodes specified in the triangulation tri
.
Extended help
This function works in several stages:
- First,
find_all_points_to_delete
is used to identify all points in the exterior faces. - Once all the points to delete have been found, all the associated triangles are found using
find_all_triangles_to_delete
, taking care for any incorrectly identified triangles and points. - Once the correct set of triangles to delete has been found, they are deleted using
delete_all_exterior_triangles!
.
DelaunayTriangulation.get_surrounding_polygon
— Functionget_surrounding_polygon(cache::TriangulationCache) -> Vector{Vertex}
Returns the polygon surrounding the triangulation stored in cache
.
get_surrounding_polygon(vor::VoronoiTessellation, i) -> Vector{Vertex}
Gets the polygon surrounding the generator with index i
in vor
.
You shouldn't need to use this, see get_polygon
instead.
get_surrounding_polygon(tri::Triangulation, u; skip_ghost_vertices=false) -> Vector
Returns the counter-clockwise sequence of neighbours of u
in tri
.
Arguments
tri::Triangulation
:Triangulation
.u
: The vertex.
Keyword Arguments
skip_ghost_vertices=false
: Whether to skip ghost vertices in the returned polygon.
Outputs
S
: The surrounding polygon. This will not be circular, meaningS[begin] ≠ S[end]
. In caseu
is an exterior ghost vertex, the returned polygon is a clockwise list of vertices for the associated boundary curve. If you do not have ghost triangles and you try to get the surrounding polygon of a ghost vertex, then this function may return an invalid polygon.
DelaunayTriangulation.delete_point!
— Functiondelete_point!(c::RepresentativeCoordinates, p)
Treating c
as an arithmetic average, updates the coordinates of c
to exclude p
.
delete_point!(tri::Triangulation, vertex; kwargs...)
Deletes the vertex
of tri
, retriangulating the cavity formed by the surrounding polygon of vertex
using triangulate_convex
.
It is not possible to delete vertices that are on the boundary, are ghost vertices, or adjoin a segment of tri
. See also check_delete_point_args
.
This function will not actually delete the corresponding coordinates from get_points(tri)
, nor will it remove the associated weight from get_weights(tri)
.
Arguments
tri::Triangulation
:Triangulation
.vertex
: The vertex to delete.
Keyword Arguments
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.store_event_history=Val(false)
: Whether to store the event history of the triangulation from deleting the point.event_history=nothing
: The event history of the triangulation from deleting the point. Only updated ifstore_event_history
is true, in which case it needs to be anInsertionEventHistory
object.rng::Random.AbstractRNG=Random.default_rng()
: The random number generator to use for the triangulation.
DelaunayTriangulation.delete_triangle!
— Functiondelete_triangle!(V, T...)
delete_triangle!(V, i, j, k)
Delete the triangles T...
from the collection of triangles V
.
Examples
julia> using DelaunayTriangulation
julia> V = Set(((1, 2, 3), (4, 5, 6), (7, 8, 9), (10, 11, 12), (13, 14, 15)))
Set{Tuple{Int64, Int64, Int64}} with 5 elements:
(7, 8, 9)
(10, 11, 12)
(4, 5, 6)
(13, 14, 15)
(1, 2, 3)
julia> delete_triangle!(V, (6, 4, 5))
Set{Tuple{Int64, Int64, Int64}} with 4 elements:
(7, 8, 9)
(10, 11, 12)
(13, 14, 15)
(1, 2, 3)
julia> delete_triangle!(V, (10, 11, 12), (1, 2, 3))
Set{Tuple{Int64, Int64, Int64}} with 2 elements:
(7, 8, 9)
(13, 14, 15)
julia> delete_triangle!(V, 8, 9, 7)
Set{Tuple{Int64, Int64, Int64}} with 1 element:
(13, 14, 15)
DelaunayTriangulation.flip_edge!
— Functionflip_edge!(tri::Triangulation, i, j, store_event_history=Val(false), event_history=nothing)
flip_edge!(tri::Triangulation, i, j, k, ℓ, store_event_history=Val(false), event_history=nothing)
Flips the edge between vertices i
and j
in tri
.
Arguments
tri::Triangulation
: TheTriangulation
.i
: The first vertex of the edge to flip.j
: The second vertex of the edge to flip.k
: The vertexk = get_adjacent(tri, j, i)
. This is only used in the second method.ℓ
: The vertexℓ = get_adjacent(tri, i, j)
. This is only used in the second method.store_event_history=Val(false)
: Whether to store the event history of the flip.event_history=nothing
: The event history. Only updated ifstore_event_history
is true, in which case it needs to be anInsertionEventHistory
object. This storage is done usingstore_flip_edge_history!
.
Outputs
There is no output, as tri
is updated in-place.
If (i, j, k, ℓ)
, where ℓ = get_adjacent(tri, i, j)
and k = get_adjacent(tri, j, i)
, is not a convex quadrilateral, then this edge flip will make the triangulation non-planar.
DelaunayTriangulation.legalise_edge!
— Functionlegalise_edge!(tri::Triangulation, i, j, r, store_event_history=Val(false), event_history=nothing; predicates::AbstractPredicateKernel=AdaptiveKernel())
Legalises the edge (i, j)
and other neighbouring edges in tri
if they are illegal, assuming the vertex r
was just added into a triangle that contains (i, j)
. flip_edge!
is used.
See also is_legal
.
Arguments
tri::Triangulation
: TheTriangulation
.i
: The first vertex of the edge to legalise.j
: The second vertex of the edge to legalise.r
: The vertex that was just added into a triangle that contains(i, j)
.store_event_history=Val(false)
: Whether to store the event history of the flip.event_history=nothing
: The event history. Only updated ifstore_event_history
is true, in which case it needs to be anInsertionEventHistory
object.
Keyword Arguments
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
There is no output, as tri
is updated in-place.
Edge flipping can lead to event_history
having triangles both in event_history.added_triangles
and event_history.deleted_triangles
. To get around this, we only store in these fields the triangles necessary to allow undo_insertion!
to work, so that at a triangle that might have appeared in both will only appear in one.
DelaunayTriangulation.lock_convex_hull!
— Functionlock_convex_hull!(tri::Triangulation; rng=Random.default_rng(), predicates::AbstractPredicateKernel=AdaptiveKernel())
Locks the convex hull of the unconstrained triangulation tri
so that it is now treated as a constrained triangulation with boundary given by its convex hull.
The random number generator (used inside add_segment!
can be provided with the rng
keyword argument, and similarly for predicates
.
If an edge is encountered along the convex hull that contains a segment from tri.interior_segments
, then this edge will be deleted from tri.interior_segments
; this will be undone from unlock_convex_hull!
, possibly splitting the segments in case they were split before unlocking.
DelaunayTriangulation.unlock_convex_hull!
— Functionunlock_convex_hull!(tri::Triangulation; reconstruct=false)
Unlocks the convex hull of the constrained triangulation tri
so that it is now treated as an unconstrained triangulation, assuming that it was locked using lock_convex_hull!
. If reconstruct = true
, then the convex hull of tri
will be reconstructed from the boundary nodes of tri
. This is useful if, for example, you have split some of the boundary edges during mesh refinement.
DelaunayTriangulation.split_edge!
— Functionsplit_edge!(tree::BoundaryRTree, i, j, r)
Splits the diametral bounding box associated with (i, j)
into two new boxes associated with the diametral circles of (i, r)
and (j, r)
.
split_edge!(enricher::BoundaryEnricher, i, j, r, update_boundary_nodes = Val(true), update_segments = Val(true), is_interior = is_segment(enricher, i, j))
Updates the fields of enricher
after splitting an edge (i, j)
at the r
th vertex. The update_boundary_nodes
argument can be used to avoid inserting an additional boundary node when boundary_nodes
was already updated somewhere else (e.g., we need this for mesh refinement which already updates the boundary_nodes
which is aliased with the same field in the enricher). The same point goes for update_segments
which can be used to avoid inserting an additional segment when segments
was already updated somewhere else. The is_interior
argument can be used to specify whether the edge is an interior segment or a boundary edge.
See also split_boundary_edge!
and split_interior_segment!
.
split_edge!(tri::Triangulation, i, j, r, store_event_history=Val(false), event_history=nothing)
Splits the edge (i, j)
in tri
at the vertex r
. For the triangulation to be valid after this splitting, it is assumed that r
is collinear with, or at least very close to collinear with, the edge (i, j)
.
See also legalise_split_edge!
and complete_split_edge_and_legalise!
.
Arguments
tri::Triangulation
: TheTriangulation
.i
: The first vertex of the edge to split.j
: The second vertex of the edge to split.r
: The vertex to split the edge at.store_event_history=Val(false)
: Whether to store the event history of the flip.event_history=nothing
: The event history. Only updated ifstore_event_history
is true, in which case it needs to be anInsertionEventHistory
object.
Outputs
There is no output, as tri
is updated in-place.
The triangulation will only be updated as if (i, j)
has been split rather than also (j, i)
. You will need to call split_edge!
again with (j, i)
if you want to split that edge as well.
DelaunayTriangulation.complete_split_edge_and_legalise!
— Functioncomplete_split_edge_and_legalise!(tri::Triangulation, i, j, r, store_event_history=Val(false), event_history=nothing; predicates::AbstractPredicateKernel=AdaptiveKernel())
Given a triangulation tri
, an edge (i, j)
, and a point r
, splits both (i, j)
and (j, i)
at r
using split_edge!
and then subsequently legalises the new edges with legalise_split_edge!
.
Arguments
tri::Triangulation
: TheTriangulation
.i
: The first vertex of the edge to split.j
: The second vertex of the edge to split.r
: The vertex to split the edge at.store_event_history=Val(false)
: Whether to store the event history of the flip.event_history=nothing
: The event history. Only updated ifstore_event_history
is true, in which case it needs to be anInsertionEventHistory
object.
Keyword Arguments
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
There is no output, as tri
is updated in-place.
DelaunayTriangulation.split_triangle!
— Functionsplit_triangle!(tri::Triangulation, args::RefinementArguments, T) -> Certificate
Splits a bad triangle T
of tri
to improve its quality.
Arguments
tri::Triangulation
: TheTriangulation
to split a triangle of.args::RefinementArguments
: TheRefinementArguments
for the refinement.T
: The triangle to split.
Output
cert
: ACertificate
indicating whether the split was successful or not. In particular, returns one of:Cert.SuccessfulInsertion
: The triangle was split successfully.Cert.EncroachmentFailure
: The triangle was not split successfully as the newly inserted point encroached upon a segment.Cert.PrecisionFailure
: The triangle was not split successfully due to precision issues.
split_triangle!(tri::Triangulation, i, j, k, r)
Splits the triangle (i, j, k)
at the vertex r
, assumed to be inside the triangle.
See also legalise_split_triangle!
and complete_split_triangle_and_legalise!
.
Arguments
tri::Triangulation
: TheTriangulation
.i
: The first vertex of the triangle.j
: The second vertex of the triangle.k
: The third vertex of the triangle.r
: The vertex to split the triangle at.
Outputs
There is no output, but tri
will be updated so that it now contains the triangles (i, j, r)
, (j, k, r)
, and (k, i, r)
.
DelaunayTriangulation.complete_split_triangle_and_legalise!
— Functioncomplete_split_triangle_and_legalise!(tri::Triangulation, i, j, k, r)
Splits the triangle (i, j, k)
at the vertex r
, assumed to be inside the triangle, and legalises the newly added edges in tri
.
Arguments
tri::Triangulation
: TheTriangulation
.i
: The first vertex of the triangle.j
: The second vertex of the triangle.k
: The third vertex of the triangle.r
: The vertex to split the triangle at.
Outputs
There is no output, as tri
is updated in-place.