Topology-Preserving Simplification for Adjacent Zones

When two geofences share a border — adjacent delivery zones, abutting congestion-charge districts, neighbouring wildlife exclusion areas — simplifying each polygon on its own quietly rips their shared edge apart, leaving strips of no-man’s-land and double-counted overlap between them. This page expands the calibration work in Polygon Simplification for High-Throughput Streams and the broader spatial index lookup architecture, addressing the failure that the single-polygon comparison of Douglas-Peucker vs Visvalingam-Whyatt cannot see: what happens along the boundary two zones are supposed to share.

Concept and specification #

The root cause is that per-polygon simplification treats a shared edge twice, once from each polygon’s vertex ordering, and the two passes make different removal decisions. Zone A keeps vertex $v$; Zone B, walking the same edge in reverse with a different chord context, drops it. The two simplified borders no longer coincide: where A pulled inward you get a gap, where B pushed outward you get an overlap. A point in the gap belongs to neither zone; a point in the overlap belongs to both. In a geofencing system that decides billing, access, or dispatch, both are correctness bugs.

The reason the two passes disagree is not randomness but the direction-dependence of the underlying algorithm. Douglas-Peucker anchors its recursion at the two endpoints of a span and measures perpendicular distance from that chord; Zone A and Zone B present the shared edge with opposite winding, so their chords, split points, and recursion trees differ. Visvalingam-Whyatt is only marginally better here — the neighbour-triangle area of an interior shared vertex is identical from either side, but as soon as the two zones reach different vertex budgets the removal order diverges and the borders separate again. No amount of parameter tuning on the independent passes closes the seam, because the seam is created by running the decision twice rather than once. This is the structural argument for lifting simplification off the polygon and onto the boundary it shares.

The fix is to simplify the topology, not the polygons. Decompose the coverage into arcs — maximal chains of vertices shared by the same set of polygons — deduplicate them so each shared boundary exists exactly once, simplify each arc a single time, and reassemble the polygons from the simplified arcs. Because both neighbours reference the identical simplified arc, they cannot diverge. This is the model TopoJSON encodes on disk and the one a planar-graph builder reconstructs in memory.

Quantify the damage of the naive approach with two integrals over the boundary strip. If $B_A$ and $B_B$ are the simplified boundaries of adjacent zones and $S$ the true shared edge, the gap and overlap areas are

$$G={\iint }_{\Omega }1[,p\notin A\wedge p\notin B\wedge p\in \text{hull}(S),],dp,O={\iint }_{\Omega }1[,p\in A\wedge p\in B,],dp$$

Topology-preserving simplification drives both to exactly zero because $B_A$ and $B_B$ are, by construction, the same polyline along $S$.

Step-by-step implementation #

Prerequisites: Python 3.11+, topojson>=1.6, shapely>=2.0, geopandas>=0.14. Input is a set of adjacent zone polygons that share borders; output is the same set with borders simplified consistently.

1. Build a topology from the zone set. The topojson library detects shared arcs and junctions across all polygons in one pass, deduplicating every common border.

import geopandas as gpd
import topojson as tp

def build_topology(zones: gpd.GeoDataFrame) -> tp.Topology:
    """Extract shared arcs across all adjacent zones (junctions detected globally)."""
    # prequantize=False keeps full precision; shared_coords stitches touching borders
    return tp.Topology(zones, prequantize=False, shared_coords=True)

2. Simplify on the shared arcs, not the polygons. Because each border arc exists once, simplifying the topology touches each shared edge a single time.

def simplify_topology(topo: tp.Topology, eps_m: float) -> gpd.GeoDataFrame:
    """Simplify each arc once; both neighbours inherit the identical result."""
    simplified = topo.toposimplify(
        epsilon=eps_m,          # meters in a projected CRS
        simplify_algorithm="dp",  # DP applied per-arc, shared between zones
    )
    return simplified.to_gdf()  # reassemble polygons from simplified arcs

Gotcha: toposimplify only shares an arc if the neighbouring borders are coordinate-identical. If two zones were digitized separately with near-but-not-equal shared vertices, snap them to a common grid first (shapely.set_precision) or the arc detection misses the join and you are back to per-polygon behavior.

3. Snap near-coincident borders before topology extraction. Real-world zone data rarely shares exact coordinates; a grid snap makes the shared edge detectable.

from shapely import set_precision

def snap_borders(zones: gpd.GeoDataFrame, grid_m: float = 0.5) -> gpd.GeoDataFrame:
    """Round coordinates to a common grid so touching borders become identical."""
    zones = zones.copy()
    zones["geometry"] = zones.geometry.apply(lambda g: set_precision(g, grid_m))
    return zones                                       # now shared_coords can stitch them

4. Validate the seal before indexing. Confirm the simplified coverage has no gap and no overlap by unioning and differencing.

from shapely.ops import unary_union

def measure_seal(zones: gpd.GeoDataFrame) -> tuple[float, float]:
    """Return (gap_area, overlap_area) in m² across the whole coverage."""
    merged = unary_union(zones.geometry.tolist())      # dissolve into coverage
    covered = merged.area
    summed = float(zones.geometry.area.sum())          # sum of individual zone areas
    overlap = max(summed - covered, 0.0)               # double-counted area
    hull_gap = merged.convex_hull.area - covered       # coarse gap proxy inside hull
    return hull_gap, overlap

Benchmark and verification #

Take four adjacent municipal zones (~3,800 vertices total, ~6 km of shared borders), simplify with $\epsilon =5$ m both ways, and measure gap and overlap area plus sliver count.

import time

def bench_seal(zones: gpd.GeoDataFrame, eps_m: float) -> dict[str, float]:
    t0 = time.perf_counter()
    topo = build_topology(snap_borders(zones))
    out = simplify_topology(topo, eps_m)
    gap, overlap = measure_seal(out)
    return {"gap_m2": gap, "overlap_m2": overlap, "ms": (time.perf_counter() - t0) * 1e3}

Representative figures at $\epsilon =5$ m:

Topology-preserving simplification costs ~2.6x the wall time (junction detection and arc reassembly) but reduces gap and overlap area to exactly zero and eliminates every sliver and misclassification — the naive path’s 63 misclassified points out of 10k are the billing/access bugs this technique exists to remove. The extra cost is a one-time ingestion expense paid off the hot path, not a per-evaluation tax: once the arcs are simplified and the polygons reassembled, the runtime containment test sees ordinary rings and pays exactly the same reduced PIP latency the Douglas-Peucker vs Visvalingam-Whyatt comparison measures. Because junction detection dominates that 21 ms, batch the whole adjacency set through one topology build rather than paying it per polygon — a coverage of 400 zones amortizes junction detection across all shared borders in a single pass, so per-zone cost stays flat as the set grows.

Verify by asserting gap_m2 == 0 and overlap_m2 == 0 before promoting a batch, and shadow-diff containment answers against the un-simplified coverage. Track the arc count before and after the snap step as a regression guard: if snapping did not reduce the arc total, the shared borders were never stitched and every seam is still being cut twice regardless of what the gap assertion reports on a coarse sample.

Failure modes and edge cases #

• Sliver polygons. The classic artifact: a long, thin gap or overlap strip along a torn edge. unary_union on naively simplified zones leaves them as separate tiny geometries; count polygons whose area-to-perimeter ratio falls below a threshold and treat any nonzero count as a topology break.
• Missed shared arcs from coordinate drift. If neighbouring borders are not coordinate-identical, shared_coords never stitches them and every edge is simplified twice. Always snap to a common grid first, and assert the arc count dropped after snapping.
• Three-zone junctions. Where three or more zones meet at a point, the junction must be preserved as an arc endpoint or simplification can pull one border across the meeting point. topojson marks junctions automatically, but a manual planar-graph builder must detect degree-≥3 vertices and pin them.
• Over-snapping collapses geometry. Too coarse a snap grid merges distinct vertices and can degenerate a thin zone. Keep the grid well below the smallest real feature width, and validate each geometry with is_valid after snapping.
• Holes and enclaves. A zone fully contained in another (an enclave) shares its entire boundary; ensure interior rings are routed through the same arc-sharing path or the enclave border tears just like an external one.

Related #

• Polygon Simplification for High-Throughput Streams — parent overview of simplification latency budgets and topology repair.
• Douglas-Peucker vs Visvalingam-Whyatt for Geofence Boundaries — the per-polygon algorithms this page makes topology-aware.
• Handling Polygon Overlaps in Quadtree Partitions — how the index copes when zones do overlap.