Algorithm Internals
Here we give some functions that are used for some of the algorithms in this package.
Unconstrained Triangulations
Here are some of the internal functions used for the Bowyer-Watson algorithm, i.e. for the computation of an unconstrained triangulation.
DelaunayTriangulation.postprocess_triangulate!
— Functionpostprocess_triangulate!(tri; delete_ghosts=false, delete_empty_features=true, recompute_representative_points=true)
Postprocesses the triangulation tri
after it has been constructed using triangulate
. This includes:
- Deleting ghost triangles using
delete_ghost_triangles!
ifdelete_ghosts
istrue
. - Clearing empty features using
clear_empty_features!
ifdelete_empty_features
istrue
. - Recomputing the representative points using
compute_representative_points!
ifrecompute_representative_points
istrue
.
DelaunayTriangulation.add_point_bowyer_watson!
— Methodadd_point_bowyer_watson!(tri::Triangulation, new_point, initial_search_point::I, rng::Random.AbstractRNG=Random.default_rng(), update_representative_point=true, store_event_history=Val(false), event_history=nothing, peek::P=Val(false), predicates::AbstractPredicateKernel=AdaptiveKernel()) -> Triangle
Adds new_point
into tri
.
Arguments
tri
: The triangulation.new_point
: The point to insert.initial_search_point::I
: The vertex to start the point location withfind_triangle
at. Seeget_initial_search_point
.rng::Random.AbstractRNG
: The random number generator to use.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.update_representative_point=true
: Iftrue
, then the representative point is updated. Seeupdate_centroid_after_addition!
.store_event_history=Val(false)
: Iftrue
, then the event history from the insertion is stored.event_history=nothing
: The event history to store the event history in. Should be anInsertionEventHistory
ifstore_event_history
istrue
, andfalse
otherwise.peek=Val(false)
: Whether to actually addnew_point
intotri
, or just record intoevent_history
all the changes that would occur from its insertion.
Output
V
: The triangle intri
containingnew_point
.
Extended help
This function works as follows:
- First, the triangle containing the new point,
V
, is found usingfind_triangle
. - Once the triangle is found, we call into
add_point_bowyer_watson_and_process_after_found_triangle
to properly insert the point. - Inside
add_point_bowyer_watson_and_process_after_found_triangle
, we first call intoadd_point_bowyer_watson_after_found_triangle
to add the point into the cavity. We then call intoadd_point_bowyer_watson_onto_segment
to make any changes necessary incase the triangulation is constrained andnew_point
lies on a segment, since the depth-first search of the triangles containingnew_point
in its circumcenter must be performed on each side of the segment thatnew_point
lies on.
The function add_point_bowyer_watson_dig_cavities!
is the main workhorse of this function from add_point_bowyer_watson_after_found_triangle
. See its docstring for the details.
DelaunayTriangulation.add_point_bowyer_watson_dig_cavities!
— Methodadd_point_bowyer_watson_dig_cavities!(tri::Triangulation, new_point::N, V, q, flag, update_representative_point=true, store_event_history=Val(false), event_history=nothing, peek::F=Val(false), predicates::AbstractPredicateKernel=AdaptiveKernel()) where {N,F}
Deletes all the triangles in tri
whose circumcircle contains new_point
. This leaves behind a polygonal cavity, whose boundary edges are then connected to new_point
, restoring the Delaunay property from new_point
's insertion.
Arguments
tri
: TheTriangulation
.new_point::N
: The point to insert.V
: The triangle intri
containingnew_point
.q
: The point to insert.flag
: The position ofq
relative toV
. Seepoint_position_relative_to_triangle
.update_representative_point=true
: Iftrue
, then the representative point is updated. Seeupdate_centroid_after_addition!
.store_event_history=Val(false)
: Iftrue
, then the event history from the insertion is stored.event_history=nothing
: The event history to store the event history in. Should be anInsertionEventHistory
ifstore_event_history
istrue
, andfalse
otherwise.peek=Val(false)
: Whether to actually addnew_point
intotri
, or just record intoevent_history
all the changes that would occur from its insertion.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
There are no changes, but tri
is updated in-place.
Extended help
This function works as follows:
- To dig the cavity, we call
dig_cavity!
on each edge ofV
, stepping towards the adjacent triangles to excavate the cavity recursively. - Once the cavity has been excavated, extra care is needed in case
is_on(flag)
, meaningnew_point
is on one of the edges ofV
. In particular, extra care is needed if:is_on(flag) && (is_boundary_triangle(tri, V) || is_ghost_triangle(V) && !is_boundary_node(tri, new_point)[1])
. The need for this check is in casenew_point
is on a boundary edge already exists, since we need to fix the associated ghost edges. For example, a boundary edge(i, k)
might have been split into the edges(i, j)
and(j, k)
, which requires that the ghost triangle(k, i, g)
be split into(j, i, g)
and(k, j, g)
, whereg
is the ghost vertex. This part of the function will fix this case. The need foris_ghost_triangle(V) && !is_boundary_node(tri, new_point)[1]
is in case the ghost edges were already correctly added. Nothing happens if the edge ofV
thatnew_point
is on is not the boundary edge.
DelaunayTriangulation.dig_cavity!
— Methoddig_cavity!(tri::Triangulation, r, i, j, ℓ, flag, V, store_event_history=Val(false), event_history=nothing, peek::F=Val(false), predicates::AbstractPredicateKernel=AdaptiveKernel()) where {F}
Excavates the cavity in tri
through the edge (i, j)
, stepping towards the adjacent triangles to excavate the cavity recursively, eliminating all triangles containing r
in their circumcircle.
Arguments
tri
: TheTriangulation
.r
: The new point being inserted.i
: The first vertex of the edge(i, j)
.j
: The second vertex of the edge(i, j)
.ℓ
: The vertex adjacent to(j, i)
, so that the triangle being stepped into is(j, i, ℓ)
.flag
: The position ofr
relative toV
. Seepoint_position_relative_to_triangle
.V
: The triangle intri
containingr
.store_event_history=Val(false)
: Iftrue
, then the event history from the insertion is stored.event_history=nothing
: The event history to store the event history in. Should be anInsertionEventHistory
ifstore_event_history
istrue
, andfalse
otherwise.peek=Val(false)
: Whether to actually addnew_point
intotri
, or just record intoevent_history
all the changes that would occur from its insertion.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
There are no changes, but tri
is updated in-place.
Extended help
This function works as follows:
- First, we check if
ℓ
is 0, which would imply that the triangle(j, i, ℓ)
doesn't exist as it has already been deleted. If this is the check, we return and stop digging. - If the edge
(i, j)
is not a segment,r
is inside of the circumcenter of(j, i, ℓ)
, andℓ
is not a ghost vertex, then we can step forward into the next triangles. In particular, we delete the triangle(j, i, ℓ)
fromtri
and then calldig_cavity!
again on two edges of(j, i, ℓ)
other than(j, i)
. - If we did not do step 2, then this means that we are on the edge of the polygonal cavity; this also covers the case that
ℓ
is a ghost vertex, i.e. we are at the boundary of the triangulation. There are two cases to consider here. Firstly, if(i, j)
is not a segment, we just need to add the triangle(r, i, j)
intotri
to connectr
to this edge of the polygonal cavity. If(i, j)
is a segment, then the situation is more complicated. In particular,r
being on an edge ofV
might imply that we are going to add a degenerate triangle(r, i, j)
intotri
, and so this needs to be avoided. So, we check ifis_on(flag) && contains_segment(tri, i, j)
and, if the edge thatr
is on is(i, j)
, we add the triangle(r, i, j)
. Otherwise, we do nothing.
DelaunayTriangulation.enter_cavity
— Functionenter_cavity(tri::Triangulation, r, i, j, ℓ, predicates::AbstractPredicateKernel=AdaptiveKernel()) -> Bool
Determines whether to enter the cavity in tri
through the edge (i, j)
when inserting r
into the triangulation.
Arguments
tri
: TheTriangulation
.r
: The new point being inserted.i
: The first vertex of the edge(i, j)
.j
: The second vertex of the edge(i, j)
.ℓ
: The vertex adjacent to(j, i)
, so that the triangle being stepped into is(j, i, ℓ)
.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
true
if the cavity should be entered, andfalse
otherwise. See alsodig_cavity!
andpoint_position_relative_to_circumcircle
.
DelaunayTriangulation.get_initial_search_point
— Methodget_initial_search_point(tri::Triangulation, num_points, new_point, insertion_order, num_sample_rule::F, rng, try_last_inserted_point) where {F} -> Vertex
For a given iteration of the Bowyer-Watson algorithm, finds the point to start the point location with find_triangle
at.
Arguments
tri
: The triangulation.num_points
: The number of points currently in the triangulation.new_point
: The point to insert.insertion_order
: The insertion order of the points. Seeget_insertion_order
.num_sample_rule::F
: The rule to use to determine the number of points to sample. Seedefault_num_samples
for the default.rng::Random.AbstractRNG
: The random number generator to use.try_last_inserted_point
: Iftrue
, then the last inserted point is also considered as the start point.
Output
initial_search_point
: The vertex to start the point location withfind_triangle
at.
DelaunayTriangulation.get_initial_triangle
— Functionget_initial_triangle(tri::Triangulation, insertion_order, itr=0) -> Triangle
Gets the initial triangle for the Bowyer-Watson algorithm.
Arguments
tri
: The triangulation.insertion_order
: The insertion order of the points. Seeget_insertion_order
.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.itr=0
: To avoid issues with degenerate triangles and infinite loops, this counts the number of timesinsertion_order
had to be shifted usingcircshift!
to find an initial non-degenerate triangle.
Output
initial_triangle
: The initial triangle.
DelaunayTriangulation.get_insertion_order
— Methodget_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.initialise_bowyer_watson!
— Functioninitialise_bowyer_watson!(tri::Triangulation, insertion_order, predicates::AbstractPredicateKernel=AdaptiveKernel()) -> Triangulation
Initialises the Bowyer-Watson algorithm.
Arguments
tri
: The triangulation.insertion_order
: The insertion order of the points. Seeget_insertion_order
.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
tri
is updated in place to contain the initial triangle from which the Bowyer-Watson algorithm starts.
DelaunayTriangulation.unconstrained_triangulation!
— Methodunconstrained_triangulation!(tri::Triangulation; kwargs...)
Computes the unconstrained Delaunay triangulation of the points in tri
.
Arguments
tri
: The triangulation.
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.randomise=true
: Iftrue
, then the insertion order is randomised. Otherwise, the insertion order is the same as the order of the points.skip_points=()
: The vertices to skip.num_sample_rule::M=default_num_samples
: The rule to use to determine the number of points to sample. Seedefault_num_samples
for the default.rng::Random.AbstractRNG=Random.default_rng()
: The random number generator to use.insertion_order=get_insertion_order(tri, randomise, skip_points, rng)
: The insertion order of the points. Seeget_insertion_order
.
Outputs
There is no output, but tri
is updated in-place.
Triangulation Rectangular Domains
Here are some of the internal functions used for triangulate_rectangle
.
DelaunayTriangulation.get_lattice_triangles
— Functionget_lattice_triangles(nx, ny, Ts, V)
Computes the triangles defining a lattice with nx
and ny
points in the x
- and y
-directions, respectively.
Arguments
nx
: The number ofx
points in the lattice.ny
: The number ofy
points in the lattice.Ts
: The type to use for representing a collection of triangles.V
: The type to use for representing an individual triangle.
Outputs
T
: The collection of triangles.sub2ind
: A map that takes cartesian indices(i, j)
into the associated linear index along the lattice. SeeLinearIndices
.
DelaunayTriangulation.get_lattice_points
— Functionget_lattice_points(a, b, c, d, nx, ny, sub2ind)
Returns the points on a lattice [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.sub2ind
: The map returned fromget_lattice_triangles
.
Outputs
points
: The points on the lattice, wherepoints[sub2ind[CartesianIndex(i, j)]]
is the point at thei
thx
point and thej
thy
point,
DelaunayTriangulation.get_lattice_boundary
— Functionget_lattice_boundary(nx, ny, sub2ind, single_boundary, IntegerType)
Returns the boundary nodes defining a lattice.
Arguments
nx
: The number ofx
points in the lattice.ny
: The number ofy
points in the lattice.sub2ind
: The map returned fromget_lattice_triangles
.single_boundary
: 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.IntegerType
: The type used for representing vertices.
Outputs
boundary_nodes
: The boundary nodes, returned according tosingle_boundary
as described above.
Triangulating Convex Polygons
Here are some functions used in triangulate_convex
for triangulating a convex polygon.
DelaunayTriangulation.add_point_convex_triangulation!
— Functionadd_point_convex_triangulation!(tri::Triangulation, u, v, w, S, predicates::AbstractPredicateKernel=AdaptiveKernel())
Adds the point u
into the triangulation tri
.
Arguments
tri::Triangulation
: TheTriangulation
.u
: The vertex to add.v
: The vertex next tou
.w
: The vertex previous tou
.S
: The set of vertices of the polygon.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 modified in-place.
Extended help
This function forms part of Chew's algorithm for triangulating a convex polygon. There are some important points to make.
- Firstly, checking that
x = get_adjacent(tri, w, v)
is needed to prevent the algorithm from exiting the polygon.
This is important in case this algorithm is used as part of delete_point!
. When you are just triangulating a convex polygon by itself, this checked is the same as checking edge_exists(tri, w, v)
.
- For this method to be efficient, the set
x ∈ S
must beO(1)
time. This is why we use aSet
type forS
. - The algorithm is recursive, recursively digging further through the polygon to find non-Delaunay edges to adjoins with
u
.
DelaunayTriangulation.delete_vertices_in_random_order!
— Methoddelete_vertices_in_random_order!(list::Triangulation, tri::ShuffledPolygonLinkedList, rng, predicates::AbstractPredicateKernel)
Deletes the vertices of the polygon represented by list
in random order, done by switching the pointers of the linked list. Only three vertices will survive. If these these three vertices are collinear, then the deletion is attempted again after reshuffling the vertices.
Arguments
list::ShuffledPolygonLinkedList
: The linked list representing the polygon to be deleted.tri::Triangulation
: TheTriangulation
.rng::Random.AbstractRNG
: The random number generator used to shuffle the vertices ofS
.predicates::AbstractPredicateKernel
: 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 list
is modified in-place.
DelaunayTriangulation.postprocess_triangulate_convex!
— Methodpostprocess_triangulate_convex!(tri::Triangulation, S; delete_ghosts, delete_empty_features)
Postprocesses the completed triangulation tri
of the convex polygon S
.
Arguments
tri::Triangulation
: TheTriangulation
.S
: The vertices of the convex polygon, as intriangulate_convex
.
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.
Output
There are no output, as tri
is modified in-place. The postprocessing that is done is:
- The convex hull of
tri
is set toS
. - The ghost triangles are deleted if
delete_ghosts=true
. - The empty features are deleted if
delete_empty_features=true
. - The representative points are set to the centroid of
S
.
DelaunayTriangulation.triangulate_convex!
— Methodtriangulate_convex!(tri::Triangulation, S; rng::Random.AbstractRNG=Random.default_rng(), predicates::AbstractPredicateKernel=AdaptiveKernel())
Triangulates the convex polygon S
in-place into tri
.
Arguments
tri::Triangulation
: The triangulation to be modified.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
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.
Outputs
There is no output, as tri
is updated in-place. This function does not do any post-processing, e.g. deleting any ghost triangles. This is done by triangulate_convex
or postprocess_triangulate_convex!
.
DelaunayTriangulation.triangulate_convex
— Methodtriangulate_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.
Constrained Triangulations
There are many functions to list for the computation of a constrained triangulation.
DelaunayTriangulation.add_segment!
— Methodadd_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_segment_to_list!
— Methodadd_segment_to_list!(tri::Triangulation, e)
Adds e
to get_interior_segments(tri)
and get_all_segments(tri)
if it, or reverse_edge(e)
, is not already in the sets.
Arguments
tri::Triangulation
: TheTriangulation
.e
: The edge to add.
Outputs
There is no output, but tri
will be updated so that e
is in get_interior_segments(tri)
and get_all_segments(tri)
.
DelaunayTriangulation.fix_edge_order_after_rotation!
— Methodfix_edge_order_after_rotation!(tri::Triangulation, segment, e)
Fixes the edge order in get_interior_segments(tri)
after segment
was rotated by optimise_edge_order
.
Arguments
tri::Triangulation
: TheTriangulation
.segment
: The segment that was arranged.e
: The arranged segment fromoptimise_edge_order
.
Outputs
There is no output, but tri
will be updated so that e
is in get_interior_segments(tri)
instead of segment
.
DelaunayTriangulation.optimise_edge_order
— Methodoptimise_edge_order(tri::Triangulation, segment) -> Edge
Optimises the orientation of segment
for inserting it into the triangulation.
Arguments
tri::Triangulation
: TheTriangulation
.segment
: The segment to arrange.
Outputs
e
: Ifsegment
is a boundary edge, thene = segment
, Otherwise,e = sort_edge_by_degree(tri, segment)
so thatinitial(e)
has the smaller degree of the two vertices.
DelaunayTriangulation.add_new_triangles!
— Methodadd_new_triangles!(tri_original::Triangulation, tris)
Adds the triangles from tris
to tri_original
.
DelaunayTriangulation.add_point_cavity_cdt!
— Functionadd_point_cavity_cdt!(tri::Triangulation, u, v, w, marked_vertices)
Adds a point to the cavity V
left behind when deleting triangles intersected in a triangulation by an edge, updating tri
to do so.
Arguments
tri::Triangulation
: TheTriangulation
to update.u
: The vertex to add.v
: The vertex along the polygon that is next tou
.w
: The vertex along the polygon that is previous tou
.marked_vertices
: Cache for marking vertices to re-triangulate during the triangulation. This gets mutated.
Outputs
There is no output, but tri
is updated in-place, as is marked_vertices
if necessary.
DelaunayTriangulation.connect_segments!
— Methodconnect_segments!(segments)
Connects the ordered vector of segments
so that the endpoints all connect, preserving order.
Example
julia> using DelaunayTriangulation
julia> segments = [(7, 12), (12, 17), (17, 22), (32, 37), (37, 42), (42, 47)];
julia> DelaunayTriangulation.connect_segments!(segments)
7-element Vector{Tuple{Int64, Int64}}:
(7, 12)
(12, 17)
(17, 22)
(22, 32)
(32, 37)
(37, 42)
(42, 47)
DelaunayTriangulation.constrained_triangulation!
— Methodconstrained_triangulation!(tri::Triangulation, segments, boundary_nodes, predicates::AbstractPredicateKernel, full_polygon_hierarchy; rng=Random.default_rng(), delete_holes=true) -> Triangulation
Creates a constrained triangulation from tri
by adding segments
and boundary_nodes
to it. This will be in-place, but a new triangulation is returned to accommodate the changed types.
Arguments
tri::Triangulation
: TheTriangulation
.segments
: The interior segments to add to the triangulation.boundary_nodes
: The boundary nodes to add to the triangulation.predicates::AbstractPredicateKernel
: Method to use for computing predicates. Can be one ofFastKernel
,ExactKernel
, andAdaptiveKernel
. See the documentation for a further discussion of these methods.full_polygon_hierarchy
: ThePolygonHierarchy
defining the boundary. This will get copied into the existing polygon hierarchy.
Keyword Arguments
rng=Random.default_rng()
: The random number generator to use.delete_holes=true
: Whether to delete holes in the triangulation. Seedelete_holes!
.
Outputs
new_tri
: The new triangulation, now containingsegments
in theinterior_segments
field andboundary_nodes
in theboundary_nodes
field, and with the updatedPolygonHierarchy
. See alsoremake_triangulation_with_constraints
andreplace_ghost_vertex_information
.
DelaunayTriangulation.delete_intersected_triangles!
— Methoddelete_intersected_triangles!(tri, triangles)
Deletes the triangles in triangles
from tri
.
DelaunayTriangulation.delete_polygon_vertices_in_random_order!
— Methoddelete_polygon_vertices_in_random_order!(list::ShuffledPolygonLinkedList, tri::Triangulation, u, v, rng::Random.AbstractRNG, predicates::AbstractPredicateKernel)
Deletes vertices from the polygon defined by list
in a random order.
Arguments
list::ShuffledPolygonLinkedList
: The linked list of polygon vertices.tri::Triangulation
: TheTriangulation
.u
,v
: The vertices of the segment(u, v)
that was inserted in order to define the polygonV = list.S
.rng::Random.AbstractRNG
: The random number generator to use.predicates::AbstractPredicateKernel
: 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, but list
is updated in-place.
DelaunayTriangulation.extend_segments!
— Methodextend_segments!(segments, segment)
Given an ordered vector of segments
, ensures that they also represent the replacement of segment
. In particular, suppose segments
represents the sequence of edges
---(i₁, i₂)---(i₂, i₃)---(⋯, ⋯)---(iₙ₋₁, iₙ)---
and segment
is (i₀, iₙ₊₁)
. Then the extended sequence becomes
---(i₀, i₁)---(i₁, i₂)---(i₂, i₃)---(⋯, ⋯)---(iₙ₋₁, iₙ)---(iₙ, iₙ₊₁)---
Example
julia> using DelaunayTriangulation
julia> segments = [(2, 7), (7, 12), (12, 49)];
julia> segment = (1, 68);
julia> DelaunayTriangulation.extend_segments!(segments, segment)
5-element Vector{Tuple{Int64, Int64}}:
(1, 2)
(2, 7)
(7, 12)
(12, 49)
(49, 68)
DelaunayTriangulation.fix_segments!
— Methodfix_segments!(segments, bad_indices)
Fixes the overlapping segments in segments
, referred to via bad_indices
, by connecting consecutive edges where needed.
Arguments
segments
: The segments to fix.bad_indices
: The indices of the segments to fix.
Outputs
There are no outputs as segments
is updated in-place.
Example
julia> using DelaunayTriangulation
julia> segments = [(2, 15), (2, 28), (2, 41)]; # the edges all start with 2, so they are not actual segments in the triangulation, and so must be fixed
julia> bad_indices = [1, 2, 3];
julia> DelaunayTriangulation.fix_segments!(segments, bad_indices)
3-element Vector{Tuple{Int64, Int64}}:
(2, 15)
(15, 28)
(28, 41)
julia> segments = [(2, 7), (2, 12), (12, 17), (2, 22), (2, 27), (2, 32), (32, 37), (2, 42), (42, 47)];
julia> bad_indices = [2, 4, 5, 6, 8]
5-element Vector{Int64}:
2
4
5
6
8
julia> DelaunayTriangulation.fix_segments!(segments, bad_indices)
9-element Vector{Tuple{Int64, Int64}}:
(2, 7)
(7, 12)
(12, 17)
(17, 22)
(22, 27)
(27, 32)
(32, 37)
(37, 42)
(42, 47)
DelaunayTriangulation.is_vertex_closer_than_neighbours
— Methodis_vertex_closer_than_neighbours([predicates::AbstractPredicateKernel=AdaptiveKernel(),] tri::Triangulation, u, v, jᵢ, jᵢ₋₁, jᵢ₊₁) -> Bool
is_vertex_closer_than_neighbours([predicates::AbstractPredicateKernel=AdaptiveKernel(),] tri::Triangulation, list::ShuffledPolygonLinkedList, u, v, j) -> Bool
Tests if the vertex jᵢ
is closer to the line (u, v)
than its neighbours jᵢ₋₁
and jᵢ₊₁
, assuming all these vertices are to the left of the line.
See also point_closest_to_line
.
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.tri::Triangulation
: TheTriangulation
.u
,v
: The vertices of the line.jᵢ
,jᵢ₋₁
,jᵢ₊₁
: The vertices to compare.
The second method extracts these latter two vertices using the doubly-linked list
of vertices.
Outputs
flag
: Whetherjᵢ
is closer to the line thanjᵢ₋₁
andjᵢ₊₁
.
DelaunayTriangulation.locate_intersecting_triangles
— Functionlocate_intersecting_triangles(tri::Triangulation, e, rotate=Val(true), rng::Random.AbstractRNG=Random.default_rng(), predicates::AbstractPredicateKernel=AdaptiveKernel()) -> (Vector, Vector, Vector, Vector)
Find all the triangles intersected by an edge e
.
See also find_triangle
.
Arguments
tri::Triangulation
: TheTriangulation
.e
: The edge going through the triangulation.rotate=Val(true)
: Whether to rotate the edge so that the minimum degree vertex ofe
is first.rng::Random.AbstractRNG=Random.default_rng()
: The random number generator to use.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
intersecting_triangles
: The intersected triangles.collinear_segments
: Segments that are collinear withe
.left_vertices
: The vertices of the intersected triangles that are left ofe
.right_vertices
: The vertices of the intersected triangles that are right ofe
.
DelaunayTriangulation.merge_segments
— Methodmerge_segments(tri::Triangulation, ghost_vertex_map) -> Edges
Creates a set of edges that merge all the boundary nodes in tri
as well as its interior segments into a single collection.
Arguments
tri::Triangulation
: TheTriangulation
.ghost_vertex_map
: The ghost vertex map to use.
Outputs
all_segments
: The set of edges that merge all the boundary nodes intri
as well as its interior segments into a single collection, with type equal to that ofget_interior_segments(tri)
's.
DelaunayTriangulation.prepare_vertex_linked_list
— Methodprepare_vertex_linked_list(V) -> ShuffledPolygonLinkedList
Given a list of polygon vertices V
, returns the doubly-linked list of polygon vertices.
Arguments
V
: The list of polygon vertices.
Outputs
list::ShuffledPolygonLinkedList
: AShuffledPolygonLinkedList
. Inlist
,prev[begin]
,prev[end]
,next[begin]
, andnext[end]
are all0
as areshuffled_indices[begin]
andshuffled_indices[end]
. Moreover,shuffled_indices
will not have been shuffled yet.
DelaunayTriangulation.process_collinear_segments!
— Methodprocess_intersecting_triangles!(tri::Triangulation, e, collinear_segments; predicates::AbstractPredicateKernel=AdaptiveKernel(), rng::Random.AbstractRNG=Random.default_rng()) -> Bool
Given segments in collinear_segments
that are collinear with an edge e
, updates tri
so that this edge e
is instead split so that it is instead represented by collinear_segments
. These new segments will be placed into the triangulation using add_segment!
.
The predicates::AbstractPredicateKernel
argument defines the method for computing predicates. Can be one of FastKernel
, ExactKernel
, and AdaptiveKernel
. See the documentation for a further discussion of these methods.
See also connect_segments!
, extend_segments!
, split_segment!
and split_boundary_edge_at_collinear_segments!
.
DelaunayTriangulation.remake_triangulation_with_constraints
— Methodremake_triangulation_with_constraints(tri::Triangulation, segments, boundary_nodes) -> (Dict, Dict, Triangulation)
Remakes the triangulation tri
so that it contains segments
and boundary_nodes
in its fields.
See also replace_ghost_vertex_information
.
Arguments
tri::Triangulation
: The triangulation to remake.segments
: The segments to add to the triangulation.boundary_nodes
: The boundary nodes to add to the triangulation.
Outputs
new_ghost_vertex_map
: The new ghost vertex map. This will not yet be added to the triangulation.new_ghost_vertex_ranges
: The new ghost vertex ranges. This will not yet be added to the triangulation.new_tri::Triangulation
: The new triangulation, now containingboundary_nodes
in theboundary_nodes
field andsegments
in theinterior_segments
field.
DelaunayTriangulation.replace_ghost_vertex_information
— Functionreplace_ghost_vertex_information(tri::Triangulation, ghost_vertex_map, ghost_vertex_ranges, polygon_hierarchy) -> Triangulation
Replaces the ghost vertex information in tri
with ghost_vertex_map
and ghost_vertex_ranges
, using the results from remake_triangulation_with_constraints
.
Arguments
tri::Triangulation
: The triangulation to remake.ghost_vertex_map
: The ghost vertex map to add to the triangulation.ghost_vertex_ranges
: The ghost vertex ranges to add to the triangulation.polygon_hierarchy
: The polygon hierarchy to add to the triangulation.
Outputs
new_tri::Triangulation
: The new triangulation, now containingghost_vertex_map
in theghost_vertex_map
field andghost_vertex_ranges
in theghost_vertex_ranges
field.
DelaunayTriangulation.retriangulate_fan!
— Methodretriangulate_fan!(tri::Triangulation, tri_fan::Triangulation, fan, fan_triangles; predicates::AbstractPredicateKernel=AdaptiveKernel(), rng::Random.AbstractRNG=Random.default_rng())
Given a sorted set of vertices fan
in a fan of triangles associated with fan_triangles
, retriangulates the fan, updating tri
to do so and using tri_fan
as a temporary triangulation. (This implements Lines 17–19 and Line 28 of the algorithms in this paper.)
The predicates
argument defines the method for computing predicates. Can be one of FastKernel
, ExactKernel
, and AdaptiveKernel
. See the documentation for a further discussion of these methods.
DelaunayTriangulation.select_random_vertex
— Methodselect_random_vertex(tri::Triangulation, list::ShuffledPolygonLinkedList, u, v, range, rng) -> Vertex
Selects a random vertex that is not closer to the line (u, v)
than both of its neighbours.
Arguments
tri::Triangulation
: TheTriangulation
.list::ShuffledPolygonLinkedList
: The linked list of polygon vertices.u
,v
: The vertices of the line.range
: The range of indices of the vertices to select from.rng::Random.AbstractRNG
: The random number generator to use.predicates::AbstractPredicateKernel
: Method to use for computing predicates. Can be one ofFastKernel
,ExactKernel
, andAdaptiveKernel
. See the documentation for a further discussion of these methods.
Outputs
j
: The selected vertex.
DelaunayTriangulation.setup_cavity_cdt
— Methodsetup_cavity_cdt(tri::Triangulation, V; rng::Random.AbstractRNG=Random.default_rng(), predicates::AbstractPredicateKernel=AdaptiveKernel()) -> ShuffledPolygonLinkedList
Prepares the linked list required for triangulating a cavity excavated by segment insertion in a constrained triangulation.
See also prepare_vertex_linked_list
and delete_polygon_vertices_in_random_order!
.
Arguments
tri::Triangulation
: TheTriangulation
.V
: The list of polygon vertices, given as a counter-clockwise list of vertices, defining the cavity.
Keyword Arguments
rng::Random.AbstractRNG=Random.default_rng()
: The random number generator to use.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
list::ShuffledPolygonLinkedList
: The linked list of polygon vertices representing the cavity.
DelaunayTriangulation.sort_fan!
— Methodsort_fan!(fan, fan_triangles, tri::Triangulation)
Given a set of triangles in a fan, fan_triangles
, associated with some triangulation tri
, places all the triangle vertices and sorts them counter-clockwise, placing the results into fan
.
DelaunayTriangulation.split_marked_vertices!
— Methodsplit_marked_vertices!(fan_triangles, tri::Triangulation, marked_vertices)
Given a set of marked_vertices
indicating a crossed triangle (like in Figure 9 of this paper), finds all triangles whose three vertices are all in marked_vertices
and places them into fan_triangles
.
DelaunayTriangulation.split_segment!
— Methodsplit_segment!(tri::Triangulation, segment, collinear_segments)
split_segment!(segments, segment, collinear_segments)
Splits segment
at the segments in collinear_segments
, which are assumed to be collinear with segment
.
Arguments
tri::Triangulation
: TheTriangulation
.segments
: The underlying set of segments. This isget_interior_segments(tri)
iftri
is aTriangulation
.segment
: The segment to split.collinear_segments
: The segments that are collinear withsegment
.
Outputs
There is no output, as segments
is updated in-place.
Example
julia> using DelaunayTriangulation
julia> segments = Set(((2, 3), (3, 5), (10, 12)))
Set{Tuple{Int64, Int64}} with 3 elements:
(2, 3)
(3, 5)
(10, 12)
julia> collinear_segments = [(2, 10), (11, 15), (2, 3)]
3-element Vector{Tuple{Int64, Int64}}:
(2, 10)
(11, 15)
(2, 3)
julia> segment = (3, 5)
(3, 5)
julia> DelaunayTriangulation.split_segment!(segments, segment, collinear_segments)
Set{Tuple{Int64, Int64}} with 4 elements:
(2, 10)
(2, 3)
(11, 15)
(10, 12)
DelaunayTriangulation.triangulate_cavity_cdt!
— Methodtriangulate_cavity_cdt!(tri::Triangulation, V, marked_vertices; rng::Random.AbstractRNG=Random.default_rng(), predicates::AbstractPredicateKernel=AdaptiveKernel())
Triangulates the cavity V
left behind when deleting triangles intersected in a triangulation by an edge, updating tri
to do so.
Arguments
tri::Triangulation
: TheTriangulation
to update. This should be an empty triangulation.V
: The list of polygon vertices, given as a counter-clockwise list of vertices, defining the cavity.tri_fan::Triangulation
: TheTriangulation
to use for the fan of triangles to be re-triangulated. This should be an empty triangulation.marked_vertices
: Cache for marking vertices to re-triangulate during the triangulation.fan_triangles
: A cache used for sorting and identifying triangles in a fan for retriangulation.
Keyword Arguments
rng::Random.AbstractRNG=Random.default_rng()
: The random number generator to use orsetup_cavity_cdt
.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, but tri
is updated in-place.
DelaunayTriangulation.each_boundary_edge
— Methodeach_boundary_edge(tri::Triangulation) -> KeySet
Returns an iterator over the boundary edges of tri
, in no specific order.
DelaunayTriangulation.get_boundary_edge_map
— Methodget_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.split_boundary_edge_map!
— Methodsplit_boundary_edge_map!(boundary_edge_map, boundary_nodes, pos)
After splitting an edge starting at pos
on the boundary, updates the boundary_edge_map
to reflect the new boundary edges. See split_boundary_edge!
.
DelaunayTriangulation.contains_boundary_edge
— Methodcontains_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.delete_boundary_node!
— Methoddelete_boundary_node!(tri::Triangulation, pos)
Deletes the boundary node at the specified position pos
in tri
.
Arguments
tri::Triangulation
: TheTriangulation
.pos
: The position to delete the node at, given as aTuple
so thatdelete_boundary_node!(tri, pos)
is the same asdeleteat!(get_boundary_nodes(tri, pos[1]), pos[2])
.
Outputs
There are no outputs, but the boundary nodes of tri
are updated in-place.
DelaunayTriangulation.get_all_boundary_nodes
— Methodget_all_boundary_nodes(tri::Triangulation) -> Set{Vertex}
Returns the set of all boundary vertices in tri
, in no specific order.
DelaunayTriangulation.get_boundary_nodes
— Methodget_boundary_nodes(tri, mnℓ...)
Given a triangulation tri
, returns the specified component of the boundary nodes. There are several forms for the methods:
get_boundary_nodes(tri, m)
: Iftri
has multiple curves, this returns them
th curve. Iftri
has multiple sections, this returns them
th section. Otherwise, this returns them
th boundary node.get_boundary_nodes(tri, m, n)
: Iftri
has multiple curves, this returns then
th section of them
th curve. Otherwise, iftri
has multiple sections, this returns then
th boundary node of them
th section.get_boundary_nodes(tri, (m, n))
: This is equivalent toget_boundary_nodes(tri, m, n)
.get_boundary_nodes(tri::A, ::A)
: This just returnsboundary_nodes
.
DelaunayTriangulation.get_left_boundary_node
— Methodget_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_right_boundary_node
— Methodget_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.insert_boundary_node!
— Methodinsert_boundary_node!(tri::Triangulation, pos, node)
Inserts a boundary node into tri
at the specified position.
Arguments
tri::Triangulation
: TheTriangulation
.pos
: The position to insert the node at, given as aTuple
so thatinsert_boundary_node!(tri, pos, node)
is the same asinsert!(get_boundary_nodes(tri, pos[1]), pos[2], node)
.node
: The node to insert.
Outputs
There are no outputs, but the boundary nodes of tri
are updated in-place.
DelaunayTriangulation.is_positively_oriented
— Methodis_positively_oriented(tri::Triangulation, curve_index) -> Bool
Tests if the curve_index
th curve in tri
is positively oriented, returning true
if so and false
otherwise.
DelaunayTriangulation.merge_boundary_edge!
— Methodmerge_boundary_edge!(tri::Triangulation, ij, node)
merge_boundary_edge!(tri::Triangulation, i, j, node)
Merges the edges (i, node)
and (node, j)
into a single edge (i, j)
, i.e. does the inverse operation to split_boundary_edge!
.
DelaunayTriangulation.num_curves
— Methodnum_curves(tri::Triangulation) -> Integer
Returns the number of curves in tri
.
DelaunayTriangulation.num_sections
— Methodnum_sections(tri::Triangulation) -> Integer
Assuming tri
only has one curve, returns the number of sections in tri
.
DelaunayTriangulation.split_boundary_edge!
— Methodsplit_boundary_edge!(tri::Triangulation, ij, node)
split_boundary_edge!(tri::Triangulation, i, j, node)
Splits the boundary edge edge
in tri
at the edge (i, j)
.
See also merge_boundary_edge!
.
DelaunayTriangulation.split_boundary_edge_at_collinear_segments!
— Methodsplit_boundary_edge_at_collinear_segments!(tri::Triangulation, collinear_segments)
Splits a boundary edge into pieces defined by collinear_segments
. In particular, if r = collinear_segments
and
u = initial(r[1])
v = terminal(r[end]),
then the boundary edge is (u, v)
and the edges are split so that all segments in collinear_segments
appear instead.
DelaunayTriangulation.is_disjoint
— Methodis_disjoint(tri::Triangulation) -> Bool
Returns true
if tri
has disjoint exterior boundary curves, and false
otherwise.
DelaunayTriangulation.is_exterior_curve
— Methodis_exterior_curve(tri::Triangulation, curve_index) -> Bool
Returns true
if the curve_index
th curve in tri
is an exterior curve, and false
otherwise.
DelaunayTriangulation.is_interior_curve
— Methodis_interior_curve(tri::Triangulation, curve_index) -> Bool
Returns true
if the curve_index
th curve in tri
is an interior curve, and false
otherwise.
DelaunayTriangulation.num_exterior_curves
— Methodnum_exterior_curves(tri::Triangulation) -> Integer
Returns the number of exterior curves in tri
.
DelaunayTriangulation.get_curve_index
— Methodget_curve_index(tri::Triangulation, ℓ) -> Integer
Returns the curve index corresponding to the ghost vertex ℓ
in tri
.
DelaunayTriangulation.get_section_index
— Methodget_section_index(tri::Triangulation, ℓ) -> Integer
Returns the section index corresponding to the ghost vertex ℓ
in tri
.
DelaunayTriangulation.is_exterior_ghost_vertex
— Methodis_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
— Methodis_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.map_ghost_vertex
— Methodmap_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.all_ghost_vertices
— Methodall_ghost_vertices(tri::Triangulation) -> KeySet
Returns the set of all ghost vertices in tri
.
DelaunayTriangulation.get_ghost_vertex_range
— Methodget_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.compute_representative_points!
— Methodcompute_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.empty_representative_points!
— Methodempty_representative_points!(tri::Triangulation)
Empties the Dict
of representative points of tri
.
DelaunayTriangulation.get_representative_point_coordinates
— Methodget_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.new_representative_point!
— Methodnew_representative_point!(tri::Triangulation, curve_index)
Creates a new representative point for the curve_index
th curve in tri
.
DelaunayTriangulation.reset_representative_points!
— Methodreset_representative_points!(tri::Triangulation)
Resets each representative point of tri
to the origin.
DelaunayTriangulation.update_centroid_after_addition!
— Methodupdate_centroid_after_addition!(tri::Triangulation, curve_index, p)
Updates the centroid of the curve_index
th curve in tri
after the addition of the point p
.
DelaunayTriangulation.update_centroid_after_deletion!
— Methodupdate_centroid_after_deletion!(tri::Triangulation, curve_index, p)
Updates the centroid of the curve_index
th curve in tri
after the deletion of the point p
.
DelaunayTriangulation.contains_segment
— Methodcontains_segment(tri::Triangulation, ij) -> Bool
contains_segment(tri::Triangulation, i, j) -> Bool
Returns true
if (i, j)
is a segment in tri
, and false
otherwise. Both (i, j)
and (j, i)
are checked.
DelaunayTriangulation.each_segment
— Methodeach_segment(tri::Triangulation) -> Edges
Returns an iterator over all segments in tri
. This includes both interior and boundary segments. If you only want interior segments, then see get_interior_segments
,
Weighted Triangulations
DelaunayTriangulation.add_weight!
— Methodadd_weight!(weights, w)
Pushes the weight w
into weights
. The default definition for this is push!(weights, w)
.
DelaunayTriangulation.add_weight!
— Methodadd_weight!(tri::Triangulation, w)
Pushes the weight w
into the weights of tri
.
DelaunayTriangulation.get_distance_to_witness_plane
— Method" getdistancetowitnessplane([kernel::AbstractPredicateKernel = AdaptiveKernel(), ] tri::Triangulation, i, V; cache = nothing) -> Number
Computes the distance between the lifted companion of the vertex i
and the witness plane to the triangle V
. If V
is a ghost triangle and i
is not on its solid edge, then the distance is -Inf
if it is below the ghost triangle's witness plane and Inf
if it is above. If V
is a ghost triangle and i
is on its solid edge, then the distance returned is the distance associated with the solid triangle adjoining V
.
In general, the distance is positive if the lifted vertex is above the witness plane, negative if it is below, and zero if it is on the plane.
The kernel
argument determines how this result is computed, and should be one of ExactKernel
, FastKernel
, and AdaptiveKernel
(the default). See the documentation for more information about these choices.
The cache
keyword argument is passed to [point_position_relative_to_circumcircle
]. Please see the documentation for that function for more information.
See also point_position_relative_to_witness_plane
and get_distance_to_plane
.
DelaunayTriangulation.get_lifted_point
— Methodget_lifted_point(p, w) -> NTuple{3, Number}
Returns the lifted companion of the point p
, in particular (x, y, x^2 + y^2 - w)
, where (x, y)
is p
.
DelaunayTriangulation.get_lifted_point
— Methodget_lifted_point(tri::Triangulation, i) -> NTuple{3, Number}
Returns the lifted companion of the i
th vertex of tri
, in particular (x, y, x^2 + y^2 - w)
, where w
is the i
th weight of tri
and (x, y)
is the i
th point of tri
.
DelaunayTriangulation.get_power_distance
— Methodget_power_distance(tri::Triangulation, i, j) -> Number
Returns the power distance between vertices i
and j
, defined by ||pᵢ - pⱼ||^2 - wᵢ - wⱼ
, where wᵢ
and wⱼ
are the respective weights.
DelaunayTriangulation.get_weight
— Methodget_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))
.
DelaunayTriangulation.get_weight
— Methodget_weight(tri::Triangulation, i) -> Number
Gets the i
th weight from tri
.
DelaunayTriangulation.get_weighted_nearest_neighbour
— Functionget_weighted_nearest_neighbour(tri::Triangulation, i, j = rand(each_solid_vertex(tri))) -> Vertex
Using a greedy search, finds the closest vertex in tri
to the vertex i
(which might not already be in tri
), measuring distance in lifted space (i.e., using the power distance - see get_power_distance
). The search starts from the vertex j
which should be in tri
.
DelaunayTriangulation.is_submerged
— Functionis_submerged([kernel::AbstractPredicateKernel = AdaptiveKernel(), ] tri::Triangulation, i; cache = nothing) -> Bool
is_submerged([kernel::AbstractPredicateKernel = AdaptiveKernel(), ] tri::Triangulation, i, V; cache = nothing) -> Bool
Returns true
if the vertex i
is submerged in tri
and false
otherwise. In the second method, V
is a triangle containing tri
.
The kernel
argument determines how this result is computed, and should be one of ExactKernel
, FastKernel
, and AdaptiveKernel
(the default). See the documentation for more information about these choices.
The cache
keyword argument is passed to [point_position_relative_to_circumcircle
]. Please see the documentation for that function for more information.
DelaunayTriangulation.is_weighted
— Methodis_weighted(weights) -> Bool
Returns true
if weights
represents a set of weights
that are not all zero, and false
otherwise. Note that even for vectors like zeros(n)
, this will return true
; by default, false
is returned only for weights = ZeroWeight()
.
Mesh Refinement
Here are some functions involved with mesh refinement.
DelaunayTriangulation._split_subsegment_curve_bounded!
— Method_split_subsegment_curve_bounded!(tri::Triangulation, args::RefinementArguments, e)
Splits a subsegment e
of tri
at a position determined by split_subcurve!
for curve-bounded domains. See split_subsegment!
. See also _split_subsegment_curve_bounded_standard!
and _split_subsegment_curve_bounded_small_angle!
, as well as the original functions _split_subsegment_curve_standard!
and _split_subcurve_complex!
, respectively, used during boundary enrichment.
DelaunayTriangulation._split_subsegment_piecewise_linear!
— Method_split_subsegment_piecewise_linear!(tri::Triangulation, args::RefinementArguments, e)
Splits a subsegment e
of tri
at a position determined by compute_split_position
for piecewise linear domains. See split_subsegment!
.
DelaunayTriangulation.assess_added_triangles!
— Methodassess_added_triangles!(args::RefinementArguments, tri::Triangulation)
Assesses the quality of all triangles in args.events.added_triangles
according to assess_triangle_quality
, and enqueues any bad quality triangles into args.queue
.
DelaunayTriangulation.assess_triangle_quality
— Methodassess_triangle_quality(tri::Triangulation, args::RefinementArguments, T) -> Float64, Bool
Assesses the quality of a triangle T
of tri
according to the RefinementArguments
.
Arguments
tri::Triangulation
: TheTriangulation
.args::RefinementArguments
: TheRefinementArguments
.T
: The triangle.
Output
ρ
: The radius-edge ratio of the triangle.flag
: Whether the triangle is bad quality.
A triangle is bad quality if it does not meet the area constraints, violates the custom constraint, or if it is skinny but neither seditious or nestled.
DelaunayTriangulation.balanced_power_of_two_quarternary_split
— Methodbalanced_power_of_two_quarternary_split(ℓ) -> Float
Returns the value of s ∈ [0, ℓ]
that gives the most balanced quarternary split of the segment pq
, so s
is a power-of-two and s ∈ [ℓ / 4, ℓ / 2]
.
DelaunayTriangulation.balanced_power_of_two_ternary_split
— Methodbalanced_power_of_two_ternary_split(ℓ) -> Float64
Returns the value of s ∈ [0, ℓ]
that gives the most balanced ternary split of the segment pq
, so s
is a power-of-two and s ∈ [ℓ / 3, 2ℓ / 3]
.
DelaunayTriangulation.check_for_invisible_steiner_point
— Methodcheck_for_invisible_steiner_point(tri::Triangulation, V, T, flag, c) -> Point, Triangle
Determines if the Steiner point c
's insertion will not affect the quality of T
, and if so instead changes c
to be T
's centroid.
Arguments
tri::Triangulation
: TheTriangulation
to split a triangle of.V
: The triangle that the Steiner point is in.T
: The triangle that the Steiner point is from.flag
: ACertificate
which isCert.On
if the Steiner point is on the boundary ofV
,Cert.Outside
if the Steiner point is outside ofV
, andCert.Inside
if the Steiner point is inside ofV
.c
: The Steiner point.
Output
c′
: The Steiner point to use instead ofc
, which isT
's centroid ifc
is not suitable.V′
: The triangle that the Steiner point is in, which isT
ifc
is not suitable.
DelaunayTriangulation.check_for_steiner_point_on_segment
— Methodcheck_for_steiner_point_on_segment(tri::Triangulation, V, V′, new_point, flag, predicates::AbstractPredicateKernel) -> Bool
Checks if the Steiner point with vertex new_point
is on a segment. If so, then its vertex is pushed into the offcenter-split list from args
, indicating that it should no longer be regarded as a free vertex (see is_free
).
Arguments
tri::Triangulation
: TheTriangulation
.V
: The triangle that the Steiner point was originally in prior tocheck_for_invisible_steiner_point
.V′
: The triangle that the Steiner point is in.new_point
: The vertex associated with the Steiner point.flag
: ACertificate
which isCert.On
if the Steiner point is on the boundary ofV
,Cert.Outside
if the Steiner point is outside ofV
, andCert.Inside
if the Steiner point is inside ofV
.predicates::AbstractPredicateKernel
: Method to use for computing predicates. Can be one ofFastKernel
,ExactKernel
, andAdaptiveKernel
. See the documentation for a further discussion of these methods.
Output
onflag
: Whether the Steiner point is on a segment or not.
DelaunayTriangulation.check_seditious_precision
— Methodcheck_seditious_precision(ℓrp, ℓrq) -> Bool
Checks if there are precision issues related to the seditiousness of a triangle, returning true
if so and false
otherwise.
DelaunayTriangulation.check_split_subsegment_precision
— Methodcheck_split_subsegment_precision(mx, my, p, q) -> Bool
Checks if there are precision issues related to the computed split position (mx, my)
of a segment (p, q)
, returning true
if so and false
otherwise.
DelaunayTriangulation.check_steiner_point_precision
— Methodcheck_steiner_point_precision(tri::Triangulation, T, c) -> Bool
Checks if the Steiner point c
of a triangle T
of tri
can be computed without precision issues, returning true
if there are precision issues and false
otherwise.
DelaunayTriangulation.compute_concentric_shell_quarternary_split_position
— Methodcompute_concentric_shell_quarternary_split_position(p, q) -> Float64
Returns the value of t ∈ [0, 1]
that gives the most balanced quarternary split of the segment pq
, so that one of the segments has a power-of-two length between 1/4 and 1/2 of the length of pq
.
DelaunayTriangulation.compute_concentric_shell_ternary_split_position
— Methodcompute_concentric_shell_ternary_split_position(p, q) -> Float64
Returns the value of t ∈ [0, 1]
that gives the most balanced ternary split of the segment pq
, so that one of the segments has a power-of-two length and both segments have lengths between 1/3 and 2/3 of the length of pq
.
DelaunayTriangulation.compute_split_position
— Methodcompute_split_position(tri::Triangulation, args::RefinementArguments, e) -> NTuple{2, Float}
Computes the position to split a segment e
of tri
at in split_subsegment!
.
Arguments
tri::Triangulation
: TheTriangulation
to split a segment of.args::RefinementArguments
: TheRefinementArguments
for the refinement.e
: The segment to split.
Output
mx, my
: The position to split the segment at.
This point is computed according to a set of rules:
- If
e
is not a subsegment, meaning it is an input segment, then its midpoint is returned. - If
e
is a subsegment and the segment adjoins two other distinct segments (one for each vertex) at an acute angle, as determined bysegment_vertices_adjoin_other_segments_at_acute_angle
, then the point is returned so thate
can be split such that one of the new subsegments has a power-of-two length between 1/4 and 1/2 of the length ofe
, computed usingcompute_concentric_shell_quarternary_split_position
. - If
e
is a subsegment and the segment adjoins one other segment at an acute angle, as determined bysegment_vertices_adjoin_other_segments_at_acute_angle
, then the point is returned so thate
can be split such that one of the new subsegments has a power-of-two length between 1/3 and 2/3 of the length ofe
, computed usingcompute_concentric_shell_ternary_split_position
. - Otherwise, the midpoint is returned.
DelaunayTriangulation.delete_free_vertices_around_subsegment!
— Methoddelete_free_vertices_around_subsegment!(tri::Triangulation, args::RefinementArguments, e)
Deletes all free vertices (i.e., Steiner points) contained in the diametral circle of e
.
Arguments
tri::Triangulation
: TheTriangulation
.args::RefinementArguments
: TheRefinementArguments
.e
: The segment.
Output
The free vertices are deleted from tri
in-place.
DelaunayTriangulation.encroaches_upon
— Methodencroaches_upon(p, q, r, args::RefinementArguments) -> Bool
Determines if a point r
encroaches upon a segment pq
according to the RefinementArguments
.
See also point_position_relative_to_diametral_circle
and point_position_relative_to_diametral_lens
.
DelaunayTriangulation.enqueue_all_bad_triangles!
— Methodenqueue_all_bad_triangles!(args::RefinementArguments, tri::Triangulation)
Enqueues all bad triangles in the triangulation into args.queue
.
DelaunayTriangulation.enqueue_all_encroached_segments!
— Methodenqueue_all_encroached_segments!(args::RefinementArguments, tri::Triangulation)
Enqueues all encroached segments in the triangulation into args.queue
.
DelaunayTriangulation.enqueue_newly_encroached_segments!
— Methodenqueue_newly_encroached_segments!(args::RefinementArguments, tri::Triangulation) -> Bool
Enqueues all segments that are newly encroached upon after a point insertion into the triangulation into args.queue
.
Arguments
args::RefinementArguments
: TheRefinementArguments
for the refinement.tri::Triangulation
: TheTriangulation
to enqueue newly encroached segments of.
Output
any_encroached
: Whether any segments were newly encroached upon.
DelaunayTriangulation.finalise!
— Methodfinalise!(tri::Triangulation, args::RefinementArguments)
Finalises the triangulation after refinement, e.g. by deleting ghost triangles and unlocking the convex hull if needed.
DelaunayTriangulation.get_init_for_steiner_point
— Methodget_init_for_steiner_point(tri::Triangulation, T) -> Vertex
Gets the initial vertex to start the search for the Steiner point of a triangle T
of tri
in get_steiner_point
. The initial vertex is chosen so that it is opposite the longest edge.
DelaunayTriangulation.get_steiner_point
— Methodget_steiner_point(tri::Triangulation, args::RefinementArguments, T) -> Certificate, Point
Computes the Steiner point for a triangle T
of tri
to improve its quality in split_triangle!
.
Arguments
tri::Triangulation
: TheTriangulation
to split a triangle of.args::RefinementArguments
: TheRefinementArguments
for the refinement.T
: The triangle to split.
Output
precision_flag
: ACertificate
which isCert.PrecisionFailure
if the Steiner point could not be computed due to precision issues, andCert.None
otherwise.c
: The Steiner point. Ifis_precision_failure(precision_flag)
, then this is just an arbitrary point ofT
to ensure type stability.
DelaunayTriangulation.is_encroached
— Methodis_encroached(tri::Triangulation, args::RefinementArguments, edge) -> Bool
Determines if a segment edge
of tri
is encroached upon.
See also encroaches_upon
.
DelaunayTriangulation.is_triangle_nestled
— Methodis_triangle_nestled(tri::Triangulation, T, idx) -> Bool
Determines if a triangle T
of tri
is nestled in the corner of a small input angle.
See also is_triangle_seditious
.
Arguments
tri::Triangulation
: TheTriangulation
.T
: The triangle.idx
: The index of the smallest edge of the triangle, so that1
meansuv
is the smallest edge,2
meansvw
is the smallest edge, and3
meanswu
is the smallest edge.
Output
flag
: Whether the triangle is nestled in the corner of a small input angle.
A triangle is nestled in the corner of a small input angle if it is nestled in the corner of a small input angle and the shortest edge is seditious. The size of the angle is not checked by this function, and is instead determined by assess_triangle_quality
.
DelaunayTriangulation.is_triangle_seditious
— Methodis_triangle_seditious(tri::Triangulation, args, u, v, w, smallest_idx) -> Bool
Determines if a triangle uvw
of tri
is seditious according to the RefinementArguments
.
See also is_triangle_nestled
.
Arguments
tri::Triangulation
: TheTriangulation
.args
: TheRefinementArguments
.u, v, w
: The vertices of the triangle.smallest_idx
: The index of the smallest edge of the triangle, so that1
meansuv
is the smallest edge,2
meansvw
is the smallest edge, and3
meanswu
is the smallest edge.
Output
flag
: Whether the triangle is seditious.
A triangle is seditious if it is nestled in the corner of a small input angle, or if it is nestled in the corner of a small input angle and the shortest edge is seditious. Here, 'small' is defined by args.constraints.seditious_angle
.
DelaunayTriangulation.locate_steiner_point
— Methodlocate_steiner_point(tri::Triangulation, args::RefinementArguments, T, c) -> Triangle, Cert
Locates the Steiner point c
of a triangle T
of tri
in get_steiner_point
. The Steiner point is located by walking from the initial vertex init
to c
using find_triangle
.
Arguments
tri::Triangulation
: TheTriangulation
to split a triangle of.args::RefinementArguments
: TheRefinementArguments
for the refinement.T
: The triangle that the Steiner point is from.c
: The Steiner point.
Output
V
: The triangle that the Steiner point is in.flag
: ACertificate
which isCert.On
if the Steiner point is on the boundary ofV
,Cert.Outside
if the Steiner point is outside ofV
, andCert.Inside
if the Steiner point is inside ofV
.
DelaunayTriangulation.refine!
— Methodrefine!(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.refine_itr!
— Methodrefine_itr!(tri::Triangulation, args::RefinementArguments)
Performs a single iteration of the refinement algorithm.
Arguments
tri::Triangulation
: TheTriangulation
to refine.args::RefinementArguments
: TheRefinementArguments
for the refinement.
Output
The triangulation is refined in-place.
DelaunayTriangulation.segment_vertices_adjoin_other_segments_at_acute_angle
— Functionsegment_vertices_adjoin_other_segments_at_acute_angle(tri::Triangulation, e) -> Int, Vertex
Determines if the vertices of a segment e
of tri
adjoin other segments at an acute angle.
Arguments
tri::Triangulation
: TheTriangulation
.e
: The segment.
Output
num_adjoin
: The number of vertices ofe
that adjoin other segments at an acute angle.adjoin_vert
: The vertex ofe
that adjoins another segment at an acute angle ifnum_adjoin == 1
, and∅
otherwise.
DelaunayTriangulation.split_all_encroached_segments!
— Methodsplit_all_encroached_segments!(tri::Triangulation, args::RefinementArguments)
Splits all encroached segments of tri
according to split_subsegment!
until no more encroached segments exist in args.queue
.
DelaunayTriangulation.split_subsegment!
— Methodsplit_subsegment!(tri::Triangulation, args::RefinementArguments, e)
Splits a subsegment e
of tri
at a position determined by compute_split_position
; for curve-bounded domains, the position is determined by split_subcurve!
. After the split, assess_triangle_quality
is used to find any new bad quality triangles. Before splitting, all free vertices in the segment's diametral circle are deleted using delete_free_vertices_around_subsegment!
.
DelaunayTriangulation.split_triangle!
— Methodsplit_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.
DelaunayTriangulation.is_encroachmentfailure
— Functionis_encroachmentfailure(cert::Certificate) -> Bool
Returns true
if cert
is EncroachmentFailure
, and false
otherwise.
DelaunayTriangulation.is_successfulinsertion
— Functionis_successfulinsertion(cert::Certificate) -> Bool
Returns true
if cert
is SuccessfulInsertion
, and false
otherwise.
DelaunayTriangulation.is_failedinsertion
— Functionis_failedinsertion(cert::Certificate) -> Bool
Returns true
if cert
is FailedInsertion
, and false
otherwise.
DelaunayTriangulation.is_precisionfailure
— Functionis_precisionfailure(cert::Certificate) -> Bool
Returns true
if cert
is PrecisionFailure
, and false
otherwise.
RTree
The RTrees we work with during boundary enrichment have several functions associated with them.
DelaunayTriangulation.EnlargementValues
— TypeEnlargementValues
Type for representing the values used in the minimisation of enlargement.
Fields
enlargement::Float64
: The enlargement of the bounding box of the child being compared with.idx::Int
: The index of the child being compared with.area::Float64
: The area of the child being compared with.bounding_box::BoundingBox
: The bounding box being compared with for enlargement.
Constructor
EnlargementValues(enlargement, idx, area, bounding_box)
EnlargementValues(bounding_box)
Base.append!
— Methodappend!(node::AbstractNode, child)
Appends child
to node
's children. Also updates node
's bounding box.
Base.delete!
— Methoddelete!(tree::RTree, id_bounding_box::DiametralBoundingBox)
Deletes id_bounding_box
from tree
.
Base.insert!
— Methodinsert!(node::AbstractNode, child, tree::RTree) -> Bool
Inserts child
into node
in tree
. Returns true
if the tree
's bounding boxes had to be adjusted and false
otherwise.
Base.insert!
— Methodinsert!(tree::RTree, bounding_box[, level = 1]) -> Bool
Inserts bounding_box
into tree
. Returns true
if the tree
's bounding boxes had to be adjusted and false
otherwise.
Base.iterate
— Methoditerate(itr::RTreeIntersectionIterator, state...)
Iterate over the next state of itr
to find more intersections with the bounding box in RTreeIntersectionIterator
.
DelaunayTriangulation.collapse_after_deletion!
— Methodcollapse_after_deletion!(node::AbstractNode, tree::RTree, detached)
Condenses tree
after a deletion of one of node
's children. The detached
argument will contain the nodes that were detached from tree
during the condensing process.
DelaunayTriangulation.detach!
— Methoddetach!(node::AbstractNode, idx) -> Bool
Detaches the idx
th child of node
. Returns true
if the node
's bounding box had to be adjusted and false
otherwise.
DelaunayTriangulation.find_bounding_box
— Methodfind_bounding_box(tree::RTree, id_bounding_box::DiametralBoundingBox) -> Tuple{Leaf{Branch}, Int}
Returns the leaf node and the index in the leaf node's children that id_bounding_box
is associated with.
DelaunayTriangulation.find_subtree
— Methodfind_subtree(tree, bounding_box, level) -> Union{Branch,Leaf{Branch}}
Returns the subtree of tree
at level
that bounding_box
should be inserted into.
DelaunayTriangulation.get_intersections
— Methodget_intersections(tree::RTree, bounding_box::BoundingBox; cache_id=1) -> RTreeIntersectionIterator
Returns an RTreeIntersectionIterator
over the elements in tree
that intersect with bounding_box
. cache_id
must be 1
or 2
, and determines what cache to use for the intersection query.
DelaunayTriangulation.get_intersections
— Methodget_intersections(tree::RTree, point::NTuple{2,<:Number}; cache_id=1) -> RTreeIntersectionIterator
Returns an RTreeIntersectionIterator
over the elements in tree
that intersect with point
, representing point
as a BoundingBox
with zero width and height centered at point
. cache_id
must be 1
or 2
, and determines what cache to use for the intersection query.
DelaunayTriangulation.get_next_child
— Methodget_next_child(node::AbstractNode, start_idx, need_tests, itr::RTreeIntersectionIterator) -> Int, QueryResult
Returns the index of the next child of node
that intersects with the bounding box in itr
and the QueryResult
of the intersection.
DelaunayTriangulation.insert_detached!
— Methodinsert_detached!(tree::RTree, detached)
Given the detached
nodes from collapse_after_deletion!
, inserts them back into tree
.
DelaunayTriangulation.minimise_enlargement
— Methodminimise_enlargement(node, bounding_box) -> EnlargementValues
Returns the EnlargementValues
associated with the child of node
that minimises the enlargement, where enlargement is defined as the difference between the area of bounding_box
and the area of the child's bounding box.
DelaunayTriangulation.overflow_insert!
— Methodoverflow_insert!(node, child, tree) -> Bool
Inserts child
into node
in tree
when node
is full. Returns true
.
DelaunayTriangulation.replace!
— Methodreplace!(node::AbstractNode, left, right, original_bounding_box, tree::RTree)
Replaces the node
in tree
with left
and right
. Returns true
if the tree
's bounding boxes had to be adjusted and false
otherwise.
DelaunayTriangulation.split!
— Methodsplit!(node::AbstractNode, tree::RTree) -> Branch, Branch
Splits node
into two other nodes using a linear splitting rule. Returns the two new nodes.
DelaunayTriangulation.split_seeds
— Methodsplit_seeds(node::AbstractNode) -> NTuple{2, Int}
Returns the indices of two children in node
used to initiate the split in split!
.
DelaunayTriangulation.test_intersection
— Methodtest_intersection(node::AbstractNode, itr::RTreeIntersectionIterator) -> QueryResult
Tests whether node
intersects with the bounding box in itr
, returning a QueryResult
.
DelaunayTriangulation.update_bounding_box!
— Methodupdate_bounding_box!(node, idx, original_bounding_box, tree)
Updates the bounding box of node
to be the union of its children's bounding boxes. If node
has a parent, updates the bounding box of the parent as well if needed.
DelaunayTriangulation.update_bounding_box!
— Methodupdate_bounding_box!(node)
Updates the bounding box of node
to be the union of its children's bounding boxes.
DelaunayTriangulation.update_enlargement
— Methodupdate_enlargement(values::EnlargementValues, child, idx) -> EnlargementValues
Compare the enlargement state in values
with child
and return the updated values
if necessary.
Point Location
Here are some functions related to point location.
DelaunayTriangulation.brute_force_search
— Methodbrute_force_search(tri::Triangulation, q; itr = each_triangle(tri), predicates::AbstractPredicateKernel=AdaptiveKernel())
Searches for the triangle containing the point q
by brute force. An exception will be raised if no triangle contains the point.
See also find_triangle
.
Arguments
tri::Triangulation
: TheTriangulation
.q
: The point to be located.
Keyword Arguments
itr = each_triangle(tri)
: The iterator over the triangles of the 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.
Output
V
: The triangle containing the pointq
.
DelaunayTriangulation.brute_force_search_enclosing_circumcircle
— Functionbrute_force_search_enclosing_circumcircle(tri::Triangulation, i, predicates::AbstractPredicateKernel=AdaptiveKernel(); cache = nothing) -> Triangle
Searches for a triangle in tri
containing the vertex i
in its circumcircle using brute force. If tri
is a weighted Delaunay triangulation, the triangle returned instead has the lifted vertex i
below its witness plane. If no such triangle exists, (0, 0, 0)
is returned. You can control the method used for computing predicates via the predicates
argument.
The cache
argument is passed to [point_position_relative_to_circumcircle
] and should be one of
nothing
: No cache is used.get_incircle_cache(tri)
: The cache stored insidetri
.AdaptivePredicates.incircleadapt_cache(number_type(tri))
: Compute a new cache.
The cache is only needed if an AdaptiveKernel()
is used.
DelaunayTriangulation.check_for_intersections_with_adjacent_boundary_edges
— Methodcheck_for_intersections_with_adjacent_boundary_edges(tri::Triangulation, k, q, ghost_vertex=𝒢, predicates::AbstractPredicateKernel=AdaptiveKernel()) -> (Certificate, Certificate, Vertex, Certificate, Certificate)
Given a boundary vertex k
, find a triangle adjacent to k
to locate a triangle or edge containing q
.
See also search_down_adjacent_boundary_edges
, which uses this function to determine an initial direction to search along a straight boundary in case q
is collinear with it.
Arguments
tri::Triangulation
: TheTriangulation
.k
: The boundary vertex to start from.q
: The query point.ghost_vertex=𝒢
: The ghost vertex corresponding to the boundary thatk
resides on.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
direction_cert
: The direction ofq
relative to the vertexk
along the boundary, given as aCertificate
Left
,Right
, orOutside
. Ifis_outside(direction_cert)
, thenq
is not collinear with either of the adjacent boundary edges.q_pos_cert
: The position ofq
relative to the vertexk
along the boundary, given as aCertificate
Left
,Right
,On
,Outside
, orDegenerate
. This is similar todirection_cert
in that it will beOutside
wheneverdirection_cert
is, but this certificate can also beOn
to indicate that not only isq
in the direction given bydirection_cert
, but it is directly on the edge in that direction. Ifis_degnerate(q_pos_cert)
, thenq = get_point(tri, next_vertex)
.next_vertex
: The next vertex along the boundary in the direction ofq
, ork
ifq
is not collinear with either of the adjacent boundary edges.right_cert
: TheCertificate
for the position ofq
relative to the boundary edge right ofk
.left_cert
: TheCertificate
for the position ofq
relative to the boundary edge left ofk
.
DelaunayTriangulation.check_for_intersections_with_interior_edges_adjacent_to_boundary_vertex
— Methodcheck_for_intersections_with_interior_edges_adjacent_to_boundary_vertex(tri::Triangulation, k, q, right_cert, left_cert, predicates::AbstractPredicateKernel=AdaptiveKernel(), store_history=Val(false), history=nothing, ghost_vertex=𝒢) -> (Bool, Vertex, Vertex, Certificate, Certificate)
Checks for intersections between the line pq
, where p = get_point(tri, k)
, and the edges neighbouring p
, assuming k
is a boundary node. This function should only be used after using check_for_intersections_with_adjacent_boundary_edges
.
Arguments
tri::Triangulation
: TheTriangulation
.k
: The boundary vertex to start from.q
: The query point.right_cert
: TheCertificate
for the position ofq
relative to the boundary edge right ofk
, coming fromcheck_for_intersections_with_adjacent_boundary_edges
.left_cert
: TheCertificate
for the position ofq
relative to the boundary edge left ofk
, coming fromcheck_for_intersections_with_adjacent_boundary_edges
.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_history=Val(false)
: Whether to store the history of the algorithm.history=nothing
: The history of the algorithm. Ifstore_history
, then this should be aPointLocationHistory
object.ghost_vertex=𝒢
: The ghost vertex corresponding to the boundary thatk
resides on.
Outputs
The output takes the form (i, j, edge_cert, triangle_cert)
. Rather than defining each output individually, here are the possible froms of the output:
(i, j, Single, Outside)
: The linepq
intersects the edgepᵢpⱼ
and(j, i, k)
is a positively oriented triangle so thatpᵢ
is left ofpq
andpⱼ
is right ofpq
.(i, j, None, Inside)
: The pointq
is inside the positively oriented triangle(i, j, k)
.(0, 0, None, Outside)
: The pointq
is outside of the triangulation.(i, j, On, Inside)
: The pointq
is on the edgepᵢpⱼ
, and thus inside the positively oriented triangle(i, j, k)
.(i, j, Right, Outside)
:The point
qis collinear with the edge
pᵢpⱼ`, but is off of it and further into the triangulation.
This function assumes that the geometry is convex.
Extended help
This function works in two stages. Firstly, using check_for_intersections_with_single_interior_edge_adjacent_to_boundary_vertex
, we check for the intersection of pq
with the edges neighbouring the vertex k
, rotating counter-clockwise until we find an intersection or reach the other side of the boundary, starting from the first edge counter-clockwise away from the boundary edge right of the vertex k
. By keeping track of the positions of pq
relative to the current vertex and the previous, we can identify when an intersection is found. If no intersection is found before reaching the boundary edge left of k
, then check_for_intersections_with_triangle_left_to_boundary_vertex
is used to check the remaining triangle.
DelaunayTriangulation.compare_distance
— Methodcompare_distance(current_dist, current_idx, pts, i, qx, qy) -> (Number, Vertex)
Computes the minimum of the distance between the i
th point of pts
and (qx, qy)
and current_dist
.
Arguments
current_dist
: The current value for the distance to the point(qx, qy)
.current_idx
: The point ofpts
corresponding to the distancecurrent_dist
.pts
: The point set.i
: The vertex to compare withcurrent_idx
.qx
: The x-coordinate of the query point.qy
: The y-coordinate of the query point.
Outputs
current_dist
: The minimum of the distance between thei
th point ofpts
and(qx, qy)
andcurrent_dist
.current_idx
: The point ofpts
corresponding to the distancecurrent_dist
, which will be eitheri
orcurrent_idx
.
DelaunayTriangulation.concavity_protection_check
— Functionconcavity_protection_check(tri::Triangulation, concavity_protection, V, q, predicates::AbstractPredicateKernel=AdaptiveKernel()) -> Bool
Check whether the find_triangle
algorithm needs to restart. This is only needed if tri
is not convex.
Arguments
tri::Triangulation
: TheTriangulation
.concavity_protection
: Whether this check is needed.V
: The triangle that the algorithm has found.q
: The query point.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
need_to_restart
: Whether the algorithm needs to restart. Will also befalse
ifconcavity_protection
.
Extended help
This function uses dist
to determine whether the query point q
is inside or outside of the polygon defined by the triangulation, and also checks the position of q
relative to V
via point_position_relative_to_triangle
. If q
is outside of this triangle, then need_to_restart = true
. If q
is inside this triangle, then issues can still arise due to overlapping ghost triangles from the non-convexity. Thus, depending on the result from dist
and whether V
is a ghost triangle, need_to_restart
will be set accordingly.
DelaunayTriangulation.default_num_samples
— Methoddefault_num_samples(n) -> Integer
Given a number of points n
, returns ∛n
rounded up to the nearest integer. This is the default number of samples used in the jump-and-march algorithm.
DelaunayTriangulation.exterior_find_triangle
— Functionexterior_find_triangle(tri::Triangulation, k, q, ghost_vertex=𝒢, predicates::AbstractPredicateKernel=AdaptiveKernel()) -> (Vertex, Vertex)
Given a query point q
outside of the triangulation, find the ghost triangle containing q
.
Arguments
tri
: TheTriangulation
.k
: The exterior ghost vertex to start from.q
: The query point.ghost_vertex=𝒢
: The ghost vertex corresponding to the boundary thatk
resides on.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
i
: The first vertex of the edge on the boundary adjoining the positively oriented ghost triangle.j
: The second vertex of the edge on the boundary adjoining the positively oriented ghost triangle.
This function assumes that the geometry is convex. If the geometry is not convex, the returned value may not be correct and should be checked separately. Additionally, if q
is actually inside the triangulation, then the result is meaningless.
Extended help
This function works by first finding the position of q
relative to pₘp
, where pₘ
is the representative point for the ghost_vertex
and p = get_point(tri, k)
. Depending on this position, we rotate counter-clockwise if q
is left of the line using exterior_find_triangle_rotate_left
and clockwise otherwise using exterior_find_triangle_rotate_right
. By keeping track of the current position of q
and its position relative to the next ghost edge, we can identify when q
resides inside a ghost triangle.
DelaunayTriangulation.find_triangle
— Methodfind_triangle(tri, q; kwargs...) -> Triangle[, Bool]
Find the triangle in the triangulation tri
containing the query point q
using the jump-and-march algorithm.
For this function to work best, the triangulation should have ghost triangles, which you can add using add_ghost_triangles!
in case tri
does not already have them. Without ghost triangles, the function may not be able to find the correct triangle containing q
.
For the variables defined below, you may want to refer to the extended help which also gives some warnings and notes.
Arguments
tri
: TheTriangulation
.q
: The query point.
Keyword Arguments
predicates::AbstractPredicateKernel=ExactKernel()
: 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 vertices to consider as possible starting points for the algorithm.m=default_num_samples(num_vertices(point_indices))
: The number of samples to use when selecting the initial point.try_points=()
: A list of points to try as the initial point in addition to them
sampled.rng=Random.default_rng()
: The random number generator to use.k=select_initial_point(tri, q; point_indices, m, try_points, rng)
: The initial point to start the algorithm from. Seeselect_initial_point
.store_history=Val(false)
: Whether to store the history of the algorithm.history=nothing
: The history of the algorithm. Ifstore_history
, then this should be aPointLocationHistory
object.maxiters=2 + num_exterior_curves(tri) - num_solid_vertices(tri) + num_solid_edges(tri)
: The maximum number of iterations to perform before restarting the algorithm withrestart_find_triangle
.concavity_protection=false
: Whether to use concavity protection. Seeconcavity_protection_check
. This is only needed if your triangulation is not convex.use_barriers::Val{U}=Val(false)
: Whether to stop searching beyond any segments in the triangulation.
Output
V
: The triangle containingq
, with type given bytriangle_type(tri)
.
If you have use_barriers == Val(true)
, then we also return
invisible_flag
:false
if the triangle was found without hitting a barrier, andtrue
otherwise.
Extended help
The algorithm underlying this function is complicated and broken into many parts. Here, we describe a brief overview of the algorithm, but note that the documentation contains a much more detailed description.
- Firstly, the algorithm is initialised depending on whether
k
is a boundary or an interior vertex, usinginitialise_find_triangle_boundary_vertex
orinitialise_find_triangle_interior_vertex
respectively. - From the initial triangle
(i, j, k)
chosen, we then check ifq
is one ofpᵢ
,pⱼ
, andp = pₖ
and then return according tofind_triangle_return_on_vertex
if needed. - If we do not return above, we need to step from the initial triangle towards
q
. Since we putpᵢ
andpⱼ
to the left and right of the linepq
, respectively, this means that we step until the trianglepᵢpⱼq
is no longer positively oriented. So, while the triangle is positively oriented, we step according tofind_triangle_across_triangle
. - If we have not yet returned and the triangle is no longer positively oriented, we check if the triangle is degenerate using
find_triangle_degenerate_arrangement
and reinitialise the algorithm if needed. Otherwise, we have found the triangle containingq
and return the triangle.
Here are some additional warnings and notes for the variables defined in this function.
If the algorithm restarts, then the initial point k
is selected again using select_initial_point
, and the algorithm is restarted from there. This is done if the algorithm gets stuck in a loop, or if the algorithm is not able to find the correct triangle containing q
after maxiters
iterations. For a convex geometry, maxiters
can be safely ignored, as the sequence of triangles visited is acyclic [see H. Edelsbrunner, An acyclicity theorem for cell complexes in d dimensions, Combinatorica 10 (1990) 251-260.)].
If you are using barriers, it will be your responsibility to verify that any found triangle from this function actually contains the triangle. This can be verified using the returned flag
(see below), although the point might still be on the triangle's boundary.
If you are using barriers, it is possible that the algorithm can walk past vertices of barriers. One way this can happen is if the initial search line intersects a vertex, in which case the segment might not be considered. Another way this can happen is if you start the algorithm directly on a segment vertex, in which case the algorithm can go past it (e.g. this means that it is possible that a ghost triangle might still be returned if you start the algorithm on a boundary node).
Some notes about the output:
If a barrier is hit before any initial triangle is properly identified, the returned triangle is (0, 0, 0)
; this is only possible if use_barriers == Val(true)
. Moreover, if use_barriers == Val(true)
, the final triangle may not even be valid if invisible_flag == true
(defined below).
While this function does still work for non-convex geometries, it may be significantly slower than for convex geometries, as most of the details of the algorithm assume that the geometry is convex, and so the algorithm may have to restart many times at new initial vertices k
.
DelaunayTriangulation.find_triangle_across_triangle
— Methodfind_triangle_across_triangle(tri::Triangulation, q, k, predicates, store_history, history, maxiters, cur_iter, concavity_protection, arrangement, original_k, last_changed, p, i, j, pᵢ, pⱼ) -> (Bool, Bool, Bool, Triangle, Integer, Certificate, Vertex, Vertex, Vertex, Point, Point, Integer, Integer)
Walks across the current triangle past the edge (i, j)
. progressing the find_triangle
algorithm.
Arguments
tri::Triangulation
: TheTriangulation
.q
: The query point.k
: The vertex that the algorithm started from.predicates::AbstractPredicateKernel
: Method to use for computing predicates. Can be one ofFastKernel
,ExactKernel
, andAdaptiveKernel
. See the documentation for a further discussion of these methods.store_history
: Whether to store the history of the algorithm.history
: The history of the algorithm. Ifstore_history
, then this should be aPointLocationHistory
object.maxiters
: The maximum number of iterations to perform before restarting the algorithm withrestart_find_triangle
.cur_iter
: The current iteration of the algorithm.concavity_protection
: Whether to use concavity protection. Seeconcavity_protection_check
. This is only needed if your triangulation is not convex.arrangement
: ACertificate
defining the orientation of the trianglepᵢpⱼq
.original_k
: The original vertex that the algorithm started from.last_changed
: The last vertex that was changed in the algorithm.p
: The point corresponding to the vertexoriginal_k
.i
: The first vertex of the triangle adjoiningk
to step from.j
: The second vertex of the triangle adjoiningk
to step from.pᵢ
: The point corresponding to the vertexi
.pⱼ
: The point corresponding to the vertexj
.
Outputs
restart_flag
: Whether the algorithm needs to be restarted.return_flag
: Whether the algorithm can return immediately, returningV
.reinitialise_flag
: Whether the algorithm needs to be reinitialised at a new vertexk
. (This would only be needed if!has_ghost_triangles(tri)
.)V
: The triangle stepped into.cur_iter
: The new number of iterations of the algorithm.arrangement
: ACertificate
defining the orientation of the trianglepᵢpⱼq
with the updated values ofi
andj
.k
: The new vertex that the algorithm is at.last_changed
: The last vertex that was changed in the algorithm.original_k
: The original vertex that the algorithm started from.pᵢ
: The point corresponding to the vertexi
.pⱼ
: The point corresponding to the vertexj
.i
: The first vertex of the triangle adjoiningk
to step from.j
: The second vertex of the triangle adjoiningk
to step from.
Extended help
This part of the algorithm is relatively complicated because there are many cases that need to be accounted for. Here we give a brief description of how this step works, and note that the documentation contains a much more detailed description.
- Firstly, we need to check whether
k
is an exterior ghost vertex or not. Ifk
is an exterior ghost vertex, then this means that we are stepping outside of the triangulation. Thus, we useexterior_find_triangle
to find whereq
is, starting from thelast_changed
vertex. Ifconcavity_protection = true
, thenconcavity_protection_check
is used to determine if a restart is needed, or if we can return safely. If we reach this step but!has_ghost_triangles(tri)
, then the algorithm should need to be reinitialised sinceq
should not be outside of the triangulation, and so we return withreinitialise_flag = true
. - Now we consider the case where
k
is not an exterior ghost vertex. We move forward by updating the value ofk
so thatk = get_adjacent(tri, i, j)
, and then consider wherepₖ
is relative to the linepq
.
2a. If pₖ
is to the right of pq
, then we should update j
by j = k
, ensuring that j
is always to the right of pq
.
2b. If pₖ
is to the left of pq
, then we should update i
by i = k
, ensuring that i
is always to the left of pq
.
2c. The alternative to 2a and 2b is that pₖ
is collinear with the edge of pq
, which could mean that q
is in the current triangle or it is in a triangle further away. We compute a Certificate
that determines where q
is relative to the triangle pᵢpⱼpₖ
. If q
is inside or on this triangle, then we return, restarting if necessary according to concavity_protection
and concavity_protection_check
. If we do not yet need to return, then we need to make a decision as to which of i
and j
to update, noting that we want i
to be left of pq
and j
to be right of pq
, but this is no longer unambiguous since pₖ
is collinear with pq
. We make this decision according to last_changed
: If last_changed = i
, then moving left is what caused us to find this collinear edge, and so we send k
left by letting i = k
. Otherwise, we send k
right by letting j = k
.
- Now having stepped forward, we recompute the
Certificate
for arrangement and return, settingrestart_flag = true
ifcur_iters ≥ maxiters
.
DelaunayTriangulation.find_triangle_degenerate_arrangement
— Methodfind_triangle_degenerate_arrangement(tri::Triangulation, q, k, predicates::AbstractPredicateKernel, store_history::F, history, pᵢ, pⱼ, i, j) -> Bool
Given a degenerate arrangement pᵢpⱼq
, reinitialise the jump and march algorithm.
Arguments
tri::Triangulation
: TheTriangulation
.q
: The query point.k
: The vertex that the algorithm started from.predicates::AbstractPredicateKernel
: Method to use for computing predicates. Can be one ofFastKernel
,ExactKernel
, andAdaptiveKernel
. See the documentation for a further discussion of these methods.store_history
: Whether to store the history of the algorithm.history
: The history of the algorithm. Ifstore_history
, then this should be aPointLocationHistory
object.pᵢ
: The point corresponding to the vertexi
.pⱼ
: The point corresponding to the vertexj
.i
: The first vertex of the triangle adjoiningk
to step from.j
: The second vertex of the triangle adjoiningk
to step from.
Outputs
Bool
: Whether the algorithm needs to be restarted.
DelaunayTriangulation.find_triangle_return_on_vertex
— Methodfind_triangle_return_on_vertex(tri::Triangulation, q, k, p, pᵢ, pⱼ, i, j, predicates::AbstractPredicateKernel=AdaptiveKernel()) -> (Bool, Bool, Triangle)
Check if q
is one of the vertices of the triangle (i, j, k)
and return if needed.
Arguments
tri::Triangulation
: TheTriangulation
.q
: The query point.k
: The vertexk
that the algorithm started from.p
: The point corresponding to the vertexk
.pᵢ
: The point corresponding to the vertexi
.pⱼ
: The point corresponding to the vertexj
.i
: The first vertex of the triangle adjoiningk
to start from.j
: The second vertex of the triangle adjoiningk
to start from.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
restart_flag
: Whether the algorithm needs to be restarted.return_flag
: Whether the algorithm can return immediately, returningV
.V
: The triangle(i, j, k)
.
Extended help
An extra check is made in this algorithm for the case that the point that q
is equal to is one of the points corresponding to a ghost vertex, so it may be for example that q == pᵢ
but is_ghost_vertex(i)
, in which case the algorithm would need to restart.
DelaunayTriangulation.initialise_find_triangle_boundary_vertex
— Methodinitialise_find_triangle_boundary_vertex(tri::Triangulation, q, k, predicates::AbstractPredicateKernel, store_history, history, ghost_vertex, concavity_protection) -> (Bool, Bool, Triangle, Point, Vertex, Vertex, Point, Point)
Initialise the jump-and-march algorithm for a boundary vertex k
.
Arguments
tri::Triangulation
: TheTriangulation
.q
: The query point.k
: The boundary vertex to start from.predicates::AbstractPredicateKernel
: Method to use for computing predicates. Can be one ofFastKernel
,ExactKernel
, andAdaptiveKernel
. See the documentation for a further discussion of these methods.store_history
: Whether to store the history of the algorithm.history
: The history of the algorithm. Ifstore_history
, then this should be aPointLocationHistory
object.ghost_vertex
: The ghost vertex corresponding to the boundary thatk
resides on.concavity_protection
: Whether to use concavity protection. Seeconcavity_protection_check
. This is only needed if your triangulation is not convex.
Outputs
restart_flag
: Whether the algorithm needs to be restarted.return_flag
: Whether the algorithm can return immediately, returningV
.V
: Either a triangle that can returned ifreturn_flag = true
, or some triangle that is used for type stability for this return value.p
: The point corresponding to the vertexk
, or it may beq
if the algorithm is going to be restarted orreturn_flag = true
.i
: The first vertex of the triangle adjoiningk
to start from, ork
if the algorithm is going to be restarted orreturn_flag = true
.j
: The second vertex of the triangle adjoiningk
to start from, ork
if the algorithm is going to be restarted orreturn_flag = true
.pᵢ
: The point corresponding to the vertexi
, or it may beq
if the algorithm is going to be restarted orreturn_flag = true
.pⱼ
: The point corresponding to the vertexj
, or it may beq
if the algorithm is going to be restarted orreturn_flag = true
.
Extended help
There are multiple stages to this initialisation, starting from check_for_intersections_with_adjacent_boundary_edges
.
If it is found that
q
is not outside of the triangulation, so thatq
is collinear with one of the boundary edges, then we usesearch_down_adjacent_boundary_edges
to find where to start, noting that we can return immediately ifq
is found to be on an adjacent boundary edge. Otherwise,exterior_find_triangle
can then be used to find the ghost triangle containingq
; ifconcavity_protection = true
, thenconcavity_protection_check
is used to determine if a restart is needed.If is is not found that
q
is outside of the triangulation yet based on information from the adjacent boundary edges, then we need to check the neighbouring interior edges usingcheck_for_intersections_with_interior_edges_adjacent_to_boundary_vertex
, returning early ifq
is found to be inside one of the neighbouring triangles. If the linepq
, wherep = get_point(tri, k)
, does not intersect any of the neighbouring edges and it is not inside any of the neighbouring triangles, then it must be outside of the triangulation and so we useexterior_find_triangle
to find the triangle; as before,concavity_protection_check
is used on the found ghost triangle if needed. If there is an intersection, then we return the triangle containing the intersection point that we can start the algorithm from, and its associated vertices and points.
DelaunayTriangulation.initialise_find_triangle_interior_vertex
— Methodinitialise_find_triangle_interior_vertex(tri::Triangulation, q, k, predicates::AbstractPredicateKernel, store_history::F, history, rng) -> (Bool, Point, Vertex, Vertex, Point, Point)
Initialise the jump-and-march algorithm for an interior vertex k
.
Arguments
tri::Triangulation
: TheTriangulation
.q
: The query point.k
: The interior vertex to start from.predicates::AbstractPredicateKernel
: Method to use for computing predicates. Can be one ofFastKernel
,ExactKernel
, andAdaptiveKernel
. See the documentation for a further discussion of these methods.store_history
: Whether to store the history of the algorithm.history
: The history of the algorithm. Ifstore_history
, then this should be aPointLocationHistory
object.rng
: The random number generator to use.
Outputs
restart_flag
: Whether the algorithm needs to be restarted.p
: The point corresponding to the vertexk
.i
: The first vertex of the triangle adjoiningk
to start from.j
: The second vertex of the triangle adjoiningk
to start from.pᵢ
: The point corresponding to the vertexi
.pⱼ
: The point corresponding to the vertexj
.
Extended help
This function works by simply using select_initial_triangle_interior_vertex
to find the initial triangle to start from. A check is made to see if the edge (i, j)
refers to a non-existent edge (0, 0)
, in which case the algorithm needs to be restarted.
DelaunayTriangulation.prepare_initial_edge
— Functionprepare_initial_edge(tri::Triangulation, edges, p, q, rng::Random.AbstractRNG=Random.default_rng(), predicates::AbstractPredicateKernel=AdaptiveKernel()) -> (Vertex, Vertex, Point, Point, Certificate, Certificate)
Selects a random edge from the set of edges edges
and computes the certificates for the points corresponding to the initial and terminal vertices of the edge.
Arguments
tri::Triangulation
: TheTriangulation
.edges
: The set of edges to sample from.p
: The initial point.q
: The query point.rng::Random.AbstractRNG=Random.default_rng()
: The random number generator to use.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
i
: The initial vertex of the edge.j
: The terminal vertex of the edge.pᵢ
: The point corresponding toi
.pⱼ
: The point corresponding toj
.line_cert_i
: TheCertificate
forpᵢ
's position relative to the oriented linepq
.line_cert_j
: TheCertificate
forpⱼ
's position relative to the oriented linepq
.
DelaunayTriangulation.restart_find_triangle
— Methodrestart_find_triangle(tri::Triangulation, q, store_history, history, rng, maxiters, cur_iter, concavity_protection, num_restarts, use_barriers, predicates::AbstractPredicateKernel) -> Triangle[, Bool]
Restart the find_triangle
algorithm, or use brute_force_search
to find q
if num_restarts ≥ 25
.
Arguments
tri::Triangulation
: TheTriangulation
.q
: The query point.store_history
: Whether to store the history of the algorithm.history
: The history of the algorithm. Ifstore_history
, then this should be aPointLocationHistory
object.rng
: The random number generator to use.maxiters
: The maximum number of iterations to perform before restarting the algorithm withrestart_find_triangle
.cur_iter
: The current iteration of the algorithm.concavity_protection
: Whether to use concavity protection. Seeconcavity_protection_check
. This is only needed if your triangulation is not convex.num_restarts
: The number of times the algorithm has been restarted.use_barriers
: Whether to use barriers, stopping the algorithm at any segment.predicates::AbstractPredicateKernel
: Method to use for computing predicates. Can be one ofFastKernel
,ExactKernel
, andAdaptiveKernel
. See the documentation for a further discussion of these methods.
Outputs
V
: The triangle containingq
.
In addition, if use_barriers = Val(true)
, then a second output is returned, which is a boolean indicating whether the algorithm reached a barrier (true
) or not (false
).
DelaunayTriangulation.search_down_adjacent_boundary_edges
— Methodsearch_down_adjacent_boundary_edges(tri::Triangulation, k, q, direction_cert, q_pos_cert, next_vertex, predicates::AbstractPredicateKernel=AdaptiveKernel(), store_history=Val(false), history=nothing, ghost_vertex=𝒢) -> (Bool, Certificate, Vertex, Vertex, Vertex)
Starting at the boundary vertex k
, walks down the boundary in the direction of q
until finding q
or finding that it is outside of the triangulation.
See also check_for_intersections_with_adjacent_boundary_edges
, which uses this function to determine an initial direction to search along.
Arguments
tri::Triangulation
: TheTriangulation
.k
: The boundary vertex to start from.q
: The query point.direction_cert
: The direction ofq
relative to the vertexk
along the boundary, defined fromcheck_for_intersections_with_adjacent_boundary_edges
.q_pos_cert
: The position ofq
relative to the vertexk
along the boundary, defined fromcheck_for_intersections_with_adjacent_boundary_edges
.next_vertex
: The next vertex along the boundary in the direction ofq
, defined fromcheck_for_intersections_with_adjacent_boundary_edges
.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_history=Val(false)
: Whether to store the history of the algorithm.history=nothing
: The history of the algorithm. Ifstore_history
, then this should be aPointLocationHistory
object.ghost_vertex=𝒢
: The ghost vertex corresponding to the boundary thatk
resides on.
Outputs
return_flag
: Whether to return, or throw an exception.q_pos_cert
: ACertificate
that isOn
ifq
is on the edge(u, v)
, andOutside
ifq
is outside of the triangulation.u
: Ifis_on(q_pos_cert)
, this is the first vertex of a positively oriented triangle thatq
is on, so thatq
is on the edge(u, v)
. Otherwise,(u, v, w)
is a ghost triangle close toq
.v
: Ifis_on(q_pos_cert)
, this is the second vertex of a positively oriented triangle thatq
is on, so thatq
is on the edge(u, v)
. Otherwise,(u, v, w)
is a ghost triangle close toq
.w
: Ifis_on(q_pos_cert)
, this is the third vertex of a positively oriented triangle thatq
is on, so thatq
is on the edge(u, v)
andw = get_adjacent(tri, u, v)
. Otherwise,(u, v, w)
is a ghost triangle close toq
.
This function assumes that the geometry is convex. The function will still be able to return, but is_outside(q_pos_cert)
may not necessarily mean q
is outside of the triangulation. The main function find_triangle
will have to restart the algorithm if it is found that is_outside(q_pos_cert)
was incorrect.
Extended help
This function works by stepping along vertices on the boundaries in the direction specified by direction_cert
, using search_right_down_adjacent_boundary_edges
if is_right(direction_cert)
and search_left_down_adjacent_boundary_edges
otherwise. In these functions, a while
loop is used to keep stepping until q_pos_cert
, which is updated at each iteration, changes value.
DelaunayTriangulation.select_initial_point
— Functionselect_initial_point(tri::Triangulation, q; kwargs...) -> Vertex
Selects the initial point for find_triangle
to start from.
Arguments
tri::Triangulation
: TheTriangulation
.q
: The query point. Can be either a point or a vertex - if it is a vertex, the corresponding pointget_point(tri, q)
will be used.
Keyword Arguments
point_indices=each_solid_vertex(tri)
: The indices to sample from.m=default_num_samples(num_vertices(point_indices))
: The number of samples to take. Replacement is not used, so there may be duplicates.try_points=()
: A list of points to try in addition to those randomly sampled.allow_boundary_points=!is_disjoint(tri)
: Whether to allow boundary points to be selected.rng=Random.default_rng()
: The random number generator to use.
Outputs
i
: The index of the point closest toq
out of those queried.
DelaunayTriangulation.select_initial_triangle_interior_vertex
— Methodselect_initial_triangle_interior_vertex(tri::Triangulation,
k,
q,
predicates::AbstractPredicateKernel=AdaptiveKernel(),
store_history=Val(false),
history=nothing,
rng::Random.AbstractRNG=Random.default_rng()) -> (Point, Vertex, Vertex, Point, Point)
Selects the initial triangle for find_triangle
to start from, for the case where k
is an interior vertex.
Arguments
tri::Triangulation
: TheTriangulation
.k
: The vertex to start from.q
: The query point.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_history=Val(false)
: Whether to store the history of the algorithm.history=nothing
: The history of the algorithm. Ifstore_history
, then this should be aPointLocationHistory
object.rng::Random.AbstractRNG=Random.default_rng()
: The random number generator to use.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
p
: The point corresponding tok
.i
: The initial vertex of the triangle.j
: The terminal vertex of the triangle.pᵢ
: The point corresponding toi
.pⱼ
: The point corresponding toj
.
This function assumes that the geometry is convex. To deal with this, when an infinite loop is detected we return ∅
for both i
and j
, and then let find_triangle
handle how to correct the algorithm from there.
Extended help
This part of the algorithm works by rotating around the vertex k
, looking for a triangle whose edges adjoining k
are to the left and to the right of k
. By choosing the initial edge at random via prepare_initial_edge
, and computing the position of q
relative to this initial edge, the rotation will be either clockwise or counter-clockwise, and the triangle is then found using either select_initial_triangle_clockwise
or select_initial_triangle_counterclockwise
, respectively.
In case the initial edge is collinear with the line pq
, where p = get_point(tri, q)
, then fix_initial_collinear_edge_for_interior_vertex
to find a non-collinear edge resample more edges from prepare_initial_edge
if necessary.
DelaunayTriangulation.select_random_edge
— Functionselect_random_edge(tri::Triangulation, edges, rng::Random.AbstractRNG=Random.default_rng()) -> (Vertex, Vertex, Point, Point)
Selects a random edge from the set of edges edges
.
Arguments
tri::Triangulation
: TheTriangulation
.edges
: The set of edges to sample from.rng::Random.AbstractRNG=Random.default_rng()
: The random number generator to use.
Outputs
i
: The initial vertex of the edge.j
: The terminal vertex of the edge.pᵢ
: The point corresponding toi
.pⱼ
: The point corresponding toj
.
Triangulation Operations
Some of the triangulation operations have internal functions associated with them.
DelaunayTriangulation.delete_all_exterior_triangles!
— Methoddelete_all_exterior_triangles!(tri::Triangulation, triangles)
Deletes all the triangles in the set triangles
from the triangulation tri
.
DelaunayTriangulation.delete_holes!
— Methoddelete_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.find_all_points_to_delete!
— Methodfind_all_points_to_delete!(points_to_process, tri::Triangulation, seed, all_bn)
Starting at seed
, finds more points to spread to and mark for deletion.
Arguments
points_to_process
: The current list of points marked for deletion.tri::Triangulation
: TheTriangulation
.seed
: The seed vertex to start spreading from.all_bn
: All the boundary nodes in the triangulation, obtained fromget_all_boundary_nodes
.
Outputs
There are no outputs, as points_to_process
is updated in-place.
Extended help
This function works by considering the neighbours around the vertex seed
. For each neighbouring vertex, we designate that as a new seed, and consider if it needs to be added into points_to_process
according to its distance from the triangulation computed from distance_to_polygon
. We then call find_all_points_to_delete!
recursively again on the new seed.
DelaunayTriangulation.find_all_points_to_delete
— Methodfind_all_points_to_delete(tri::Triangulation) -> Set{Int}
Returns a set of all the points that are in the exterior faces of the triangulation tri
.
Extended help
This function works by 'spreading' from some initial vertex. In particular, starting at each boundary node, we spread outwards towards adjacent vertices, recursively spreading so that all exterior points are identified with the help of find_all_points_to_delete!
.
DelaunayTriangulation.find_all_triangles_to_delete
— Methodfind_all_triangles_to_delete(tri::Triangulation, points_to_process) -> Set{Triangle}
Returns a set of all the triangles that are in the exterior faces of the triangulation tri
.
Arguments
tri::Triangulation
: TheTriangulation
.points_to_process
: The set of points that are in the exterior faces of the triangulationtri
, obtained fromfind_all_points_to_delete
.
Outputs
triangles_to_delete
: The set of triangles that are in the exterior faces of the triangulationtri
.
Extended help
This function works in two stages.
- Firstly, all the non-boundary vertices, i.e. those from
points_to_process
, are processed. For each vertexv
, the triangles adjoining it, given byget_adjacent2vertex(tri, v)
, aremarked for deletion. - Next, all the boundary vertices need to be processed and carefully analysed to determine if any other triangles need to be deleted since, for example, a triangle may be adjoining three vertices that are all boundary vertices, and it might not be obvious if it is inside or outside of the triangulation. By applying
dist
to compute the distance between the triangle's centroid and the triangulation, the triangle can be accurately marked for deletion if it is outside of the triangulation.
DelaunayTriangulation.check_delete_point_args
— Methodcheck_delete_point_args(tri::Triangulation, vertex, S) -> Bool
Checks that the vertex vertex
can be deleted from the triangulation tri
. Returns true
if so, and throws an InvalidVertexDeletionError
otherwise. This will occur if:
vertex
is a boundary node oftri
.vertex
is a ghost vertex oftri
.vertex
adjoins a segment oftri
.
DelaunayTriangulation.delete_point!
— Methoddelete_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.fix_edges_after_deletion!
— Methodfix_edges_after_deletion!(tri::Triangulation, S)
Ensures that the edges in S
surrounding a deleted vertex of tri
are correctly updated.
DelaunayTriangulation.get_surrounding_polygon
— Methodget_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.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.get_edges_for_split_edge
— Methodget_edges_for_split_edge(tri::Triangulation, i, j, r)
Returns the edges (i, j)
, (j, i)
, (i, r)
, (r, i)
, (r, j)
, and (j, r)
.
DelaunayTriangulation.legalise_split_edge!
— Functionlegalise_split_edge!(tri::Triangulation, i, j, k, r, store_event_history=Val(false), event_history=nothing; predicates::AbstractPredicateKernel=AdaptiveKernel())
Legalises the newly added edges in tri
after the edge (i, j)
was split using split_edge!
.
See also complete_split_edge_and_legalise!
.
Arguments
tri::Triangulation
: TheTriangulation
.i
: The first vertex of the edge that was split.j
: The second vertex of the edge that was split.k
: The vertex that was originally adjacent to(i, j)
.r
: The vertex that(i, j)
was split 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_edge!
— Functionsplit_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_triangle_and_legalise!
— Methodcomplete_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.
DelaunayTriangulation.legalise_split_triangle!
— Methodlegalise_split_triangle!(tri::Triangulation, i, j, k, r; predicates::AbstractPredicateKernel=AdaptiveKernel())
Legalises the newly added edges in tri
after the triangle (i, j, k)
was split using split_triangle!
.
See also 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.
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!
— Methodsplit_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._safe_get_adjacent
— Method_safe_get_adjacent(tri::Triangulation, uv) -> Vertex
This is the safe version of get_adjacent
, which is used when the triangulation has multiple sections, ensuring that the correct ghost vertex is returned in case uv
is a ghost edge.
DelaunayTriangulation.add_adjacent!
— Methodadd_adjacent!(tri::Triangulation, uv, w)
add_adjacent!(tri::Triangulation, u, v, w)
Adds the key-value pair (u, v) ⟹ w
to the adjacency map of tri
.
DelaunayTriangulation.delete_adjacent!
— Methoddelete_adjacent!(tri::Triangulation, uv)
delete_adjacent!(tri::Triangulation, u, v)
Deletes the key (u, v)
from the adjacency map of tri
.
DelaunayTriangulation.get_adjacent
— Methodget_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.
DelaunayTriangulation.add_adjacent2vertex!
— Methodadd_adjacent2vertex!(tri::Triangulation, w, uv)
add_adjacent2vertex!(tri::Triangulation, w, u, v)
Adds the edge (u, v)
into the set of edges returned by get_adjacent2vertex(tri, w)
.
DelaunayTriangulation.delete_adjacent2vertex!
— Methoddelete_adjacent2vertex!(tri::Triangulation, w, uv)
delete_adjacent2vertex!(tri::Triangulation, w, u, v)
Deletes the edge (u, v)
from the set of edges returned by get_adjacent2vertex(tri, w)
.
DelaunayTriangulation.delete_adjacent2vertex!
— Methoddelete_adjacent2vertex!(tri::Triangulation, w)
Deletes the key w
from the Adjacent2Vertex
map of tri
.
DelaunayTriangulation.get_adjacent2vertex
— Methodget_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.add_neighbour!
— Methodadd_neighbour!(tri::Triangulation, u, v...)
Adds the neighbours v...
to u
in the graph of tri
.
DelaunayTriangulation.add_vertex!
— Methodadd_vertex!(tri::Triangulation, u...)
Adds the vertices u...
into the graph of tri
.
DelaunayTriangulation.delete_ghost_vertices_from_graph!
— Methoddelete_ghost_vertices_from_graph!(tri::Triangulation)
Deletes all ghost vertices from the graph of tri
.
DelaunayTriangulation.delete_neighbour!
— Methoddelete_neighbour!(tri::Triangulation, u, v...)
Deletes the neighbours v...
from u
in the graph of tri
.
DelaunayTriangulation.delete_vertex!
— Methoddelete_vertex!(tri::Triangulation, u...)
Deletes the vertices u...
from the graph of tri
.
DelaunayTriangulation.get_neighbours
— Methodget_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_neighbours
— Methodget_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.
DelaunayTriangulation.has_ghost_vertices
— Methodhas_ghost_vertices(tri::Triangulation) -> Bool
Returns true
if tri
has ghost vertices, and false
otherwise.
DelaunayTriangulation.has_vertex
— Methodhas_vertex(tri::Triangulation, u) -> Bool
Returns true
if u
is a vertex in tri
, and false
otherwise.
DelaunayTriangulation.num_edges
— Methodnum_edges(tri::Triangulation) -> Integer
Returns the number of edges in tri
. Note that, if has_ghost_triangles(tri)
, then some of these edges will be ghost edges.
See also num_solid_edges
and num_ghost_edges
.
DelaunayTriangulation.num_ghost_edges
— Methodnum_ghost_edges(tri::Triangulation) -> Integer
Returns the number of ghost edges in tri
.
See also num_solid_edges
and num_edges
.
DelaunayTriangulation.num_ghost_vertices
— Methodnum_ghost_vertices(tri::Triangulation) -> Integer
Returns the number of ghost vertices in tri
.
See also num_solid_vertices
and num_vertices
.
DelaunayTriangulation.num_neighbours
— Methodnum_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.num_solid_edges
— Methodnum_solid_edges(tri::Triangulation) -> Integer
Returns the number of solid edges in tri
.
See also num_ghost_edges
and num_edges
.
DelaunayTriangulation.num_solid_vertices
— Methodnum_solid_vertices(tri::Triangulation) -> Integer
Returns the number of solid vertices in tri
.
See also num_ghost_vertices
and num_vertices
.
DelaunayTriangulation.num_vertices
— Methodnum_vertices(tri::Triangulation) -> Integer
Returns the number of vertices in tri
. Note that, if has_ghost_triangles(tri)
, then some of these vertices will be ghost vertices.
See also num_solid_vertices
and num_ghost_vertices
.
DelaunayTriangulation.sort_edge_by_degree
— Methodsort_edge_by_degree(tri::Triangulation, e) -> Edge
Returns the edge e
sorted so that initial(e)
has the smaller degree of the two vertices.
DelaunayTriangulation.get_point
— Methodget_point(tri::Triangulation, i) -> NTuple{2, Number}
get_point(tri::Triangulation, i...) -> NTuple{length(i), NTuple{2, Number}}
Returns the coordinates corresponding to the vertices i...
of tri
, given as a Tuple
of the form (x, y)
for each point. If i
is a ghost vertex, then the coordinates of the representative point of the curve associated with i
are returned instead.
DelaunayTriangulation.num_points
— Methodnum_points(tri::Triangulation) -> Integer
Returns the number of points in tri
.
If tri
has vertices that are not yet present in the triangulation, e.g. if you have deleted vertices or have some submerged vertices in a weighted triangulation, then the corresponding points will still be counted in this function. It is recommended that you instead consider num_vertices
, num_solid_vertices
, or num_ghost_vertices
.
DelaunayTriangulation.pop_point!
— Methodpop_point!(tri::Triangulation)
Pops the last point from the points of tri
.
DelaunayTriangulation.push_point!
— Methodpush_point!(tri::Triangulation, x, y)
push_point!(tri::Triangulation, p)
Pushes the point p = (x, y)
into the points of tri
.
DelaunayTriangulation.set_point!
— Methodset_point!(tri::Triangulation, i, x, y)
set_point!(tri::Triangulation, i, p)
Sets the i
th point of tri
to p = (x, y)
.
DelaunayTriangulation.construct_positively_oriented_triangle
— Functionconstruct_triangle(tri::Triangulation, i, j, k, predicates::AbstractPredicateKernel=AdaptiveKernel()) -> Triangle
Returns a triangle in tri
from the vertices i
, j
, and k
such that the triangle is positively oriented.
You can use the predicates
argument to determine how the orientation predicate is computed. Can be one of FastKernel
, ExactKernel
, and AdaptiveKernel
. See the documentation for more discussion on these choices.
DelaunayTriangulation.contains_triangle
— Methodcontains_triangle(tri::Triangulation, T) -> (Triangle, Bool)
contains_triangle(tri::Triangulation, i, j, k) -> (Triangle, Bool)
Tests whether tri
contains T = (i, j, k)
up to rotation, returning
V
: The rotated form ofT
that is intri
, or simplyT
ifT
is not intri
.flag
:true
ifT
is intri
, andfalse
otherwise.
DelaunayTriangulation.num_ghost_triangles
— Methodnum_ghost_triangles(tri::Triangulation) -> Integer
Returns the number of ghost triangles in tri
.
DelaunayTriangulation.num_solid_triangles
— Methodnum_solid_triangles(tri::Triangulation) -> Integer
Returns the number of solid triangles in tri
.
DelaunayTriangulation.num_triangles
— Methodnum_triangles(tri::Triangulation) -> Integer
Returns the number of triangles in tri
. Note that, if has_ghost_triangles(tri)
, then some of these triangles will be ghost triangles.
Voronoi Tessellations
Here are some functions related to the computation of unbounded Voronoi tessellations.
DelaunayTriangulation.add_edge_to_voronoi_polygon!
— Methodadd_edge_to_voronoi_polygon!(B, vorn::VoronoiTessellation, i, k, S, m, encountered_duplicate_circumcenter) -> (Vertex, Bool, Vertex)
Add the next edge to the Voronoi polygon for the point i
in the VoronoiTessellation
vorn
.
Arguments
B
: The vector of circumcenters defining the polygon.vorn
: TheVoronoiTessellation
.i
: The polygon index.k
: The vertex to add.S
: The surrounding polygon ofi
. Seeget_surrounding_polygon
.m
: The index of the next vertex inS
.encountered_duplicate_circumcenter
: Whether or not a duplicate circumcenter has been encountered.
Outputs
ci
: The index for the circumcenter of the triangle considered.encountered_duplicate_circumcenter
: Whether or not a duplicate circumcenter has been encountered.k
: The next vertex inS
after the inputk
.
DelaunayTriangulation.add_voronoi_polygon!
— Methodadd_voronoi_polygon!(vorn::VoronoiTessellation, i) -> Vector
Add the Voronoi polygon for the point i
to the VoronoiTessellation
vorn
.
Arguments
vorn
: TheVoronoiTessellation
.i
: The polygon index.
Outputs
B
: The vector of circumcenters defining the polygon. This is a circular vector, i.e.B[begin] == B[end]
.
DelaunayTriangulation.close_voronoi_polygon!
— Methodclose_voronoi_polygon!(vorn::VoronoiTessellation, B, i, encountered_duplicate_circumcenter, prev_ci)
Close the Voronoi polygon for the point i
in the VoronoiTessellation
vorn
.
Arguments
vorn
: TheVoronoiTessellation
.B
: The vector of circumcenters defining the polygon.i
: The polygon index.encountered_duplicate_circumcenter
: Whether or not a duplicate circumcenter has been encountered.prev_ci
: The previous circumcenter index.
Outputs
There are no outputs, as vorn
and B
are modified in-place.
DelaunayTriangulation.connect_circumcenters!
— Methodconnect_circumcenters!(B, ci)
Add the circumcenter index ci
to the array B
.
DelaunayTriangulation.get_next_triangle_for_voronoi_polygon
— Methodget_next_triangle_for_voronoi_polygon(vorn::VoronoiTessellation, i, k, S, m) -> (Vertex, Vertex)
Get the next triangle for the Voronoi polygon for the point i
in the VoronoiTessellation
.
Arguments
vorn
: TheVoronoiTessellation
.i
: The polygon index.k
: The vertex to add.S
: The surrounding polygon ofi
. Seeget_surrounding_polygon
.m
: The index of the next vertex inS
.
Outputs
ci
: The index for the circumcenter of the next triangle.k
: The next vertex inS
after the inputk
.
DelaunayTriangulation.initialise_voronoi_tessellation
— Methodinitialise_voronoi_tessellation(tri::Triangulation) -> VoronoiTessellation
Initialise a VoronoiTessellation
from the triangulation tri
.
Arguments
tri
: TheTriangulation
.
Output
vorn
: TheVoronoiTessellation
. This tessellation is not yet filled in, as all the polygons and other fields need to be properly defined. This simply defines all the initial objects that will be pushed into.
DelaunayTriangulation.prepare_add_voronoi_polygon
— Methodprepare_add_voronoi_polygon(vorn::VoronoiTessellation, i) -> (Vector, Vector)
Prepare to add a Voronoi polygon for the vertex i
to the Voronoi tessellation vorn
.
Arguments
vorn
: TheVoronoiTessellation
.i
: The vertex.
Outputs
S
: The surrounding polygon ofi
. Seeget_surrounding_polygon
.B
: The buffer for the circumcenters. This is an emptyVector{I}
, whereI = integer_type(tri)
.
Clipped Voronoi Tessellations
Here are some functions related to the computation of clipped Voronoi tessellations.
DelaunayTriangulation.add_all_boundary_polygons!
— Methodadd_all_boundary_polygons!(vorn::VoronoiTessellation, boundary_sites)
Add all of the boundary polygons to the Voronoi tessellation.
Arguments
vorn
: TheVoronoiTessellation
.boundary_sites
: A dictionary of boundary sites.
Outputs
There are no outputs, but the boundary polygons are added in-place.
DelaunayTriangulation.add_intersection_points!
— Methodadd_intersection_points!(vorn::VoronoiTessellation, segment_intersections) -> Integer
Adds all of the segment_intersections
into the polygon vertices of vorn
.
Arguments
vorn
: TheVoronoiTessellation
.segment_intersections
: The intersection points fromfind_all_intersections
.
Outputs
n
: The number of polygon vertices before the intersections were added.
DelaunayTriangulation.add_segment_intersection!
— Methodadd_segment_intersection!(segment_intersections, boundary_sites, intersection_point, incident_polygon::I) where {I} -> Integer
Adds the intersection_point
into the list of segment_intersections
.
Arguments
segment_intersections
: The list of segment intersections.boundary_sites
: A mapping from boundary sites to the indices of the segment intersections that are incident to the boundary site.intersection_point
: The intersection point to add.incident_polygon
: The index of the polygon that is incident to the intersection point.
Outputs
idx
: The index of the intersection point in the list of segment intersections. If the intersection point already exists in the list, then the index of the existing point is returned and used instead.
DelaunayTriangulation.add_to_intersected_edge_cache!
— Methodadd_to_intersected_edge_cache!(intersected_edge_cache, u, v, a, b)
Add the edge uv
to the list of intersected edges.
Arguments
intersected_edge_cache
: The list of intersected edges.u
: The first vertex of the edge of the Voronoi polygon intersecting the edgeab
of the boundary.v
: The second vertex of the edge of the Voronoi polygon intersecting the edgeab
of the boundary.a
: The first vertex of the edge of the boundary.b
: The second vertex of the edge of the boundary.
Outputs
There are no outputs, as intersected_edge_cache
is modified in-place.
DelaunayTriangulation.classify_intersections!
— Methodclassify_intersections!(intersected_edge_cache, left_edge_intersectors, right_edge_intersectors, current_edge_intersectors, left_edge, right_edge, current_edge)
Classify the intersections in intersected_edge_cache
into left_edge_intersectors
, right_edge_intersectors
, and current_edge_intersectors
based on whether they intersect left_edge
, right_edge
, or current_edge
, respectively.
Arguments
intersected_edge_cache
: The list of intersected edges currently being considered.left_edge_intersectors
: The set of sites that intersect the edge to the left of an edge currently being considered.right_edge_intersectors
: The set of sites that intersect the edge to the right of an edge currently being considered.current_edge_intersectors
: The set of sites that intersect the current edge being considered.left_edge
: The edge to the left ofe
on the boundary.right_edge
: The edge to the right ofe
on the boundary.current_edge
: The edge on the boundary being considered.
Outputs
There are no outputs, but left_edge_intersectors
, right_edge_intersectors
, or current_edge_intersectors
are updated all in-place depending on the type of intersection for each edge in intersected_edge_cache
.
DelaunayTriangulation.cleanup_unbounded_polygons!
— Methodcleanup_unbounded_polygons!(vorn::VoronoiTessellation)
Removes any remaining unbounded polygons from vorn
after clipping. This is needed only for power diagrams.
DelaunayTriangulation.clip_all_polygons!
— Methodclip_all_polygons!(vorn::VoronoiTessellation, n, boundary_sites, exterior_circumcenters, equal_circumcenter_mapping)
Clip all of the polygons in the Voronoi tessellation.
Arguments
vorn
: TheVoronoiTessellation
.n
: The number of vertices in the tessellation before clipping.boundary_sites
: A dictionary of boundary sites.exterior_circumcenters
: Any exterior circumcenters to be filtered out.equal_circumcenter_mapping
: A mapping from the indices of the segment intersections that are equal to the circumcenter of a site to the index of the site.
Outputs
There are no outputs, but the polygons are clipped in-place.
DelaunayTriangulation.clip_polygon!
— Methodclip_polygon!(vorn::VoronoiTessellation, n, points, polygon, new_verts, exterior_circumcenters, equal_circumcenter_mapping)
Clip the polygon polygon
by removing the vertices that are outside of the domain and adding the new vertices new_verts
to the polygon.
Arguments
vorn
: TheVoronoiTessellation
.n
: The number of vertices in the tessellation before clipping.points
: The polygon points of the tessellation.polygon
: The index of the polygon to be clipped.new_verts
: The indices of the new vertices that are added to the polygon.exterior_circumcenters
: Any exterior circumcenters to be filtered out.equal_circumcenter_mapping
: A mapping from the indices of the segment intersections that are equal to the circumcenter of a site to the index of the site.
Outputs
There are no outputs, but the polygon is clipped in-place.
DelaunayTriangulation.clip_voronoi_tessellation!
— Methodclip_voronoi_tessellation!(vorn::VoronoiTessellation; rng=Random.default_rng(), predicates::AbstractPredicateKernel=AdaptiveKernel())
Clip the Voronoi tessellation vorn
to the convex hull of the generators in vorn
.
Arguments
vorn
: TheVoronoiTessellation
.
Keyword Arguments
clip_polygon=nothing
: Ifclip=true
, then this is the polygon to clip the Voronoi tessellation to. Ifnothing
, the convex hull of the triangulation is used. The polygon should be defined as aTuple
of the form(points, boundary_nodes)
where theboundary_nodes
are vertices mapping to coordinates inpoints
, adhering to the usual conventions for defining boundaries. Must be a convex polygon.rng::Random.AbstractRNG=Random.default_rng()
: The random number generator.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 are no outputs, but the Voronoi tessellation is clipped in-place.
DelaunayTriangulation.convert_to_edge_adjoining_ghost_vertex
— Methodconvert_to_edge_adjoining_ghost_vertex(tri::Triangulation, e) -> Edge
Returns the edge e
if it is not a boundary edge, and the edge reverse(e)
if it is a boundary edge.
See also is_boundary_edge
.
DelaunayTriangulation.dequeue_and_process!
— Functiondequeue_and_process!(vorn, polygon_edge_queue, edges_to_process,
intersected_edge_cache, left_edge_intersectors, right_edge_intersectors, current_edge_intersectors,
processed_pairs, boundary_sites, segment_intersections, exterior_circumcenters, equal_circumcenter_mapping,
rng::Random.AbstractRNG=Random.default_rng(), predicates::AbstractPredicateKernel=AdaptiveKernel())
Dequeue an edge from polygon_edge_queue
and process it. If polygon_edge_queue
is empty, then we process the first edge in edges_to_process
.
Arguments
vorn
: TheVoronoiTessellation
.polygon_edge_queue
: The queue of edges that need to be processed.edges_to_process
: The edges that need to be processed.intersected_edge_cache
: A cache of intersected edges.left_edge_intersectors
: The intersection points ofleft_edge
with other edges.right_edge_intersectors
: The intersection points ofright_edge
with other edges.current_edge_intersectors
: The intersection points ofcurrent_edge
with other edges.processed_pairs
: A set of pairs of edges and polygons that have already been processed.boundary_sites
: A dictionary of boundary sites.segment_intersections
: A dictionary of segment intersections.exterior_circumcenters
: A dictionary of exterior circumcenters.equal_circumcenter_mapping
: A mapping from the indices of the segment intersections that are equal to the circumcenter of a site to the index of the site.rng::Random.AbstractRNG=Random.default_rng()
: The random number generator.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 are no outputs. Instead, the caches and queues are updated in-place.
Extended help
This function works as follows:
- Firstly, if there are no edges queued in
polygon_edge_queue
, then we enqueue the first in edge inedges_to_process
usingenqueue_new_edge!
. - We then dequeue the next edge to be processed. If the edge has already been processed, then we return early.
- If we're still here, then we process the
(edge, polygon)
pair enqueued frompolygon_edge_queue
usingprocess_polygon!
. This function checks for intersections of theedge
with thepolygon
. - Once the polygon has been processed, we then needed to classify all of the intersections using
classify_intersections!
, which determines, for each intersection, if the intersection is withedge
, or with the edge left ofedge
, or to the edge right ofedge
. - Then,
process_intersection_points!
is used to process the intersection points, enqueueing new edges when needed. - We then delete the edge from
edges_to_process
if it is in there and return.
DelaunayTriangulation.enqueue_new_edge!
— Functionenqueue_new_edge!(polygon_edge_queue, vorn::VoronoiTessellation, e)
Enqueue the edge e
of the boundary to be processed.
Arguments
polygon_edge_queue
: The queue of edges that are to be processed.vorn
: TheVoronoiTessellation
.e
: The edge to be processed.rng::AbstractRNG=Random.default_rng()
: Random number generator. Needed forget_nearest_neighbour
.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. Needed forget_nearest_neighbour
.
Outputs
There are no outputs, as polygon_edge_queue
is modified in-place.
DelaunayTriangulation.find_all_exterior_circumcenters
— Methodfind_all_exterior_circumcenters(vorn::VoronoiTessellation, clip_points, clip_vertices) -> Set{I}
Finds all the polygon vertices in vorn
that are outside of the polygon defined by (clip_points, clip_vertices)
. The return is the set of all exterior vertices.
DelaunayTriangulation.find_all_intersecting_polygons
— Methodfind_all_intersecting_polygons(vorn::VoronoiTessellation, clip_points, clip_vertices, exterior_circumcenters) -> (Set{I}, Set{I})
Classifies all polygons in vorn
into three sets depending on how they intersect the polygon defined by (clip_points, clip_vertices)
, returning (interior, intersecting, exterior)
, where
interior
are the polygons that are entirely inside the polygon defined by(clip_points, clip_vertices)
,possibly_intersecting
are the polygons that might intersect the polygon defined by(clip_points, clip_vertices)
.
The exterior_circumcenters
are the circumcenters that are outside of the domain, see find_all_exterior_circumcenters
.
DelaunayTriangulation.find_all_intersections
— Methodfind_all_intersections(vorn::VoronoiTessellation; rng=Random.default_rng(), predicates::AbstractPredicateKernel=AdaptiveKernel()) -> (Dict, Vector, Set, Dict)
Find all intersections between the edges of the Voronoi tessellation and the boundary of the polygon.
Arguments
vorn
: TheVoronoiTessellation
.
Keyword Arguments
rng::Random.AbstractRNG=Random.default_rng()
: The random number generator.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
boundary_sites
: A dictionary of boundary sites.segment_intersections
: The intersection points.exterior_circumcenters
: The circumcenters that are outside of the domain.equal_circumcenter_mapping
: A mapping from the indices of the segment intersections that are equal to the circumcenter of a site to the index of the site.
Extended help
This algorithm works as follows:
- First, using
initialise_clipping_arrays
, we initialise the arrays that we will use to store the intersections, and queue up all boundary edges for processing. - Then, starting with the first edge in
edges_to_process
, we dequeue an edge frompolygon_edge_queue
and process it viadequeue_and_process!
. - We repeat step 2 until
polygon_edge_queue
andedges_to_process
are both empty. - We then return.
DelaunayTriangulation.fix_no_intersections!
— Methodfix_no_intersections!(vorn::VoronoiTessellation; predicates::AbstractPredicateKernel=AdaptiveKernel())
In the case that the Voronoi tiles have no intersections at all with the convex hull, this function adds in the missing intersections. clip_polygon
is used for this.
DelaunayTriangulation.get_neighbouring_boundary_edges
— Methodget_neighbouring_boundary_edges(tri::Triangulation, e) -> (Edge, Edge)
Returns the two boundary edges adjacent to the boundary edge e
in the triangulation tri
.
Arguments
tri::Triangulation
: a triangulatione
: The boundary edge.
Outputs
left_e
: The left edge.right_e
: The right edge.
DelaunayTriangulation.get_shared_vertex
— Methodget_shared_vertex(e, f) -> Vertex
Returns the vertex shared by the edges e
and f
, or ∅
if they do not share a vertex.
Arguments
e
: The first edge.f
: The second edge.
Outputs
u
: The shared vertex.
Example
julia> using DelaunayTriangulation
julia> DelaunayTriangulation.get_shared_vertex((1, 3), (5, 7))
0
julia> DelaunayTriangulation.get_shared_vertex((1, 3), (3, 7))
3
julia> DelaunayTriangulation.get_shared_vertex((10, 3), (10, 5))
10
julia> DelaunayTriangulation.get_shared_vertex((9, 4), (9, 5))
9
DelaunayTriangulation.initialise_clipping_arrays
— Methodinitialise_clipping_arrays(vorn::VoronoiTessellation) -> (Set{E}, Queue{Tuple{E,I}}, Dict{I,Set{I}}, NTuple{2,F}[], Set{Tuple{E,I}}, Pair{E,E}[], Set{I}, Set{E}, Set{E}, Set{E}, Dict{I,I})
Initialise the arrays used in the clipping algorithm for the VoronoiTessellation
vorn
.
Arguments
vorn
: TheVoronoiTessellation
.
Outputs
edges_to_process
: The set of edges that are to be processed.polygon_edge_queue
: The queue of edges that are to be processed.boundary_sites
: A mapping from boundary sites to the indices of the segment intersections that are incident to the boundary site.segment_intersections
: The list of segment intersections.processed_pairs
: The set of pairs of edges and polygons that have been processed.intersected_edge_cache
: The list of intersected edges currently being considered.exterior_circumcenters
: The list of circumcenters of sites that are outside the boundary.left_edge_intersectors
: The set of sites that intersect the edge to the left of an edge currently being considered.right_edge_intersectors
: The set of sites that intersect the edge to the right of an edge currently being considered.current_edge_intersectors
: The set of sites that intersect the current edge being considered.equal_circumcenter_mapping
: A mapping from the indices of the segment intersections that are equal to the circumcenter of a site to the index of the site.
DelaunayTriangulation.is_finite_segment
— Methodis_finite_segment(u, v) -> Bool
Returns true
if the segment (u, v)
is finite, and false
otherwise.
DelaunayTriangulation.is_ray_going_in
— Methodis_ray_going_in(u, v) -> Bool
Returns true
if the ray (u, v)
is coming in from infinity, and false
otherwise.
DelaunayTriangulation.is_ray_going_out
— Methodis_ray_going_out(u, v) -> Bool
Returns true
if the ray (u, v)
is going out to infinity, and false
otherwise.
DelaunayTriangulation.is_segment_between_two_ghosts
— Methodis_segment_between_two_ghosts(u, v) -> Bool
Returns true
if the segment (u, v)
is between two ghost vertices, and false
otherwise.
DelaunayTriangulation.process_intersection_points!
— Methodprocess_intersection_points!(polygon_edge_queue, vorn, current_incident_polygon,
left_edge_intersectors, right_edge_intersectors, current_edge_intersectors,
left_edge, right_edge, current_edge, processed_pairs, segment_intersections, boundary_sites)
Process the intersection points in left_edge_intersectors
, right_edge_intersectors
, and current_edge_intersectors
and add the new edges to polygon_edge_queue
if necessary. Special care is taken to not miss any corner points.
Arguments
polygon_edge_queue
: The queue of edges that need to be processed.vorn
: TheVoronoiTessellation
.current_incident_polygon
: The index of the current polygon being processed.left_edge_intersectors
: The intersection points ofleft_edge
with other edges.right_edge_intersectors
: The intersection points ofright_edge
with other edges.current_edge_intersectors
: The intersection points ofcurrent_edge
with other edges.left_edge
: The left edge of the current polygon.right_edge
: The right edge of the current polygon.current_edge
: The current edge of the current polygon.processed_pairs
: A set of pairs of edges and polygons that have already been processed.segment_intersections
: A dictionary of segment intersections.boundary_sites
: A dictionary of boundary sites.
Outputs
There are no outputs, but the caches and queues are updated in-place.
Extended help
The rules are based on the paper "Efficient Computation of Clipped Voronoi Diagram for Mesh Generation" by Yan, Wang, Levy, and Liu. Namely, an edge that intersects a boundary edge and one adjacent to it has its shared vertex added to the queue together with the current polygon (current_incident_polygon
) being considered, and any intersections have the adjacent polygon added to the queue together with the intersecting edge. (These are not strictly the rules in the paper.)
This function works as follows:
- First, assuming that there is more than one triangle in the underlying triangulation of
vorn
, we need to considerleft_edge
andright_edge
individually. - The procedure for each edge is the same, so here we just describe the
left_edge
. If there are any intersectors with theleft_edge
, and neither of(left_edge, current_incident_polygon)
or(reverse_edge(left_edge), current_incident_polygon)
have already been processed (i.e., inprocessed_pairs
), then we enqueue(left_edge, i)
and(left_edge, j)
intopolygon_edge_queue
, wherei
andj
are the vertices ofleft_edge
which correspond to polygons. This will ensure that we can find intersections next to this polygon. - After enqueueing these pairs, we also need to protect against corner points, which we check for by considering
current_incident_polygon ∈ all_indices
, whereall_indices
are the vertices ofleft_edge
,right_edge
, andcurrent_edge
. If this is true, and if the shared vertex ofcurrent_edge
andleft_edge
is equal tocurrent_incident_polygon
, then we need to add the point generator ofcurrent_incident_polygon
as an intersection. This need comes from having to worry about corners, i.e. points where the two unbounded polygons meet and go directly left and right of a vertex so that that vertex is not considered an intersection; this point needs to be included. - Once the
left_edge
andright_edge
have been processed as above, we need to then consider all ofleft_edge
,right_edge
, andcurrent_edge
, and each of the intersections through the respective edge. This step is done regardless of whether there is a single triangle in the underlying triangulation. The procedure for each edge is the same, so let us just describe thecurrent_edge
. For each edgeuv
in thecurrent_edge_intersectors
, we need to get the polygon adjacent to that edge. Then, if(current_edge, adjacent_incident_polygon)
or(reverse_edge(current_edge), adjacent_incident_polygon)
have not been processed, we enqueue(current_edge, adjacent_incident_polygon)
. - Once the edges have all been processed as above, we return.
DelaunayTriangulation.process_polygon!
— Functionprocess_polygon!(vorn::VoronoiTessellation, e, incident_polygon, boundary_sites, segment_intersections, intersected_edge_cache, exterior_circumcenters, equal_circumcenter_mapping, predicates::AbstractPredicateKernel=AdaptiveKernel()) -> (Edge, Edge, Edge)
Processes the polygon incident_polygon
for all of its intersections based on the boundary edge e
.
Arguments
vorn::VoronoiTessellation
: TheVoronoiTessellation
.e
: The edge on the boundary being considered.incident_polygon
: The index of the polygon being considered.boundary_sites
: The mapping from the indices of the sites on the boundary to the indices of the edges on the boundary that they intersect.segment_intersections
: The list of segment intersections.intersected_edge_cache
: A cache of the edges that have been intersected by the ray fromu
tov
.exterior_circumcenters
: A list of the circumcenters of the sites that are outside the convex hull of the sites on the boundary.equal_circumcenter_mapping
: A mapping from the indices of the segment intersections that are equal to the circumcenter of a site to the index of the site.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
left_edge
: The edge to the left ofe
on the boundary.right_edge
: The edge to the right ofe
on the boundary.e
: The edge on the boundary being considered.
In addition to these outputs, the caches are also updated in-place.
Extended help
This function works as follows:
- First, for the current edge
e
, we get the edgesleft_edge
andright_edge
that neighbour it viaget_neighbouring_boundary_edges
. - For each edge of the
incident_polygon
, we need to process it depending on whether the edge(u, v)
is finite, between two ghosts, going out to infinity, or coming in from infinity. If the edge is between two ghosts, we skip the edge. For rays that go out or in to infinity, we useprocess_ray_intersection!
andprocess_ray_intersection_with_other_edges!
to process the intersection of the ray with the boundary edges. The functionprocess_ray_intersection_with_other_edges!
is needed since rays going out to infinity may have to go through other boundary edges in order to do so, e.g. at a corner it may be that it crosses two boundary edges. For finite segments,process_segment_intersection!
is used to process the intersection. We apply this function with each ofe
,left_edge
, andright_edge
to check for all intersections. - The function is done once each of the polygon edges has been considered.
DelaunayTriangulation.process_ray_intersection!
— Functionprocess_ray_intersection!(
vorn::VoronoiTessellation,
u,
v,
incident_polygon,
intersected_edge_cache,
segment_intersections,
boundary_sites,
exterior_circumcenters,
equal_circumcenter_mapping,
predicates::AbstractPredicateKernel=AdaptiveKernel()) -> Point
Process the intersection of the Voronoi polygon of the site u
with the ray emanating from the circumcenter of the site v
.
Arguments
vorn
: TheVoronoiTessellation
.u
: The index of the siteu
, given as a ghost vertex for the associated ghost triangle.v
: The index of the sitev
.incident_polygon
: The index of the Voronoi polygon of the siteu
that is incident to the ray emanating from the circumcenter of the sitev
.intersected_edge_cache
: The list of intersected edges currently being considered.segment_intersections
: The list of segment intersections.boundary_sites
: A mapping from boundary sites to the indices of the segment intersections that are incident to the boundary site.exterior_circumcenters
: The list of circumcenters of sites that are outside the boundary.equal_circumcenter_mapping
: A mapping from the indices of the segment intersections that are equal to the circumcenter of a site to the index of the site.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
p
: The coordinates of the intersection.
In addition to the point p
, add_segment_intersection!
is also updated to incorporate the new intersection point, as is add_to_intersected_edge_cache!
.
DelaunayTriangulation.process_ray_intersection_with_other_edges!
— Functionprocess_ray_intersection_with_other_edges!(vorn::VoronoiTessellation,
u,
v,
e,
left_edge,
right_edge,
r,
segment_intersections,
boundary_sites,
incident_polygon,
equal_circumcenter_mapping,
intersected_edge_cache,
predicates::AbstractPredicateKernel=AdaptiveKernel())
Process the intersection of the ray from the ghost site u
to the site v
with the edges e
, left_edge
and right_edge
.
Arguments
vorn::VoronoiTessellation
: TheVoronoiTessellation
.u
: The index of the ghost site.v
: The index of the siteu
is going to.e
: The edge on the boundary being considered.left_edge
: The edge to the left ofe
on the boundary.right_edge
: The edge to the right ofe
on the boundary.r
: The coordinates of the intersection of the ray fromu
tov
with some edge. Ifany(isnan, r)
, then the ray does not intersect any edge and we skip.segment_intersections
: The list of segment intersections.boundary_sites
: The mapping from the indices of the sites on the boundary to the indices of the edges on the boundary that they intersect.incident_polygon
: The index of the polygon that contains the intersection of the ray fromu
tov
with the boundary.equal_circumcenter_mapping
: A mapping from the indices of the segment intersections that are equal to the circumcenter of a site to the index of the site.intersected_edge_cache
: A cache of the edges that have been intersected by the ray fromu
tov
.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 are no outputs, but add_segment_intersection!
and add_to_intersected_edge_cache!
are used to update the intersection objects.
DelaunayTriangulation.process_segment_intersection!
— Functionprocess_segment_intersection!(
vorn::VoronoiTessellation,
u,
v,
e,
incident_polygon,
intersected_edge_cache,
segment_intersections,
boundary_sites,
exterior_circumcenters,
equal_circumcenter_mapping,
predicates::AbstractPredicateKernel=AdaptiveKernel()) -> Point
Process the intersection of the Voronoi polygon's edge (u, v)
with the edge e
of the boundary, returning the coordinates of the intersection and updating via add_segment_intersection!
.
Arguments
vorn
: TheVoronoiTessellation
.u
: The index of the siteu
.v
: The index of the sitev
.e
: The edgee
of the boundary.incident_polygon
: The index of the Voronoi polygon currently being considered.intersected_edge_cache
: The list of intersected edges currently being considered.segment_intersections
: The list of segment intersections.boundary_sites
: A mapping from boundary sites to the indices of the segment intersections that are incident to the boundary site.exterior_circumcenters
: The list of circumcenters of sites that are outside the boundary.equal_circumcenter_mapping
: A mapping from the indices of the segment intersections that are equal to the circumcenter of a site to the index of the site.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
p
: The coordinates of the intersection. If there is no intersection, this is(NaN, NaN)
.
In addition to the point p
, add_segment_intersection!
is also updated to incorporate the new intersection point, as is add_to_intersected_edge_cache!
.
DelaunayTriangulation.single_triangle_clip!
— Methodsingle_triangle_clip!(vorn::VoronoiTessellation; predicates::AbstractPredicateKernel=AdaptiveKernel())
In the case that vorn
is dual to a triangulation with only a single triangle, this function clips the tessellation more efficiently than the general case with find_all_intersections
.
DelaunayTriangulation._get_ray
— Function_get_ray(vorn, i, ghost_vertex, predicates::AbstractPredicateKernel=AdaptiveKernel()) -> (Point, Point)
Extracts the ray from the i
th polygon of vorn
corresponding to the ghost_vertex
, where ghost_vertex
here means that get_polygon(vorn, i)[ghost_vertex]
is a ghost vertex.
Arguments
vorn
: TheVoronoiTessellation
.i
: The index of the polygon.ghost_vertex
: The index of the ghost vertex in the polygon.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
p
: The first point of the ray.q
: A second point of the ray, so thatpq
gives the direction of the ray (which extends to infinity).
DelaunayTriangulation.clip_bounded_polygon_to_bounding_box
— Methodclip_bounded_polygon_to_bounding_box(vorn::VoronoiTessellation, i, bounding_box; predicates::AbstractPredicateKernel=AdaptiveKernel()) -> Vector{NTuple{2,Number}}
Clips the i
th polygon of vorn
to bounding_box
.
See also clip_polygon
.
Arguments
vorn
: TheVoronoiTessellation
.i
: The index of the polygon.bounding_box
: The bounding box to clip the polygon to.
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
coords
: The coordinates of the clipped polygon. This is a circular vector.
DelaunayTriangulation.clip_unbounded_polygon_to_bounding_box
— Methodclip_unbounded_polygon_to_bounding_box(vorn::VoronoiTessellation, i, bounding_box; predicates::AbsractPredicateType=AdaptiveKernel()) -> Vector{NTuple{2,Number}}
Clips the i
th polygon of vorn
to bounding_box
. The polygon is assumed to be unbounded. See also clip_polygon
.
Use the keyword arguments predicates
to determine how predicates are computed. Should be one of ExactKernel
, AdaptiveKernel
, and FastKernel
. See the documentation for more information about these choices.
DelaunayTriangulation.get_bounded_polygon_coordinates
— Methodget_bounded_polygon_coordinates(vorn::VoronoiTessellation, i, bounding_box; predicates::AbstractPredicateKernel=AdaptiveKernel()) -> Vector{NTuple{2,Number}}
Returns the coordinates of the i
th polygon of vorn
, clipped to bounding_box
.
Use the keyword arguments predicates
to determine how predicates are computed. Should be one of ExactKernel
, AdaptiveKernel
, and FastKernel
. See the documentation for more information about these choices.
DelaunayTriangulation.get_clipping_poly_structs
— Methodget_clipping_poly_structs(vorn::VoronoiTessellation, i, bounding_box) -> (Polygon, Polygon)
Returns the polygons used for clipping the i
th polygon of vorn
to bounding_box
.
See also clip_polygon
.
Arguments
vorn
: TheVoronoiTessellation
.i
: The index of the polygon.bounding_box
: The bounding box to clip the polygon to.
Outputs
poly
: The polygon to clip.clip_poly
: The polygon to clip to.
DelaunayTriangulation.get_new_polygon_indices
— Methodget_new_polygon_indices(vorn, vertices) -> (Vector{Int}, Vector{NTuple{2,Float64}}, Tuple{Int, Int})
Returns the new vertices and points of the polygon, as well as the indices of the ghost vertices in the polygon.
Arguments
vorn
: TheVoronoiTessellation
.vertices
: The vertices of the polygon.
Outputs
new_vertices
: The new vertices of the polygon. This is not a circular vector. The vertices corresponding to a ghost vertex will be given by the ghost vertex itself.new_points
: The new points of the polygon. This is not a circular vector. The points corresponding to a ghost vertex will be given by by(NaN, NaN)
.ghost_vertices
: The indices of the ghost vertices innew_vertices
.
DelaunayTriangulation.get_polygon_coordinates
— Functionget_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
: TheVoronoiTessellation
.i
: The index of the polygon.bounding_box=nothing
: The bounding box to clip the polygon to. Ifnothing
, then the polygon is not clipped. If the polygon is unbounded, thenbounding_box
must be provided.
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
coords
: The coordinates of the polygon. This is a circular vector.
DelaunayTriangulation.get_unbounded_polygon_coordinates
— Methodget_unbounded_polygon_coordinates(vorn::VoronoiTessellation, i, bounding_box; predicates::AbstractPredicateKernel=AdaptiveKernel()) -> Vector{NTuple{2,Number}}
Returns the coordinates of the i
th polygon of vorn
, clipped to bounding_box
. The polygon is assumed to be unbounded.
Use the keyword arguments predicates
to determine how predicates are computed. Should be one of ExactKernel
, AdaptiveKernel
, and FastKernel
. See the documentation for more information about these choices.
DelaunayTriangulation.grow_polygon_outside_of_box
— Functiongrow_polygon_outside_of_box(vorn::VoronoiTessellation, i, bounding_box, predicates::AbstractPredicateKernel=AdaptiveKernel()) -> (Vector{Int}, Vector{NTuple{2,Number}})
Truncates the unbounded edges of the i
th polygon of vorn
so that the line connecting the truncated unbounded edges is entirely outside of bounding_box
.
Arguments
vorn
: TheVoronoiTessellation
.i
: The index of the polygon. The polygon must be unbounded.bounding_box
: The bounding box to clip the polygon to. See alsopolygon_bounds
.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
new_vertices
: The new vertices of the polygon. This is not a circular vector.new_points
: The new points of the polygon. This is not a circular vector.
DelaunayTriangulation.is_first_ghost_vertex
— Methodis_first_ghost_vertex(cell, i) -> Bool
Assuming that the circular vector cell
is such that ghost vertices only appear next to eachother in cell
and there are only two, tests if i
is the first ghost vertex in cell
.
See also is_last_ghost_vertex
.
Arguments
cell
: The circular vector.i
: The index of the vertex incell
.
Outputs
flag
:true
ifi
is the first ghost vertex incell
, andfalse
otherwise.
DelaunayTriangulation.is_last_ghost_vertex
— Methodis_last_ghost_vertex(cell, i) -> Bool
Assuming that the circular vector cell
is such that ghost vertices only appear next to eachother in cell
, tests if i
is the last ghost vertex in cell
.
See also is_first_ghost_vertex
.
Arguments
cell
: The circular vector.i
: The index of the vertex incell
.
Outputs
flag
:true
ifi
is the last ghost vertex incell
, andfalse
otherwise.
DelaunayTriangulation.maximum_distance_to_box
— Methodmaximum_distance_to_box(a, b, c, d, p) -> Number
Computes the maximum squared distance from the point p
to the box with corners (a, c)
, (b, c)
, (b, d)
, (a, d)
.
Arguments
a
: The minimumx
-coordinate of the box.b
: The maximumx
-coordinate of the box.c
: The minimumy
-coordinate of the box.d
: The maximumy
-coordinate of the box.p
: The point.
Outputs
dist
: The maximum squared distance fromp
to the box.
DelaunayTriangulation.liang_barsky
— Methodliang_barsky(a, b, c, d, p, q) -> (Point, Point)
Applies the Liang-Barsky algorithm to find the intersection of the line segment pq
with the rectangle [a, b] × [c, d]
.
Arguments
p
: The first point of the line segment.q
: The second point of the line segment.a
: The minimum x-coordinate of the rectangle.b
: The maximum x-coordinate of the rectangle.c
: The minimum y-coordinate of the rectangle.d
: The maximum y-coordinate of the rectangle.
Output
u
: The first coordinate of the intersection, or(NaN, NaN)
if there is no intersection.v
: The second coordinate of the intersection, or(NaN, NaN)
if there is no intersection.
DelaunayTriangulation.clip_polygon
— Methodclip_polygon(vertices, points, clip_vertices, clip_points; predicates::AbstractPredicateKernel=AdaptiveKernel()) -> Vector
Clip a polygon defined by (vertices, points)
to a convex clip polygon defined by (clip_vertices, clip_points)
with the Sutherland-Hodgman algorithm. The polygons should be defined in counter-clockwise order.
Arguments
vertices
: The vertices of the polygon to be clipped.points
: The underlying point set that the vertices are defined over.clip_vertices
: The vertices of the clipping polygon.clip_points
: The underlying point set that the clipping vertices are defined over.
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.
Output
clipped_polygon
: The coordinates of the clipped polygon, given in counter-clockwise order andclipped_polygon[begin] == clipped_polygon[end]
.
Centroidal Voronoi Tessellations
Here are some functions related to the computation of centroidal Voronoi tessellations.
DelaunayTriangulation._centroidal_smooth_itr
— Function_centroidal_smooth_itr(vorn::VoronoiTessellation, points, rng, predicates::AbstractPredicateKernel=AdaptiveKernel(); kwargs...) -> (VoronoiTessellation, Number)
Performs a single iteration of the centroidal smoothing algorithm.
Arguments
vorn
: TheVoronoiTessellation
.points
: The underlying point set. This is adeepcopy
of the points of the underlying triangulation.rng
: The random number generator.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.
Keyword Arguments
clip_polygon=nothing
: Ifclip=true
, then this is the polygon to clip the Voronoi tessellation to. Ifnothing
, the convex hull of the triangulation is used. The polygon should be defined as aTuple
of the form(points, boundary_nodes)
where theboundary_nodes
are vertices mapping to coordinates inpoints
, adhering to the usual conventions for defining boundaries.kwargs...
: Extra keyword arguments passed toretriangulate
.
Outputs
vorn
: The updatedVoronoiTessellation
.max_dist
: The maximum distance moved by any generator.
DelaunayTriangulation.centroidal_smooth
— Methodcentroidal_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
vorn
: TheVoronoiTessellation
.
Keyword Arguments
maxiters=1000
: The maximum number of iterations.clip_polygon=nothing
: Ifclip=true
, then this is the polygon to clip the Voronoi tessellation to. Ifnothing
, the convex hull of the triangulation is used. The polygon should be defined as aTuple
of the form(points, boundary_nodes)
where theboundary_nodes
are vertices mapping to coordinates inpoints
, adhering to the usual conventions for defining boundaries. Must be a convex polygon.tol=default_displacement_tolerance(vorn)
: The displacement tolerance. Seedefault_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 ofFastKernel
,ExactKernel
, andAdaptiveKernel
. See the documentation for a further discussion of these methods.kwargs...
: Extra keyword arguments passed toretriangulate
.
Outputs
vorn
: The updatedVoronoiTessellation
. This is not done in-place.
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.
DelaunayTriangulation.default_displacement_tolerance
— Methoddefault_displacement_tolerance(vorn::VoronoiTessellation) -> Number
Returns the default displacement tolerance for the centroidal smoothing algorithm. The default is given by max_extent / 1e4
, where max_extent = max(width, height)
, where width
and height
are the width and height of the bounding box of the underlying point set.
Arguments
vorn
: TheVoronoiTessellation
.
Outputs
tol
: The default displacement tolerance.
DelaunayTriangulation.move_generator_to_centroid!
— Methodmove_generator_to_centroid!(points, vorn::VoronoiTessellation, generator) -> Number
Moves the generator generator
to the centroid of its Voronoi cell. Returns the distance moved.
Arguments
points
: The underlying point set. This is adeepcopy
of the points of the underlying triangulation.vorn
: TheVoronoiTessellation
.generator
: The generator to move.
Outputs
dist
: The distance moved.
Triangulating Curve-Bounded Domains
We have several functions related to the triangulation of curve-bounded domains.
DelaunayTriangulation.coarse_discretisation!
— Methodcoarse_discretisation!(points, boundary_nodes, boundary_curve; n=0)
Constructs an initial coarse discretisation of a curve-bounded domain with bonudary defines by (points, boundary_nodes, boundary_curves)
, where boundary_nodes
and boundary_curves
should come from convert_boundary_curves!
. The argument n
is the amount of times to split an edge. If non-zero, this should be a power of two (otherwise it will be rounded up to the next power of two). If it is zero, then the splitting will continue until the maximum total variation over any subcurve is less than π/2.
DelaunayTriangulation.compute_split_position
— Functioncompute_split_position(enricher::BoundaryEnricher, i, j, predicates::AbstractPredicateKernel=AdaptiveKernel()) -> (Float64, Float64, NTuple{2,Float64})
Gets the point to split the edge (i, j)
at.
Arguments
enricher::BoundaryEnricher
: The enricher.i
: The first point of the edge.j
: The second point of the edge.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
t
: The parameter value of the split point.Δθ
: The total variation of the subcurve(i, t)
. If a split was created due to a small angle, this will be set to zero.ct
: The point to split the edge at.
DelaunayTriangulation.enrich_boundary!
— Methodenrich_boundary!(enricher::BoundaryEnricher; predicates::AbstractPredicateKernel=AdaptiveKernel())
Enriches the initial boundary defined inside enricher
, implementing the algorithm of Gosselin and Ollivier-Gooch (2007). At the termination of the algorithm, all edges will contain no other points inside their diametral circles.
The predicates
argument determines how predicates are computed, and should be one of ExactKernel
, FastKernel
, and AdaptiveKernel
(the default). See the documentation for more information about these choices.
DelaunayTriangulation.has_acute_neighbouring_angles
— Methodhas_acute_neighbouring_angles(predicates::AbstractPredicateKernel, enricher::BoundaryEnricher, i, j) -> Int, Vertex
Given a boundary edge (i, j)
, tests if the neighbouring angles are acute. The first returned value is the number of angles adjoining (i, j)
that are acute (0, 1, or 2). The second returned value is the vertex that adjoins the edge (i, j)
that is acute. If there are no such angles, or if there are two, then this returned vertex is 0
.
(The purpose of this function is similar to segment_vertices_adjoin_other_segments_at_acute_angle
.)
DelaunayTriangulation.split_subcurve!
— Functionsplit_subcurve!(enricher::BoundaryEnricher, i, j, predicates::AbstractPredicateKernel=AdaptiveKernel()) -> Bool
Splits the curve associated with the edge (i, j)
into two subcurves by inserting a point r
between (i, j)
such that the total variation of the subcurve is equal on (i, r)
and (r, j)
. The returned value is a flag
that is true
if there was a precision issue, and false
otherwise.
The predicate
argument determines how predicates are computed, and should be one of ExactKernel
, FastKernel
, and AdaptiveKernel
(the default). See the documentation for more information about these choices.
DelaunayTriangulation.test_visibility
— Methodtest_visibility(predicates::AbstractPredicateKernel, enricher::BoundaryEnricher, i, j, k) -> Certificate
Tests if the vertex k
is visible from the edge (i, j)
. Returns a Certificate
which is
Invisible
: Ifk
is not visible from(i, j)
.Visible
: Ifk
is visible from(i, j)
.
For this function, k
should be inside the diametral circle of (i, j)
.
We say that k
is invisibile from (i, j)
if the edges (i, k)
or (j, k)
intersect any other boundary edges, or there is a hole between (i, j)
and k
.
This is not the same definition used in defining constrained Delaunay triangulations, where visibility means visible from ANY point on the edge instead of only from the endpoints.
DelaunayTriangulation.triangulate_curve_bounded
— Methodtriangulate_curve_bounded(points::P;
segments=nothing,
boundary_nodes=nothing,
predicates::AbstractPredicateKernel=AdaptiveKernel(),
IntegerType::Type{I}=Int,
polygonise_n=4096,
coarse_n=0,
check_arguments=true,
delete_ghosts=false,
delete_empty_features=true,
recompute_representative_points=true,
rng::Random.AbstractRNG=Random.default_rng(),
insertion_order=nothing,
kwargs...) where {P,I} -> Triangulation
Triangulates a curve-bounded domain defined by (points, segments, boundary_nodes)
. Please see triangulate
for a description of the arguments. The only differences are:
insertion_order=nothing
: This argument is ignored for curve-bounded domains.polygonise_n=4096
: For generating a high-resolution discretisation of a boundary initially for the construction of aPolygonHierarchy
, many points are needed. This number of points is defined bypolygonise_n
, and must be a power of 2 (otherwise, the next highest power of 2 is used). Seepolygonise
.coarse_n=0
: This is the number of points to use for initialising a curve-bounded domain viacoarse_discretisation!
. The defaultcoarse_n=0
means the discretisation is performed until the maximum variation over any subcurve is less thanπ/2
.skip_points
: This is still used, but it is ignored during the enrichment phase (seeenrich_boundary!
).
See also BoundaryEnricher
and enrich_boundary!
.
To refine the mesh further beyond its initial coarse discretisation, as produced from this function, please see refine!
.
DelaunayTriangulation.convert_boundary_curves!
— Methodconvert_boundary_curves!(points, boundary_nodes, IntegerType) -> (NTuple{N, AbstractParametricCurve} where N, Vector)
Converts the provided points
and boundary_nodes
into a set of boundary curves and modified boundary nodes suitable for triangulation. In particular:
- The function gets
boundary_curves
fromto_boundary_curves
. boundary_nodes
is replaced with a set of initial boundary nodes (fromget_skeleton
). These nodes come from evaluating each boundary curve att = 0
andt = 1
. In the case of a piecewise linear boundary, the vertices are copied directly. Note that not all control points of aCatmullRomSpline
(whichis_interpolating
) will be added - only those att = 0
andt = 1
.- The
points
are modified to include the new boundary nodes. If a point is already inpoints
, it is not added again.
Arguments
points
: The point set. This is modified in place with the new boundary points.boundary_nodes
: The boundary nodes to be converted. This is not modified in place.IntegerType
: The type of integer to use for the boundary nodes.
Output
boundary_curves
: The boundary curves associated withboundary_nodes
.boundary_nodes
: The modified boundary nodes.
DelaunayTriangulation.is_curve_bounded
— Functionis_curve_bounded(tri::Triangulation) -> Bool
is_curve_bounded(boundary_nodes) -> Bool
Returns true
if tri
is curve bounded, and false
otherwise; similarly for the boundary_nodes
method.
DelaunayTriangulation.to_boundary_curves
— Methodto_boundary_curves(points, boundary_nodes) -> NTuple{N, AbstractParametricCurve} where N
Returns the set of boundary curves associated with boundary_nodes
and points
.
DelaunayTriangulation.is_visible
— Functionis_visible(cert::Certificate) -> Bool
Returns true
if cert
is Visible
, and false
otherwise.
DelaunayTriangulation.is_invisible
— Functionis_invisible(cert::Certificate) -> Bool
Returns true
if cert
is Invisible
, and false
otherwise.
DelaunayTriangulation.test_visibility
— Functiontest_visibility([kernel::AbstractPredicateKernel = AdaptiveKernel(),] tri::Triangulation, q, i) -> Certificate
Tests if the vertex i
and the point q
can see each other. Here, visibility means that the line segment joining the two does not intersect any segments.
Arguments
kernel::AbstractPredicateKernel = AdaptiveKernel()
: Method to use for computing predicates. Can be one ofFastKernel
,ExactKernel
(the default), andAdaptiveKernel
. See the documentation for a further discussion of these methods.tri
: TheTriangulation
.q
: The point from which we are testing visibility.i
: The vertex we are testing visibility of.
Outputs
cert
: ACertificate
. This will beVisible
ifi
is visible fromq
, andInvisible
otherwise.
test_visibility([kernel::AbstractPredicateKernel=AdaptiveKernel(),] tri::Triangulation, u, v, i; shift=0.0, attractor=get_point(tri,i)) -> Certificate
Tests if the edge (u, v)
and the point i
can see each other. Here, visibility means that any point in the interior of (u, v)
can see i
. To test this, we only check 10
points equally spaced between u
and v
, excluding u
and v
.
Arguments
kernel::AbstractPredicateKernel=AdaptiveKernel()
: How predicates are computed. See the documentation for information on the choices betweenFastKernel
,ExactKernel
, andAdaptiveKernel
.tri
: TheTriangulation
.u
: The first vertex of the edge.v
: The second vertex of the edge.i
: The vertex we are testing visibility of.
Keyword Arguments
shift=0.0
: The amount by which to shift each point on the edge towardsattractor
, i.e. ifp
is a point on the edge, thenp .+ shift .* (attractor - p)
is the point used to test visibility rather thanp
itself.attractor=get_point(tri,i)
: Related toshift
; see above.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
cert
: ACertificate
. This will beVisible
ifi
is visible from(u, v)
, andInvisible
otherwise.
test_visibility(predicates::AbstractPredicateKernel, enricher::BoundaryEnricher, i, j, k) -> Certificate
Tests if the vertex k
is visible from the edge (i, j)
. Returns a Certificate
which is
Invisible
: Ifk
is not visible from(i, j)
.Visible
: Ifk
is visible from(i, j)
.
For this function, k
should be inside the diametral circle of (i, j)
.
We say that k
is invisibile from (i, j)
if the edges (i, k)
or (j, k)
intersect any other boundary edges, or there is a hole between (i, j)
and k
.
This is not the same definition used in defining constrained Delaunay triangulations, where visibility means visible from ANY point on the edge instead of only from the endpoints.