Skip duplicates by default (#218)

.gitignore

@@ -9,3 +9,5 @@ docs/Manifest.toml
 /.vscode
 LocalPreferences.toml
 docs/src/applications/cell_simulation.mp4
+varying_weight_power.mp4
+varying_weight.mp4

NEWS.md

@@ -1,5 +1,9 @@
 # Changelog
 
+## 1.6.4 
+
+- An error is no longer thrown for inputs with duplicate points. Instead, a warning is thrown and any duplicates are merged into the `skip_points` keyword argument. With this, `DuplicatePointsError` has been removed. To silence the new warning, use `DelaunayTriangulation.toggle_warn_on_dupes!()`.
+
 ## 1.6.1
 - Fix issue with clipping Voronoi tessellation dual to a single right-angled triangle. See [#207](https://github.com/JuliaGeometry/DelaunayTriangulation.jl/pull/207)
 

Project.toml

@@ -1,7 +1,7 @@
 name = "DelaunayTriangulation"
 uuid = "927a84f5-c5f4-47a5-9785-b46e178433df"
 authors = ["Daniel VandenHeuvel <[email protected]>"]
-version = "1.6.3"
+version = "1.6.4"
 
 [deps]
 AdaptivePredicates = "35492f91-a3bd-45ad-95db-fcad7dcfedb7"

README.md

@@ -124,7 +124,7 @@ ax = Axis(fig[2, 3]; title = "Curve-Bounded", wh...);            triplot!(ax, tr
 ax = Axis(fig[2, 4]; title = "Disjoint Curve-Bounded", wh...);   triplot!(ax, tri8)
 ax = Axis(fig[3, 1]; title = "Weighted", wh...);                 triplot!(ax, tri9)
 ax = Axis(fig[3, 2]; title = "Power Diagram", wh...);            voronoiplot!(ax, vorn10)
-ax = Axis(fig[3, 3]; title = "Generic Clipped Voronoi", wh...);  voronoiplot!(ax, vorn11)
+ax = Axis(fig[3, 3]; title = "Clipped Voronoi", wh...);          voronoiplot!(ax, vorn11, color=:white, strokewidth = 4)
 ```
 
 ![](readme.png)

docs/src/api/triangulation.md

@@ -12,6 +12,7 @@ refine!
 retriangulate 
 convert_boundary_points_to_indices
 check_args
+toggle_warn_on_dupes!
 get_points
 get_triangles
 get_boundary_nodes

docs/src/literate_tutorials/clipped_polygon.jl

@@ -36,11 +36,11 @@ clipped_vorn = voronoi(tri, clip = true, clip_polygon = clip_polygon)
 fig = Figure()
 ax1 = Axis(fig[1, 1], title = "Unclipped", width = 600, height = 400)
 ax2 = Axis(fig[1, 2], title = "Clipped", width = 600, height = 400)
-voronoiplot!(ax1, vorn, show_generators = false, colormap = :matter, strokewidth = 4)
+voronoiplot!(ax1, vorn, show_generators = false, color = :white, strokewidth = 4)
 xlims!(ax1, -2, 2)
 ylims!(ax1, -2, 2)
 lines!(ax1, [clip_points..., clip_points[begin]], color = :black, linewidth = 4, linestyle = :dash)
-voronoiplot!(ax2, clipped_vorn, show_generators = false, colormap = :matter, strokewidth = 4)
+voronoiplot!(ax2, clipped_vorn, show_generators = false, color = :white, strokewidth = 4)
 xlims!(ax2, -2, 2)
 ylims!(ax2, -2, 2)
 resize_to_layout!(fig)
@@ -61,11 +61,11 @@ clipped_vorn = voronoi(tri, clip = true, clip_polygon = clip_polygon)
 fig = Figure()
 ax1 = Axis(fig[1, 1], title = "Unclipped", width = 600, height = 400)
 ax2 = Axis(fig[1, 2], title = "Clipped", width = 600, height = 400)
-voronoiplot!(ax1, vorn, show_generators = false, colormap = :matter, strokewidth = 4)
+voronoiplot!(ax1, vorn, show_generators = false, color = :white, strokewidth = 4)
 xlims!(ax1, -2, 2)
 ylims!(ax1, -2, 2)
 lines!(ax1, [clip_points..., clip_points[begin]], color = :black, linewidth = 4, linestyle = :dash)
-voronoiplot!(ax2, clipped_vorn, show_generators = false, colormap = :matter, strokewidth = 4)
+voronoiplot!(ax2, clipped_vorn, show_generators = false, color = :white, strokewidth = 4)
 xlims!(ax2, -2, 2)
 ylims!(ax2, -2, 2)
 resize_to_layout!(fig)

src/algorithms/triangulation/check_args.jl

@@ -1,13 +1,15 @@
 """
-    check_args(points, boundary_nodes, hierarchy::PolygonHierarchy, boundary_curves = ()) -> Bool 
+    check_args(points, boundary_nodes, hierarchy::PolygonHierarchy, boundary_curves = (); skip_points = Set{Int}()) -> Bool 
 
 Check that the arguments `points` and `boundary_nodes` to [`triangulate`](@ref), and a constructed 
 [`PolygonHierarchy`](@ref) given by `hierarchy`, are valid. In particular, the function checks:
 
 - The dimension of the points. If the points are not 2D, a warning is thrown.
-- The points are all unique. If they are not, a `DuplicatePointsError` is thrown.
+- The points are all unique. If they are not, a warning is thrown and the indices of the duplicates are merged into `skip_points`.
 - There are at least three points. If there are not, an `InsufficientPointsError` is thrown.
 
+If any duplicate points are found, the indices of the duplicates are merged into `skip_points` in-place.
+
 If `boundary_nodes` are provided, meaning [`has_boundary_nodes`](@ref), then the function also checks:
 
 - If the boundary curves all connect consistently. Here, this means that each section of a boundary curve ends at the start of the next boundary section;
@@ -21,9 +23,9 @@ If `boundary_nodes` are provided, meaning [`has_boundary_nodes`](@ref), then the
     Another requirement for [`triangulate`](@ref) is that none of the boundaries intersect in their interior, which also prohibits 
     interior self-intersections. This is NOT checked. Similarly, segments should not intersect in their interior, which is not checked.
 """
-function check_args(points, boundary_nodes, hierarchy, boundary_curves = ())
+function check_args(points, boundary_nodes, hierarchy, boundary_curves = (); skip_points = Set{Int}())
     check_dimension(points)
-    has_unique_points(points)
+    has_unique_points!(skip_points, points)
     has_enough_points(points)
     has_bnd = has_boundary_nodes(boundary_nodes)
     if has_bnd
@@ -33,9 +35,6 @@ function check_args(points, boundary_nodes, hierarchy, boundary_curves = ())
     return true
 end
 
-struct DuplicatePointsError{P} <: Exception
-    points::P
-end
 struct InsufficientPointsError{P} <: Exception
     points::P
 end
@@ -52,22 +51,6 @@ struct InconsistentOrientationError{I} <: Exception
     is_sectioned::Bool
     is_curve_bounded::Bool
 end
-function Base.showerror(io::IO, err::DuplicatePointsError)
-    points = err.points
-    dup_seen = find_duplicate_points(points)
-    println(io, "DuplicatePointsError: There were duplicate points. The following points are duplicated:")
-    n = length(dup_seen)
-    ctr = 1
-    for (p, ivec) in dup_seen
-        if ctr < n
-            println(io, "    ", p, " at indices ", ivec)
-        else
-            print(io, "    ", p, " at indices ", ivec, ".")
-        end
-        ctr += 1
-    end
-    return io
-end
 function Base.showerror(io::IO, err::InsufficientPointsError)
     points = err.points
     print(io, "InsufficientPointsError: The provided point set has ", num_points(points), " points, but triangulations require at least three points.")
@@ -116,9 +99,29 @@ function check_dimension(points)
     return valid
 end
 
-function has_unique_points(points)
+function has_unique_points!(skip_points, points)
     all_unique = points_are_unique(points)
-    !all_unique && throw(DuplicatePointsError(points))
+    if !all_unique 
+        dup_seen = find_duplicate_points(points)
+        if WARN_ON_DUPES[]
+            io = IOBuffer()
+            println(io, "There were duplicate points. Only one of each duplicate will be used, and all other duplicates will be skipped. The indices of the duplicates are:")
+        end
+        for (p, ivec) in dup_seen 
+            for j in 2:lastindex(ivec)
+                push!(skip_points, ivec[j])
+            end
+            if WARN_ON_DUPES[]
+                println(io, "  ", p, " at indices ", ivec)
+            end
+        end
+        if WARN_ON_DUPES[]
+            println(io, "To suppress this warning, call `DelaunayTriangulation.toggle_warn_on_dupes!()`.")
+        end
+        if WARN_ON_DUPES[]
+            @warn String(take!(io))
+        end
+    end
     return true
 end
 

src/algorithms/triangulation/main.jl

@@ -151,10 +151,15 @@ function triangulate(
     if isnothing(full_polygon_hierarchy)
         full_polygon_hierarchy = construct_polygon_hierarchy(points, boundary_nodes; IntegerType)
     end
-    check_arguments && check_args(points, boundary_nodes, full_polygon_hierarchy, boundary_curves)
+    skip_points_set = Set{IntegerType}(skip_points)
+    n = length(skip_points_set)
+    check_arguments && check_args(points, boundary_nodes, full_polygon_hierarchy, boundary_curves; skip_points = skip_points_set)
+    if length(skip_points_set) > n 
+        setdiff!(insertion_order, skip_points_set) 
+    end
     tri = Triangulation(points; IntegerType, EdgeType, TriangleType, EdgesType, TrianglesType, weights, boundary_curves, boundary_enricher, build_cache = Val(true))
     return _triangulate!(
-        tri, segments, boundary_nodes, predicates, randomise, try_last_inserted_point, skip_points, num_sample_rule, rng, insertion_order,
+        tri, segments, boundary_nodes, predicates, randomise, try_last_inserted_point, skip_points_set, num_sample_rule, rng, insertion_order,
         recompute_representative_points, delete_holes, full_polygon_hierarchy, delete_ghosts, delete_empty_features,
     )
 end

src/public.jl

@@ -52,6 +52,7 @@
                 num_generators,
                 polygon_features,
                 toggle_inf_warn!,
+                toggle_warn_on_dupes!,
                 AbstractParametricCurve,
                 angle_between,
                 arc_length,

src/setup.jl

@@ -48,6 +48,15 @@ By default, this warning is enabled.
 """
 toggle_inf_warn!() = (INF_WARN[] = !INF_WARN[])
 
+const WARN_ON_DUPES = Ref(true)
+"""
+    toggle_warn_on_dupes!()
+
+Toggle the warning for duplicate points in the input data.
+By default, this warning is enabled.
+"""
+toggle_warn_on_dupes!() = (WARN_ON_DUPES[] = !WARN_ON_DUPES[])
+
 @eval macro $(Symbol("const"))(field)
     if VERSION >= v"1.8.0-DEV.1148"
         return Expr(:const, esc(field))

test/Project.toml

@@ -6,6 +6,7 @@ CairoMakie = "13f3f980-e62b-5c42-98c6-ff1f3baf88f0"
 ColorSchemes = "35d6a980-a343-548e-a6ea-1d62b119f2f4"
 DataStructures = "864edb3b-99cc-5e75-8d2d-829cb0a9cfe8"
 Dates = "ade2ca70-3891-5945-98fb-dc099432e06a"
+DelaunayTriangulation = "927a84f5-c5f4-47a5-9785-b46e178433df"
 DelimitedFiles = "8bb1440f-4735-579b-a4ab-409b98df4dab"
 Distributions = "31c24e10-a181-5473-b8eb-7969acd0382f"
 Documenter = "e30172f5-a6a5-5a46-863b-614d45cd2de4"

test/makie/makie.jl

@@ -302,7 +302,6 @@ end
         f
     end
 
-
     @test_reference "Voronoiplot with some custom bounding boxes may not contain all data sites.png" begin
         reseed!()
         points = [(-3.0, 7.0), (1.0, 6.0), (-1.0, 3.0), (-2.0, 4.0), (3.0, -2.0), (5.0, 5.0), (-4.0, -3.0), (3.0, 8.0)]

test/readme_example.jl

@@ -93,7 +93,7 @@ ax = Axis(fig[2, 3]; title = "Curve-Bounded", wh...);            triplot!(ax, tr
 ax = Axis(fig[2, 4]; title = "Disjoint Curve-Bounded", wh...);   triplot!(ax, tri8)
 ax = Axis(fig[3, 1]; title = "Weighted", wh...);                 triplot!(ax, tri9)
 ax = Axis(fig[3, 2]; title = "Power Diagram", wh...);            voronoiplot!(ax, vorn10)
-ax = Axis(fig[3, 3]; title = "Clipped Voronoi", wh...);          voronoiplot!(ax, vorn11)
+ax = Axis(fig[3, 3]; title = "Clipped Voronoi", wh...);          voronoiplot!(ax, vorn11, color=:white, strokewidth = 4)
 
 readme_img = joinpath(dirname(dirname(pathof(DelaunayTriangulation))), "readme.png")
 @test_reference readme_img fig by = psnr_equality(10)