Fix curve-bounded refinement with custom structs (#175)

Author: DanielVandH

.zenodo.json

@@ -19,4 +19,4 @@
 	"title": "DelaunayTriangulation.jl",
 	"upload_type": "software",
 	"language": "eng"
-}
+}

NEWS.md

@@ -1,8 +1,12 @@
 # Changelog
 
+## v1.1.4
+
+- Fixed a bug with curve-bounded refinement with custom edge(s) structs. See [#175](https://github.com/JuliaGeometry/DelaunayTriangulation.jl/pull/175).
+
 ## v1.1.3
 
-- `sort_triangle` is now public. 
+- `sort_triangle` is now public. See [#174](https://github.com/JuliaGeometry/DelaunayTriangulation.jl/pull/174).
 
 ## v1.1.2
 

Project.toml

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

src/algorithms/triangulation/triangulate_curve_bounded.jl

@@ -39,6 +39,8 @@ function triangulate_curve_bounded(
         boundary_nodes = nothing,
         predicates::AbstractPredicateKernel = AdaptiveKernel(),
         IntegerType::Type{I} = Int,
+        EdgeType::Type{E} = isnothing(segments) ? NTuple{2, IntegerType} : (edge_type ∘ typeof)(segments),
+        EdgesType::Type{Es} = isnothing(segments) ? Set{EdgeType} : typeof(segments),
         polygonise_n = 4096,
         coarse_n = 0,
         check_arguments = true,
@@ -48,12 +50,14 @@ function triangulate_curve_bounded(
         rng::Random.AbstractRNG = Random.default_rng(),
         insertion_order = nothing, # use this so that it gets ignored by the kwargs
         kwargs...,
-    ) where {P, I}
-    enricher = BoundaryEnricher(points, boundary_nodes, segments; IntegerType, n = polygonise_n, coarse_n)
+    ) where {P, I, E, Es}
+    enricher = BoundaryEnricher(points, boundary_nodes, segments; IntegerType, EdgeType, EdgesType, n = polygonise_n, coarse_n)
     return _triangulate_curve_bounded(
         points, enricher;
         predicates,
         IntegerType,
+        EdgeType,
+        EdgesType,
         check_arguments,
         delete_ghosts,
         delete_empty_features,
@@ -67,14 +71,16 @@ function _triangulate_curve_bounded(
         points::P, enricher;
         predicates::AbstractPredicateKernel = AdaptiveKernel(),
         IntegerType::Type{I} = Int,
+        EdgeType::Type{E},
+        EdgesType::Type{Es},
         check_arguments = true,
         delete_ghosts = false,
         delete_empty_features = true,
         recompute_representative_points = true,
         rng::Random.AbstractRNG = Random.default_rng(),
         insertion_order = nothing, # use this so that it gets ignored by the kwargs
         kwargs...,
-    ) where {P, I}
+    ) where {P, I, E, Es}
     check_arguments && check_args(enricher)
     enrich_boundary!(enricher; predicates)
     new_boundary_nodes = get_boundary_nodes(enricher)
@@ -84,6 +90,8 @@ function _triangulate_curve_bounded(
     tri = triangulate(
         points;
         IntegerType,
+        EdgeType,
+        EdgesType,
         segments = new_segments,
         boundary_nodes = new_boundary_nodes,
         predicates,

src/data_structures/mesh_refinement/boundary_enricher.jl

@@ -337,7 +337,7 @@ function get_minimum_edge_length(complex::SmallAngleComplex, points)
 end
 
 """
-    BoundaryEnricher{P,B,C,I,BM,S}
+    BoundaryEnricher{P,B,C,I,BM,S,E}
 
 Struct used for performing boundary enrichment on a curve-bounded boundary. 
 
@@ -367,34 +367,39 @@ field will no longer aliased with the input `boundary_nodes`, although `points`
 The argument `n` is used in [`polygonise`](@ref) for filling out the boundary temporarily in order to construct the [`PolygonHierarchy`](@ref). The argument `coarse_n` defines the initial coarse discretisation 
 through [`coarse_discretisation!`](@ref); the default `n=0` means that the coarse discretisation will be performed until the maximum total variation of a subcurve is less than π/2.
 """
-struct BoundaryEnricher{P, B, C, I, BM, S}
+struct BoundaryEnricher{P, B, C, I, BM, S, E}
     points::P
     boundary_nodes::B
     segments::S
     boundary_curves::C
     polygon_hierarchy::PolygonHierarchy{I}
-    parent_map::Dict{NTuple{2, I}, I}
+    parent_map::Dict{E, I}
     curve_index_map::Dict{I, I}
     boundary_edge_map::BM
     spatial_tree::BoundaryRTree{P}
     queue::Queue{I}
     small_angle_complexes::Dict{I, Vector{SmallAngleComplex{I}}}
 end
-function BoundaryEnricher(points::P, boundary_nodes::B, segments = nothing; IntegerType = Int, n = 4096, coarse_n = 0) where {P, B}
-    boundary_curves, new_boundary_nodes = convert_boundary_curves!(points, boundary_nodes, IntegerType)
+function BoundaryEnricher(
+        points::P, boundary_nodes::B, segments = nothing; n = 4096, coarse_n = 0,
+        IntegerType::Type{I} = Int,
+        EdgeType::Type{E} = isnothing(segments) ? NTuple{2, IntegerType} : (edge_type ∘ typeof)(segments),
+        EdgesType::Type{Es} = isnothing(segments) ? Set{EdgeType} : typeof(segments),
+    ) where {P, B, I, E, Es}
+    boundary_curves, new_boundary_nodes = convert_boundary_curves!(points, boundary_nodes, I)
     polygon_hierarchy = construct_polygon_hierarchy(points, new_boundary_nodes, boundary_curves; IntegerType, n)
-    return _construct_boundary_enricher(points, new_boundary_nodes, boundary_curves, polygon_hierarchy, segments, n, coarse_n, IntegerType)
+    return _construct_boundary_enricher(points, new_boundary_nodes, boundary_curves, polygon_hierarchy, segments, n, coarse_n, I, E, Es)
 end
-function _construct_boundary_enricher(points, boundary_nodes, boundary_curves, polygon_hierarchy, segments, n, coarse_n, ::Type{I}) where {I}
+function _construct_boundary_enricher(points, boundary_nodes, boundary_curves, polygon_hierarchy, segments, n, coarse_n, ::Type{I}, ::Type{E}, ::Type{Es}) where {I, E, Es}
     expand_bounds!(polygon_hierarchy, ε(number_type(points)))
     coarse_discretisation!(points, boundary_nodes, boundary_curves; n = coarse_n)
     boundary_edge_map = construct_boundary_edge_map(boundary_nodes, I)
-    parent_map = Dict{NTuple{2, I}, I}()
+    parent_map = Dict{E, I}()
     curve_index_map = Dict{I, I}()
     spatial_tree = BoundaryRTree(points)
     queue = Queue{I}()
     small_angle_complexes = get_small_angle_complexes(points, boundary_nodes, boundary_curves, segments; IntegerType = I)
-    _segments = isnothing(segments) ? Set{NTuple{2, I}}() : segments
+    _segments = isnothing(segments) ? Es() : segments
     enricher = BoundaryEnricher(points, boundary_nodes, _segments, boundary_curves, polygon_hierarchy, parent_map, curve_index_map, boundary_edge_map, spatial_tree, queue, small_angle_complexes)
     construct_parent_map!(enricher)
     construct_curve_index_map!(enricher)

test/operations/delete_triangle.jl

@@ -432,9 +432,9 @@ end
         DT.add_triangle!(tri, i, j, k; update_ghost_edges = true)
     end
     [
-        DT.delete_triangle!(tri, i, j, k; update_ghost_edges = true) for (i, j, k) in (
-            (1, 8, 9), (9, 8, 10), (10, 11, 5), (5, 7, 4), (4, 6, 3), (3, 6, 2), (2, 6, 1),
-        )
+    DT.delete_triangle!(tri, i, j, k; update_ghost_edges = true) for (i, j, k) in (
+        (1, 8, 9), (9, 8, 10), (10, 11, 5), (5, 7, 4), (4, 6, 3), (3, 6, 2), (2, 6, 1),
+    )
     ]
     true_T = Set{NTuple{3, Int}}(
         [

test/refinement/curve_bounded.jl

@@ -1636,3 +1636,50 @@ end
         @test tri.boundary_enricher.segments ∈ (Set(((r + 1, r + 2),)), Set(((r + 2, r + 1),)))
     end
 end
+
+@testset "Custom structs" begin
+    # Used to be a broken example
+    DT = DelaunayTriangulation
+    struct Custom2Points
+        points::Vector{NTuple{2, Float64}}
+    end
+    struct Custom2Segments
+        segments::Set{NTuple{2, Int}}
+    end
+    struct Custom2Triangles
+        triangles::Set{NTuple{3, Int}}
+    end
+    Base.eachindex(points::Custom2Points) = Base.eachindex(points.points)
+    Base.iterate(points::Custom2Points, state...) = Base.iterate(points.points, state...)
+    Base.length(points::Custom2Points) = length(points.points)
+    Base.getindex(points::Custom2Points, i) = points.points[i]
+    DT.number_type(::Type{Custom2Points}) = Float64
+    Base.iterate(triangles::Custom2Triangles, state...) = Base.iterate(triangles.triangles, state...)
+    Base.sizehint!(triangles::Custom2Triangles, n) = sizehint!(triangles.triangles, n)
+    Base.eltype(::Type{Custom2Triangles}) = NTuple{3, Int}
+    Base.length(triangles::Custom2Triangles) = length(triangles.triangles)
+    Custom2Triangles() = Custom2Triangles(Set{NTuple{3, Int}}())
+    Base.iterate(segments::Custom2Segments, state...) = Base.iterate(segments.segments, state...)
+    Base.eltype(::Type{Custom2Segments}) = NTuple{2, Int}
+    Base.rand(rng::AbstractRNG, segments::Custom2Segments) = rand(rng, segments.segments)
+    Base.length(segments::Custom2Segments) = length(segments.segments)
+    Custom2Segments() = Custom2Segments(Set{NTuple{2, Int}}())
+    DT.contains_edge(e::NTuple{2, Int}, Es::Custom2Segments) = e ∈ Es.segments
+    Base.empty!(triangles::Custom2Triangles) = empty!(triangles.triangles)
+    Base.empty!(segments::Custom2Segments) = empty!(segments.segments)
+    Base.push!(segments::Custom2Segments, e::NTuple{2, Int}) = push!(segments.segments, e)
+    Base.delete!(segments::Custom2Segments, e::NTuple{2, Int}) = delete!(segments.segments, e)
+    Base.pop!(points::Custom2Points) = pop!(points.points)
+    DT.push_point!(points::Custom2Points, x, y) = push!(points.points, (x, y))
+    Base.delete!(triangles::Custom2Triangles, T::NTuple{3, Int}) = delete!(triangles.triangles, T)
+    DT.set_point!(points::Custom2Segments, i, x, y) = points.points[i] = (x, y)
+    Base.push!(triangles::Custom2Triangles, T::NTuple{3, Int}) = push!(triangles.triangles, T)
+    points = [(0.0, 0.0), (1.0, 0.0), (1.0, 1.0), (0.0, 1.0)]
+    append!(points, [(clamp(0.5 + randn(), -2, 2), clamp(0.5 + randn(), -2, 2)) for _ in 1:100])
+    unique!(points)
+    inner_circle = CircularArc((0.5, 0.25), (0.5, 0.25), (0.5, 0.5), positive = false)
+    boundary_nodes = [[[1, 2, 3, 4, 1]], [[inner_circle]]]
+    points = Custom2Points(points)
+    tri = triangulate(points; boundary_nodes, TrianglesType = Custom2Triangles, EdgesType = Custom2Segments)
+    @test DT.validate_triangulation(tri)
+end