Point Location
DelaunayTriangulation.brute_force_search
— Functionbrute_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.find_triangle
— Functionfind_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.get_nearest_neighbour
— Functionget_nearest_neighbour(tri_or_vor, q; kwargs...)
Get the index of the nearest neighbour of q
in tri_or_vor
.
For power diagrams, distance is measured using get_power_distance
(with q
being assigned zero weight).
Arguments
tri_or_vor
: ATriangulation
orVoronoiTessellation
.q
: The point to be located.
Keyword Arguments
kwargs...
: Keyword arguments passed tofind_triangle
.
Output
i
: The index of the nearest neighbour. This is a point of the triangulation iftri_or_vor
is aTriangulation
or of a generator iftri_or_vor
is aVoronoiTessellation
.
DelaunayTriangulation.find_polygon
— Functionfind_polygon(tri::Triangulation, q) -> Integer
Given a point q
, finds the index of the polygon in the triangulation tri
that contains q
. If q
is on the boundary of the triangulation or outside the triangulation, the function returns 0
.
See also dist
and distance_to_polygon
.