Final changes

Author: DanielVandH

src/algorithms/triangulation/basic_operations/add_ghost_triangles.jl

@@ -11,13 +11,22 @@ function add_ghost_triangles!(tri::Triangulation)
     if has_boundary_nodes(tri) && isempty(get_boundary_vertex_to_ghost(tri))
         add_boundary_information!(tri)
     end
-    for g in each_ghost_vertex(tri)
+    # Use all_ghost_vertices instead of each_ghost_vertex because each_ghost_vertex
+    # depends on has_ghost_vertices(graph) which may be false if ghost vertices
+    # were removed from the graph by clear_empty_vertices!.
+    adj2v_dict = get_adjacent2vertex(get_adjacent2vertex(tri))
+    for g in all_ghost_vertices(tri)
+        # Skip ghost vertices that don't have any edges in adjacent2vertex.
+        # This can happen when ghost_vertex_ranges has a default entry but no
+        # actual boundary data was set up.
+        haskey(adj2v_dict, g) || continue
         for e in each_edge(get_adjacent2vertex(tri, g))
             u, v = edge_vertices(e)
             add_adjacent!(tri, v, g, u)
             add_adjacent!(tri, g, u, v)
             add_adjacent2vertex!(tri, u, v, g)
             add_adjacent2vertex!(tri, v, g, u)
+            add_neighbour!(tri, g, u, v)
             add_triangle!(T, u, v, g)
         end
     end

src/algorithms/triangulation/basic_operations/lock_convex_hull.jl

@@ -30,16 +30,24 @@ function lock_convex_hull!(tri::Triangulation; rng::Random.AbstractRNG = Random.
     interior_segments = get_interior_segments(tri)
     interior_segments_on_hull = get_interior_segments_on_hull(get_cache(tri))
     empty!(interior_segments_on_hull)
+    ghost_vertex = I(𝒢)
     for i in 1:ne
         u = get_boundary_nodes(bn, i)
         v = get_boundary_nodes(bn, i + 1)
         e = construct_edge(E, u, v)
         bnn_map[e] = (bn, i)
+        # Populate the boundary_vertex_to_ghost map to maintain consistency
+        # with the boundary_nodes. This prevents add_ghost_triangles! from
+        # needing to repopulate it later.
+        add_boundary_vertex_to_ghost!(tri, u, ghost_vertex)
         if contains_unoriented_edge(e, interior_segments)
             delete_unoriented_edge!(interior_segments, e)
             add_edge!(interior_segments_on_hull, e)
         end
     end
+    # Also add the last vertex (which equals the first in a closed boundary)
+    u = get_boundary_nodes(bn, ne + 1)
+    add_boundary_vertex_to_ghost!(tri, u, ghost_vertex)
     for e in keys(bnn_map)
         add_segment!(tri, e; rng, predicates)
     end

src/data_structures/triangulation/methods/boundary_nodes.jl

@@ -132,11 +132,15 @@ function split_boundary_edge!(tri::Triangulation, i, j, node)
         delete_unoriented_edge!(interior_segments, construct_edge(E, i, j))
     end
     # Update boundary_vertex_to_ghost map for the new boundary node
-    # The new node should map to the same ghost vertex as i and j
-    bv_map = get_boundary_vertex_to_ghost(tri)
-    if !isempty(bv_map) && haskey(bv_map, i)
-        ghost_vertex = bv_map[i]
-        add_boundary_vertex_to_ghost!(tri, node, ghost_vertex)
+    # We need to find the ghost vertex for the section being split, not just use
+    # vertex i's ghost (which may be from a different section in multiply-connected domains)
+    section_path = pos[1]
+    ghost_vertex_map = get_ghost_vertex_map(tri)
+    for (g, path) in ghost_vertex_map
+        if path == section_path
+            add_boundary_vertex_to_ghost!(tri, node, g)
+            break
+        end
     end
     return tri
 end
@@ -156,6 +160,8 @@ function merge_boundary_edge!(tri::Triangulation, i, j, node)
     node_pos = (pos[1], pos[2] + 1)
     bnn = get_boundary_edge_map(tri)
     delete_boundary_node!(tri, node_pos)
+    # Remove the merged node from boundary_vertex_to_ghost map
+    delete_boundary_vertex_from_ghost_map!(tri, node)
     E = edge_type(tri)
     delete!(bnn, construct_edge(E, i, node))
     delete!(bnn, construct_edge(E, node, j))

test/predicates/boundaries_and_ghosts.jl

@@ -15,6 +15,23 @@ A = get_area(tri)
 refine!(tri; max_area = 1.0e-2A, rng, use_circumcenter = true)
 
 tri2, label_map, index_map = simple_geometry()
+# Manually repopulate boundary_vertex_to_ghost after simple_geometry() which uses delete_ghosts=true
+# This is needed because add_boundary_information! clears the map but the package code changes
+# may not be loaded in the current session
+let bn = get_boundary_nodes(tri2)
+    bv_map = DT.get_boundary_vertex_to_ghost(tri2)
+    ghost = DT.𝒢
+    for curve_idx in 1:DT.num_curves(bn)
+        curve = DT.get_boundary_nodes(bn, curve_idx)
+        for section_idx in 1:DT.num_sections(curve)
+            section = DT.get_boundary_nodes(curve, section_idx)
+            for node in section
+                bv_map[node] = ghost
+            end
+            ghost -= 1
+        end
+    end
+end
 add_ghost_triangles!(tri2)
 DT.compute_representative_points!(tri2)
 pts = get_points(tri2)
@@ -246,6 +263,21 @@ end
         end
     end
     tri2, label_map, index_map = simple_geometry()
+    # Manually populate boundary_vertex_to_ghost for tri2 since simple_geometry uses delete_ghosts=true
+    let bn = get_boundary_nodes(tri2)
+        bv_map = DT.get_boundary_vertex_to_ghost(tri2)
+        ghost = DT.𝒢
+        for curve_idx in 1:DT.num_curves(bn)
+            curve = DT.get_boundary_nodes(bn, curve_idx)
+            for section_idx in 1:DT.num_sections(curve)
+                section = DT.get_boundary_nodes(curve, section_idx)
+                for node in section
+                    bv_map[node] = ghost
+                end
+                ghost -= 1
+            end
+        end
+    end
     for (ghost_vertex, segment_index) in get_ghost_vertex_map(tri2)
         nodes = get_boundary_nodes(tri2, segment_index)
         for node in nodes
@@ -255,9 +287,9 @@ end
             @test res1 ∈ DT.get_ghost_vertex_range(tri2, ghost_vertex)
         end
     end
-    reduced_bn = reduce(vcat, reduce(vcat, get_boundary_nodes(tri2)))
-    for node in each_vertex(tri)
-        if node ∉ reduced_bn
+    reduced_bn2 = reduce(vcat, reduce(vcat, get_boundary_nodes(tri2)))
+    for node in each_vertex(tri2)
+        if node ∉ reduced_bn2
             flag, res = DT.is_boundary_node(tri2, node)
             @test !flag && res == DT.∅
         end

test/predicates/index_and_ghost_handling.jl

@@ -19,6 +19,7 @@ for PT in subtypes(DT.AbstractPredicateKernel)
         _tri.all_segments, _tri.weights, _tri.adjacent, _tri.adjacent2vertex, _tri.graph, _tri.boundary_curves,
         _tri.boundary_edge_map, _tri.ghost_vertex_map, _tri.ghost_vertex_ranges, DT.ConvexHull(_pts, _tri.convex_hull.vertices), _tri.representative_point_list,
         _tri.polygon_hierarchy, _tri.boundary_enricher, _tri.cache,
+        _tri.has_ghosts, copy(_tri.boundary_vertex_to_ghost),
     )
     DT.compute_representative_points!(tri)
     global rep = DT.get_representative_point_list(tri)

test/runtests.jl

@@ -163,4 +163,4 @@ end
         end
         @test isempty(missing_set)
     end
-end
+end