cs.CG

Computing the diameter of the intersection graphs of objects is a basic problem in computational geometry. Previous works showed that the complexity of computing the diameter mainly depends on the object types: for unit disks and squares in 2D, the problem is solvable in truly subquadratic time, while for other objects, including unit segments and equilateral triangles in 2D or unit balls and axis-parallel unit cubes in 3D, there is no truly subquadratic time algorithm under the Orthogonal Vector (OV) hypothesis. We undertake a comprehensive study of computing the diameter of geometric intersection graphs for various types of objects. We discover many new irregularities, showing that the landscape is extremely nuanced: the source of hardness is a combination of the object type, the true diameter value, and how the objects intersect with each other. Our highlighted results for the 2D case include: 1. The diameter of non-degenerate, axis-aligned line segments can be computed in truly subquadratic time. Previous hardness result for line segments applies only to degenerate instances. On the other hand, for the degenerate case, we show that a truly subquadratic time algorithm exists when the true diameter is constant. 2. An almost-linear-time algorithm for unit-square graphs of constant diameter. Previous algorithms rely on succinct representation assuming bounded VC-dimension; for such a strategy $\Omega(n^{7/4})$ time is an inherent barrier. 3. An $\tilde{O}(n^{4/3})$-time algorithm to decide if the diameter of a unit-disk graph is at most 2. This improves upon the recent algorithm with running time $\tilde{O}(n^{2-1/9})$. 4. Deciding if the diameter of intersection graphs of fat triangles or line segments is at most 2 is truly subquadratic-hard under fine-grained complexity assumptions. Previous lower bounds only hold when deciding if diameter is at most 3.

Many topological data analysis (TDA) pipelines compute large collections of persistence diagrams, yet vectorizations and kernel methods discard the rank-induced implication relations among persistence intervals that are essential for faithful structural comparison and interpretability. We introduce higher-order persistence diagrams, a recursive construction in which containment relations among persistence intervals define higher-order persistence intervals. This construction performs comparison and aggregation directly on persistence diagrams and preserves interval-level structure. We use harmonic analysis to reduce frequency-space evaluations of aggregated diagrams to zeta transforms. This reduction avoids explicit construction of higher-order diagrams and replaces quadratic pair enumeration with nearly linear-time evaluation. Experiments on random network models show substantial speedups over explicit aggregation. Anonymized code is available at https://anonymous.4open.science/r/higher-order-persistence-8201.

Vertical decomposition is a widely used general technique for decomposing the cells of arrangements of semi-algebraic sets in ${{\mathbb R}}^d$ into constant-complexity subcells. In this paper, we settle in the affirmative a few long-standing open problems involving the vertical decomposition of substructures of arrangements for $d = 3, 4$. For example, we obtain sharp bounds on the complexity of the vertical decomposition of the complement of the union of a set of semi-algebraic regions of constant complexity in ${{\mathbb R}}^3$, and of the minimization diagram of a set of trivariate functions. These results lead to efficient algorithms for a variety of problems involving vertical decompositions, including algorithms for constructing the decompositions themselves and for constructing $(1/r)$-cuttings of substructures of arrangements. They also lead to a data structure for answering point-enclosure queries amid semi-algebraic sets in ${{\mathbb R}}^3$ and ${{\mathbb R}}^4$.

Persistence-based topological optimization deforms a point cloud $X \subset \mathbb{R}^d$ by minimizing objectives of the form $L(X) = \ell(\mathrm{Dgm}(X))$, where $\mathrm{Dgm}(X)$ is a persistence diagram. In practice, optimization is limited by two coupled issues: persistent homology is typically computed on subsamples, and the resulting topological gradients are highly sparse, with only a few anchor points receiving nonzero updates. Motivated by diffeomorphic interpolation, which extends sparse gradients to smooth ambient vector fields via Reproducing Kernel Hilbert Space (RKHS) interpolation, we propose a more scalable pipeline that improves both subsampling and gradient extension. We introduce subsampling via random slicing, a lightweight scheme that promotes iteration-wise geometric coverage and mitigates density bias. We further replace the costly kernel solve with a fast Nadaraya-Watson (NW) Gaussian convolution, producing a globally defined smooth update field at a fraction of the computational cost, while being more suited for topological optimization tasks. We provide theoretical guarantees for NW smoothing, including anchor approximation bounds and global Lipschitz estimates. Experiments in $2$D and $3$D show that combining random slicing with NW smoothing yields consistent speedups and improved objective values over other baselines on common persistence losses.

We study the problem of computing a shortest tour that visits a sequence of $k$ polygons $P_1,\dots, P_k$ with a total number of $n$ vertices. A tour is an oriented curve such that there exist points $p_i\in P_i$ for all $i$ where $p_i$ appears not after $p_{i+1}$. In a seminal paper, Dror, Efrat, Lubiw and Mitchell (STOC 2003) considered the problem under $L_2$ distance, and gave $\widetilde O(nk)$ and $\widetilde O(nk^2)$ algorithms for disjoint and intersecting convex polygons, respectively. In this paper, we consider the orthogonal setting (with orthogonal polygons and Manhattan distance) and obtain the following results: - a truly subquadratic $\widetilde O(n^{2-\frac{1}{48}})$ algorithm when consecutive polygons in the sequence are disjoint; - an $\widetilde O(n)$ algorithm for ortho-convex polygons when consecutive polygons are disjoint; - an $O(n)$ algorithm for axis-aligned rectangles; - $\widetilde O(n^2)$ and $\widetilde O(n^{1.5}k^2)$ algorithms without restrictions. Our algorithms build on a wide range of techniques, including additively weighted Voronoi diagrams, rectangle decompositions, persistent data structures, and dynamic distance oracles for weighted planar graphs.

In the imprecise geometry model, the input is an imprecise point set, which is a family of regions $F = (R_1, \ldots,R_n)$, where for each $R_i$ one may retrieve the true point $p_i \in R_i$. By preprocessing $F$, we can construct the output, in our case the Pareto front, on $P$ faster. We efficiently construct the Pareto front of an imprecise point set in the plane. Efficiency is interpreted in two ways: minimizing (i) the number of retrievals, and (ii) the computation time used to determine the set of regions that must be retrieved and to construct the Pareto front. We present an algorithm to construct the Pareto front for possibly overlapping rectangles that is \emph{instance-optimal} with respect to the number of retrievals, meaning that for every fixed input $(F, P)$, there is no algorithm that retrieves asymptotically fewer regions to compute the output. This is a strong algorithmic quality, as it means that our algorithm is competitive even to clairvoyant algorithms which know a correct guess of the output and only have to verify its correctness. In terms of algorithmic running time, instance-optimality is provably unobtainable. We instead present an algorithm which is within a $\log n$-factor of instance-optimality. This generalizes earlier results to overlapping input regions, at only a minor cost in running time. For unit squares, we present an algorithm that is not only instance-optimal in the number of retrievals, but also \emph{universally} optimal in terms of running time, meaning that for any fixed set of regions $F$, no algorithm has a better worst-case running time for all possible point sets $P$. This is the first universally optimal algorithm for overlapping planar input. Compared to previous work, this result improves the degree of overlap, the preprocessing time, the number of retrievals, and the running time.

Voronoi diagrams, and their more general weighted counterpart, power diagrams, are fundamental geometric constructs with wide-ranging applications. Recently, they have gained renewed attention in mesh-based neural rendering. Despite being extensively studied, the construction of 3D Voronoi diagrams for large-scale point sets remains computationally expensive, limiting their adoption in large-scale applications. Existing CPU-based approaches typically rely on computing its dual, the Delaunay tetrahedralization, but are prohibitively slow for large diagrams, while GPU-based methods either struggle to scale efficiently to large point sets or assume homogeneous point distributions. The weighted case, power diagrams, is even less explored in this context. Existing approaches are typically tailored to the application at hand, assuming homogeneous point distributions and small weight variations, making them unsuitable for general use in more complex heterogeneous data. In this paper, we present a highly parallelizable GPU algorithm for the fast construction of large-scale 3D Voronoi and power diagrams. Our approach constructs each convex cell from a weighted 3D point by progressively clipping an initial cell volume against bisecting planes induced by candidate neighboring points. To efficiently identify candidate neighbors under arbitrary spatial distributions, we introduce a culling criterion based on directional geometric bounds of the evolving cell, combined with a hierarchical best-first traversal of bounding volumes. We achieve performance on par with state-of-the-art Delaunay tetrahedralization methods on small and moderate problem sizes, while exhibiting robust scalability to large point sets and diverse spatial distributions. Moreover, our method naturally generalizes to power diagrams without additional assumptions. See https://research.zenseact.com/publications/paragram .

Continuous Dynamic Time Warping (CDTW) is a robust similarity measure for polygonal curves that has recently found a variety of applications. Despite its practical use, not much is known about the algorithmic complexity of computing it in 2D, especially when one requires either an exact solution or strong approximation guarantees. We fill this gap by introducing a $5$-approximation algorithm with running time $O(n^5)$ under the 1-norm. This is the first constant-factor approximation for 2D CDTW with polynomial running time. We extend our algorithm to all polygonal norms on $\mathbb{R}^2$, which we subsequently use in order to achieve a $(5+\varepsilon)$-approximation with time complexity $O(n^5 / \varepsilon^{1/2})$ for CDTW in 2D under any fixed norm. The latter result in particular includes the usual Euclidean 2-norm.

We study the computational complexity of exact cardinality-constrained minimum Riesz $s$-energy subset selection in finite metric spaces: given $n$ points, select $k<n$ points of minimum Riesz $s$-energy. The objective sums inverse-power pair interactions and therefore promotes well-separated subsets; as $s$ becomes large, it increasingly approaches a bottleneck criterion governed by the closest selected pair, linking it to minimum pairwise distance (MPD). Building on the general-metric NP-hardness result of Pereverdieva et al. (2025), we prove that NP-hardness persists for point sets in the Euclidean plane when $s$ is part of the input. In contrast, finite ultrametric spaces form an exact tractable regime: on rooted binary ultrametric trees with $n$ leaves, an optimal size-$k$ subset can be computed by dynamic programming in $O(nk^2)$ time. We also discuss the ordered one-dimensional Euclidean case, where the classical MPD objective admits simple dynamic programming, but the additive Riesz energy does not appear to allow the same state compression. Finally, we explain why one natural route to fixed-$s$ Euclidean hardness does not close: Fowler-style 3SAT gadgets, together with zeta-function bounds for far-field interactions, show why this approach still requires an exponent depending on $k$. Together, these results provide a compact complexity landscape for a natural diversity or dispersion objective, distinguishing Euclidean hardness, ultrametric tractability, and the ordered one-dimensional case.

Computing the convex hull of a planar $n$-point set $P$ is one of the most fundamental problems in computational geometry. It has an $\Omega(n \log n)$ lower bound in the algebraic computation tree model, and many convex hull algorithms match this bound. Classical results show that, under special input assumptions, sub-$O(n \log n)$ algorithms are possible. For instance, when the points are given in lexicographic or angular order, the convex hull can be computed in linear time. Even under the weaker assumption that the sequence of points corresponds to the ordered vertices of a simple polygonal chain, linear-time algorithms exist. This naturally raises the question: can the convex hull of a point set be computed in sub-$O(n \log n)$ time under weaker input assumptions? We answer this positively. Under the promise that the input sequence contains the convex hull as a subsequence, we give a deterministic $O(n \sqrt{\log n})$-time algorithm to compute the convex hull of $P$. With randomisation, we achieve expected running time $O(n \log^{\varepsilon} n)$ for any constant $\varepsilon > 0$. We find this surprising, as points not on the convex hull may behave adversarially toward our convex hull construction algorithm. Yet the promise that \emph{only} the hull points are sorted suffices for $o(n \log n)$-time algorithms. Finally, we show that this promise is tight: if it is even slightly broken, i.e., allowing just one hull point to appear out of order, we prove an adversarial $\Omega(n \log n)$-time lower bound. Consequently, the promise cannot be verified with fewer than $\Omega(n \log n)$ comparisons. This also negatively resolves an open problem of L\"{o}ffler and Raichel, who conjectured sub-$O(n \log n)$-time algorithms for computing the convex hull of a supersequence containing the hull as a subsequence.

Given a simple polygon $P$ with $n$ vertices, we consider the problem of constructing a data structure for visibility queries: for any query point $q \in P$, compute the visibility polygon of $q$ in $P$. To obtain $O(\log n + k)$ query time, where $k$ is the size of the visibility polygon of $q$, the previous best result requires $O(n^3)$ space. In this paper, we propose a new data structure that uses $O(n^{2+\epsilon})$ space, for any $\epsilon > 0$, while achieving the same query time. If only $O(n^2)$ space is available, the best known result provides $O(\log^2 n + k)$ query time. We improve this to $O(\log n \log \log n + k)$ time. When restricted to $o(n^2)$ space, the only previously known approach, aside from the $O(n)$-time algorithm that computes the visibility polygon without preprocessing, is an $O(n)$-space data structure that supports $O(k \log n)$-time queries. We construct a data structure using $O(n \log n)$ space that answers visibility queries in $O(n^{1/2+\epsilon} + k)$ time. In addition, for the special case in which $q$ lies on the boundary of $P$, we build a data structure of $O(n \log n)$ space supporting $O(\log^2 n + k)$ query time; alternatively, we achieve $O(\log n + k)$ query time using $O(n^{1+\epsilon})$ space. To achieve our results, we propose a new method for decomposing simple polygons, which may be of independent interest.

We answer an open problem in the \emph{American Mathematical Monthly} about covering large triangles. Given a triangle $T$ of any triangular shape with a selected side length between $n \in \mathbb{N}$ and $n+1$, Baek and Lee proved that $T$ could not be covered with $n^2+1$ homothetic unit triangles (with the selected side of length 1). Letting $T_{n+d}$ denote a triangle with selected side length $n + d$ with $d \in (0, 1)$, Baek and Lee extended their proof to establish upper bounds for $d$ above which a $T_{n+d}$ cannot be covered with $n^2+2$ or $n^2+3$ homothetic unit triangles. Then, they showed that these bounds are tight based on analyses of a method by Conway and Soifer for the $n^2+2$ case and their own method for the $n^2+3$ case. Baek and Lee stated as an open problem the need to find tight upper bounds for the $n^2 + k$ cases for $4 \le k \le 2n$. We extend the Baek and Lee proof to establish upper bounds for those higher cases, and we show the upper bounds are tight by presenting two new triangle covering methods for the odd and even cases of $k$ that meet the upper bounds, as well as an optimal consolidated method that uses whichever of the two will cover a given $T_{n+d}$ with the fewest homothetic unit triangles.

k-nearest neighbor (k-NN) search is a fundamental primitive in geometry processing and computer graphics. While spatial partitioning structures such as kd-trees are standard, they are often manifold-blind, failing to exploit the intrinsic low-dimensional structure of points sampled from 2-manifolds. Recent advances in dynamic programming-based nearest neighbor search (DP-NNS) leverage incrementally constructed Voronoi diagrams to accelerate queries, where each site p maintains a list of successors that progressively refine its Voronoi cell. However, DP-NNS is restricted to single nearest neighbor (k=1) searches, precluding their adoption in applications that require local neighborhood statistics. In this paper, we generalize the DP-NNS framework to support arbitrary k-NN queries for manifold-aligned data. Our approach is founded on the geometric observation that if p_i is the nearest neighbor of a query q in P, then the second nearest neighbor of q must reside either within the prefix set P_{1:i-1} = {p_1, \dots, p_{i-1}} or within p_i's successor list. By recursively extending this principle, we introduce Manifold k-NN, a recursive algorithmic scheme that significantly outperforms conventional kd-trees for manifold-aligned data. Our method achieves a 1\times--10\times speedup in volume-to-surface query scenarios and inherently supports dynamic prefix queries -- enabling k-NN searches within any subset P_{1:m} (m \leq n) with zero overhead. Furthermore, we extend the framework to support point deletion via local Delaunay updates, providing a complete suite of dynamic operations for point set modification. Comprehensive experiments on diverse geometric datasets demonstrate the efficiency and broad applicability of our approach for modern graphics pipelines. Source code is available at https://github.com/sssomeone/manifold-knn.

We give a greedy sweepline algorithm for a Jordan curve and prove that it is maximal in the sense of [1]. Our proof uses K\H{o}nig's lemma.

We study the Witness Set problem, a natural dual to the classical Art Gallery problem. In the Witness Set problem, we are given a polygon $P$ and an integer $k$ as input, and the objective is to determine whether $P$ has a witness set of size at least $k$. A point set $X$ in $P$ is called a witness set if every point in $P$ is visible from at most one point in $X$. For simple polygons, we show that Witness Set lies in both $NP$ and $XP$. This stands in sharp contrast to its dual, the Art Gallery problem, which was recently shown to be $\exists \mathbb{R}$-complete by Abrahamsen et al. and is therefore neither in $NP$ nor admits a polynomial-size discretization unless $NP=\exists \mathbb{R}$. In contrast, we prove that Witness Set for simple polygons admits a finite discretization of size $n^{f(k)}$ for some function $f$. For comparison, even for simple polygons, Efrat and Har-Peled gave an algorithm for Art Gallery running in time $n^{O(k)}$ using tools from real algebraic geometry, and it appears difficult to obtain such algorithms without this machinery. On the other hand, our approach for Witness Set is purely combinatorial and relies on discretization, leading to an $n^{f(k)}$-time algorithm. Although Amit et al. claimed more than fifteen years ago that Witness Set is $NP$-hard, no proof or reference was provided. We show that the discrete version of the Witness Set problem - where the witness set must be chosen from a given finite point set $Q$ (instead of allowing witnesses to be chosen anywhere in the polygon), referred to as Discrete Witness Set - is $NP$-complete, even when the input is restricted to rectilinear polygons with holes. However, for simple polygons, Discrete Witness Set admits a polynomial-time algorithm by Das et al. Thus, it remains an open question whether the Witness Set problem is $NP$-hard.

This paper establishes a rigorous geometrical framework for spherical origami, origami using spherical sheets based on spherical geometry. Two settings are treated: origami restricted to the unit sphere ($\mathbb{S}^2$), and three-dimensional folding of spherical sheets in space. For origami on $\mathbb{S}^2$, the definitions of Euclidean origami are systematically extended to the spherical setting, and all seven Huzita--Justin axioms are shown to admit explicit equations in spherical geometry. For three-dimensional folding, equidistant curves are introduced as fold curves, replacing geodesics and enabling a richer family of folds. The framework is validated by successfully constructing computer graphics of spherical origami birds, demonstrating both the theoretical completeness and practical utility of the proposed approach.

Let $S$ be a set of $n$ points in $\mathbb{R}^2$. Our goal is to preprocess $S$ to efficiently compute the smallest enclosing disk of the points in $S$ that lie inside an axis-aligned query rectangle. Previous data structures for this problem achieve a query time of $O(\log^6 n)$ with $O(n \log^2 n)$ preprocessing time and space by lifting the points to 3D, dualizing them into polyhedra, and searching through their intersections. We present a significantly simpler approach, solely based on 2D geometric structures, specifically 2D farthest-point Voronoi diagrams. Our approach achieves a deterministic query time of $O(\log^4 n)$ and, via randomization, an expected query time of $O(\log^{5/2} n \log\log n)$ with the same preprocessing bounds.

We consider upward-planar layered drawings of directed graphs, i.e., crossing-free drawings in which each edge is drawn as a y-monotone curve going upward from its tail to its head, and the y-coordinates of the vertices are integers. The span of an edge in such a drawing is the absolute difference between the y-coordinates of its endpoints, and the span of the drawing is the maximum span of any edge. The span of an upward-planar graph is the minimum span over all its upward-planar drawings. We study the problem of determining the span of upward-planar graphs and provide both combinatorial and algorithmic results. On the combinatorial side, we present upper and lower bounds for the span of directed trees. On the algorithmic side, we show that the problem of determining the span of an upward-planar graph is NP-complete already for directed trees and for biconnected single-source graphs. Moreover, we give efficient algorithms for several graph families with a bounded number of sources, including st-planar graphs and graphs where the planar or upward-planar embedding is prescribed. Furthermore, we show that the problem is fixed-parameter tractable with respect to the vertex cover number and the treedepth plus the span.

This paper studies exact fixed-cardinality Solow--Polasky diversity subset selection on ordered finite $\ell_1$ point sets, with monotone biobjective Pareto fronts and their higher-dimensional staircase analogues as central applications. Solow--Polasky diversity was introduced in biodiversity conservation, whereas the same inverse-matrix expression appears in metric geometry as magnitude: for a finite metric space $(X,d)$ with exponential similarity matrix $Z_{ij}=e^{-q d(x_i,x_j)}$, the quantity $\1^\top Z^{-1}\1$ is the magnitude of the scaled finite metric space $(X,qd)$ whenever the weighting is defined by the inverse matrix. Thus, in this finite exponential-kernel setting, Solow--Polasky diversity and magnitude are mathematically the same object viewed through different motivations. Building on the linear-chain magnitude formula of Leinster and Willerton, the paper gives a detailed proof of the scaled consecutive-gap identity $ \SP(X)=1+\sum_r \tanh(qg_r/2),$ where the $g_r$ are the gaps between consecutive selected points. It then proves an exact Bellman-recursion theorem for maximizing this value over all subsets of a prescribed cardinality, yielding an $O(kn^2)$ dynamic program for an ordered $n$-point candidate set and subset size $k$. Finally, the paper proves ordered $\ell_1$ reductions showing that the same algorithm applies to monotone biobjective Pareto-front approximations and, more generally, to finite coordinatewise monotone $\ell_1$ staircases in $\R^d$. These are precisely the ordered $\ell_1$ chains for which the $\ell_1$-distance becomes a line metric along the chosen order, so the one-dimensional dynamic program applies without modification. Keywords: Solow--Polasky diversity; magnitude; metric geometry; dynamic programming; ordered points; $\ell_1$ geometry; Pareto-front approximation.

We prove the Jordan curve theorem by generalizing the sweepline algorithm for trapezoidal decomposition of a polygon. Our proof uses Zorn's lemma (or, equivalently the axiom of choice). Though several proofs have been given for the Jordan curve theorem by various authors, ours is the {\bf first algorithmic proof} of Jordan curve theorem using computational geometry. Further, with some preparation, the proof can be taught as part of an undergraduate discrete mathematics course, where till now the proof of this theorem was considered inaccessible.

In the (Nesting) Bird Box Problem we are given a polygonal domain P and a number k and we want to know if there is a set B of k points inside P such that no two points in B can see each other. The underlying idea is that each point represents a birdhouse and many birds only use a birdhouse if there is no other occupied birdhouse in its vicinity. We say two points a,b see each other if the open segment ab intersects neither the exterior of P nor any vertex of P. We show that the Nesting Bird Box problem is ER-complete. The complexity class ER can be defined by the set of problems that are polynomial time equivalent to finding a solution to the equation $p(x) = 0$, with $x\in R^n$ and $p\in $Z[X_1,...,X_n]$. The proof builds on the techniques developed in the original ER-completeness proof of the Art Gallery problem. However our proof is significantly shorter for two reasons. First, we can use recently developed tools that were not available at the time. Second, we consider polygonal domains with holes instead of simple polygons.

We study detection of collapse in high-dimensional point clouds, where mass concentrates near a lower-dimensional set relative to a non-collapsed geometry. We propose persistent homology-based test statistics under two well-studied filtrations, with cutoffs calibrated under a broad set of non-collapsed reference models. We benchmark power across three alternative collapse mechanisms (linear/spectral, nonlinear-support, and contamination/heterogeneity) and distill the results into a mechanism map guiding the choice of filtration and statistic.

Symmetry is an implicit objective in structural form-finding that often reconciles efficiency and aesthetics. This paper identifies the symmetry of polyhedral diagrams in three-dimensional graphic statics (3DGS) as point groups and formulates them as constraints, enabling the optimization and manipulation of polyhedral diagrams that preserve such symmetry. 3DGS has been an efficient and effective tool for the form-finding of funicular structures. However, when modifying complex diagrams for design exploration or optimization, one can easily break the symmetry of the reciprocal design input, rendering the result undesirable for practical use. To address this problem, this paper investigates symmetry transformations and introduces point groups, an abstract algebra tool commonly used in crystallography to represent the symmetry and equivalence between a network of atoms (points with labels). It then discusses the hierarchy of symmetry in the geometry types of a polyhedral diagram, and proposes the constraint of symmetry through edge lengths. Based on the crystal symmetry search algorithm by spglib and pymatgen, a fast fingerprinting algorithm is developed to identify the point group of a polyhedral diagram and sort equivalent edges into sets. Finally, the paper shows that the necessary and sufficient condition for preserving the point group symmetry is that each set of edges has the same length. This constraint is compatible with the algebraic formulation of 3DGS and effectively preserves symmetry while reducing the dimension of the solution space. The method is implemented in the PolyFrame 2 plug-in for Rhino and Grasshopper.

We maintain a $(1+\varepsilon)$-spanner over the disk intersection graph of a dynamic set of disks. We restrict all disks to have their diameter in $[4,\Psi]$ for some fixed and known $\Psi$. The resulting $(1+\varepsilon)$-spanner has size $O(n \varepsilon^{-2} \log \Psi \log (\varepsilon^{-1}))$, where $n$ is the present number of disks. We develop a novel use of persistent data structures to dynamically maintain our $(1+\varepsilon)$-spanner. Our approach requires $O(\varepsilon^{-2} n \log^4 n \log \Psi)$ space and has an $O( \left( \frac{\Psi}{\varepsilon} \right)^2 \log^4 n \log^2 \Psi \log^2 (\varepsilon^{-1}))$ expected amortised update time. For constant $\varepsilon$ and $\Psi$, this spanner has near-linear size, uses near-linear space and has polylogarithmic update time. Furthermore, we observe that for any $\varepsilon < 1$, our spanner also serves as a connectivity data structure. With a slight adaptation of our techniques, this leads to better bounds for dynamically supporting connectivity queries in a disk intersection graph. In particular, we improve the space usage when compared to the dynamic data structure of (Baumann et al., DCG'24), replacing the linear dependency on $\Psi$ by a polylogarithmic dependency. Finally, we generalise our results to $d$-dimensional hypercubes.

The Fr\'echet distance is a popular distance measure between trajectories or curves in space, or between walks in graphs. We study computing the Fr\'echet distance between walks in the $d$-dimensional grid graphs, i.e. $\mathbb{Z}^d$ where points share an edge if they differ by one in one coordinate. We give an algorithm, that for two simple paths on $n$ vertices, $(1+\varepsilon)$-approximates the Fr\'echet distance in time $\widetilde{O}((\frac{n}{\varepsilon})^{2-2/d} +n)$. We complement this by a near-matching fine-grained lower bound: for constant dimensions $d \geq 3$, there is no $O((\varepsilon^{2/d}(\frac{n}{\varepsilon})^{2-2/d})^{1-\delta})$ algorithm for any $\delta>0$ unless the Orthogonal Vector Hypothesis fails. Thus, our results are tight up to a factor $\varepsilon^{2/d}$ and $\log(n)$-factors. We extend our results to imbalanced lower and upper bounds, where the curves have $n$ and $m$ vertices respectively, and also obtain near-tight bounds. Driemel, Har-Peled and Wenk [DCG'12] studied \emph{realistic assumptions} for curves to speed up Fr\'echet distance computation. One of these assumptions is $\lambda$-low density and they can compute a $(1+\varepsilon)$-approximation between $\lambda$-low dense curves in time $\widetilde{O}( \varepsilon^{-2} \lambda^2 n^{2(1-1/d)})$. By adapting our lower bound, we show that their algorithm has a tight dependency on $n$ and a tight dependency on $\varepsilon$ as $d$ goes to infinity. A gap remains in terms of $\lambda$.

We study the common intersection of arrangements of double-wedges. We consider arrangements where double-wedges may be either bowties (which do not contain a vertical line) or hourglasses (which contain a vertical line), in contrast to earlier studies that focused on arrangements of only bowties. This generalization changes the setting drastically, in particular, with respect to all arguments involving the point-line duality. Namely, a point in the intersection of all double-wedges is equivalent to a line that stabs a set of segments $\mathcal{S}$ (corresponding to the bowties) while it avoids a different set of segments $\mathcal{A}$ (corresponding to the complement of the hourglasses). We show that in this general setting, the intersection of $n$ double-wedges may consist of $\Omega(n^2)$ interior-disjoint regions. Further, we discuss Gallai-type results for arrangements of segments and anti-segments, and we provide algorithms for computing the intersection of such arrangements with worst-case optimal running time. Finally, we also prove that we can find a single intersection point in almost optimal running time, assuming that 3SUM admits no truly subquadratic-time algorithm.

We study SINGLE-SOURCE SHORTEST PATH (SSSP) on unweighted intersection graphs whose node set corresponds to a set of $n$ constant-complexity objects in the plane. We prove SSSP can be solved in $O(U(n)\ \mathrm{polylog}\,n)$ expected time for any class of objects whose union complexity is $U(n)$. In particular, we obtain an $O(n 2^{\alpha(n)}\log^2 n)$ algorithm for arbitrary pseudodisks, and an $O(\lambda_{s+2}(n)2^{O(\log^* n)} \log^2 n)$ algorithm for locally fat objects. This significantly extends the class of objects for which SSSP can be solved in $O(n\ \mathrm{polylog}\, n)$ time: so far, $O(n\ \mathrm{polylog}\, n)$ SSSP algorithms were not even known for pseudodisks that are fat and convex and similarly-sized. Our second result concerns the DIAMETER problem, which asks for the maximum distance between any two nodes in a graph. Even for intersection graphs, near-quadratic algorithms are difficult to obtain, and the $O(n^2\ \mathrm{polylog}\, n)$ running time that follows from our SSSP algorithm is the first near-quadratic running time for such general classes of intersection graphs. Obtaining subquadratic running time is even more challenging. We prove that the diameter of a set of arbitrary pseudodisks can be computed almost exactly, namely up to an additive error of 2, in $\tilde{O}(n^{2-1/14})$ expected time. This generalizes and speeds up a recent algorithm by Chang, Gao, and Le~(SoCG 2024) that works for similarly-sized disks (or similarly-sized pseudodisks that are fat and satisfy a strong monotonicity assumption) and runs in $\tilde{O}(n^{2-1/18})$ time. To this end, we develop a so-called star-based $r$-clustering for intersection graphs of pseudodisks, which is interesting in its own right. Our star-based $r$-clustering can also be used to obtain an almost exact distance oracle for pseudodisks that uses $O(n^{2-1/13})$ storage and has $O(1)$ query time.

Consider a stream of $n$ random points (say, from the unit square) arriving one by one, where a player has to make an irreversible immediate decision for each arriving point whether to pick it. The player has to pick a single point, and the payoff is the area of the cell of the picked point, in the final Voronoi diagram of \emph{all} the points. We show that there is a simple strategy so that with probability $\geq 1 - \tilde O(1/\sqrt{n})$, the player's payoff is only a constant factor smaller than the optimal choice (i.e., the one made by the prophet). This competitiveness is somewhat surprising, as this payoff is larger by a factor of $\Theta( \log n)$ than the average payoff.

We study the exact counting problem for all lattice rectangles contained in the square $[0,n)\times[0,n)$, including non-axis-parallel ones. Starting from the standard parametrization by a primitive direction $(u,v)$ and two side lengths, we derive several exact algorithms: the classical $O(n^2)$ sweep, decompositions of complexity $O(n^{3/2}\log n)$ and $O(n^{4/3}\log n)$, a ten-moment weighted-floor-sum reduction of complexity $O(n\log^3 n)$, and a divisor-layer algorithm with the complexity $O(n\log^2 n)$. We also give an all-values algorithm that computes $F(1),\ldots,F(N)$ in $O(N^{3/2})$ arithmetic operations. The main idea behind the near-linear one-value algorithms is to reduce the geometric summation to constant-size families of weighted floor sums closed under Euclidean-style affine and reciprocal transformations. Besides the exact algorithmic results, we derive a two-term asymptotic expansion, $F(n)=\frac{4\log 2-1}{\pi^2}n^4\log n+B\,n^4+o(n^4)$ with the explicit formula for $B$, which provides an independent consistency check for the large-$n$ numerical data produced by the algorithms.

Origami structures, particularly Miura-ori patterns, offer unique capabilities for surface approximation and deployable designs. In this study, a constrained mapping optimization algorithm is designed for designing surface-aligned Miura-ori via a narrow band approximation of the input surface. The Miura-fold, embedded in the narrow band, is parameterized to a planar domain, and a mapping is computed on the parameter pattern by optimizing certain energy terms and constraints. Extensive experiments are conducted, showing the significance and flexibility of our methods.

We prove that, for every fixed $\theta_0>0$, selecting a subset of prescribed cardinality that maximizes the Solow--Polasky diversity indicator is NP-hard for finite point sets in $\mathbb{R}^2$ with the Euclidean metric, and therefore also for finite point sets in $\mathbb{R}^d$ for every fixed dimension $d \ge 2$. This strictly strengthens our earlier NP-hardness result for general metric spaces by showing that hardness persists under the severe geometric restriction to the Euclidean plane. At the same time, the Euclidean proof technique is different from the conceptually easier earlier argument for arbitrary metric spaces, and that general metric-space construction does not directly translate to the Euclidean setting. In the earlier proof one can use an exact construction tailored to arbitrary metrics, essentially exploiting a two-distance structure. In contrast, such an exact realization is unavailable in fixed-dimensional Euclidean space, so the present reduction requires a genuinely geometric argument. Our Euclidean proof is based on two distance thresholds, which allow us to separate yes-instances from no-instances by robust inequalities rather than by the exact construction used in the general metric setting. The main technical ingredient is a bounded-box comparison lemma for the nonlinear objective $\mathbf{1}^{\top}Z^{-1}\mathbf{1}$, where $Z_{ij}=e^{-\theta_0 d(x_i,x_j)}$. This lemma controls the effect of perturbations in the pairwise distances well enough to transfer the gap created by the reduction. The reduction is from \emph{Geometric Unit-Disk Independent Set}. We present the main argument in geometric form for finite subsets of $\mathbb{R}^2$, with an appendix supplying the bit-complexity details needed for polynomial-time reducibility.

We study a family of local depth-based corrections to maxmin landmark selection for lazy witness persistence. Starting from maxmin seeds, we partition the cloud into nearest-seed cells and replace or move each seed toward a deep representative of its cell. The principal implemented variant, \emph{support-weighted partial recentering}, scales the amount of movement by cell support. The contributions are both mathematical and algorithmic. On the mathematical side, we prove local geometric guarantees for these corrections: a convex-core robustness lemma derived from halfspace depth, a $2r$ cover bound for subset recentering, and projected cover bounds for the implemented partial-recentering rules. On the algorithmic side, we identify a practically effective variant through a layered empirical study consisting of planar synthetic benchmarks, a parameter-sensitivity study, and an MPEG-7 silhouette benchmark, together with a modest three-dimensional torus extension. The main planar experiments show that support-weighted partial recentering gives a consistent geometric improvement over maxmin while preserving the thresholded $H_1$ summary used in the study. The three-dimensional experiment shows the same geometric tendency but only mixed topological behavior. The paper should therefore be read as a controlled study of a local depth-based alternative to maxmin, rather than as a global witness-approximation theorem or a claim of uniform empirical superiority.

The minimum convex cover problem seeks to cover a polygon $P$ with the fewest convex polygons that lie within $P$. This problem is $\exists\mathbb R$-complete, and the best previously known algorithm, due to Eidenbenz and Widmayer (2001), achieves an $O(\log n)$-approximation in $O(n^{29} \log n)$ time, where $n$ is the complexity of $P$. In this work we present a novel approach that preserves the $O(\log n)$ approximation guarantee while significantly reducing the running time. By discretizing the problem and formulating it as a set cover problem, we focus on efficiently finding a convex polygon that covers the largest number of uncovered regions, in each iteration of the greedy algorithm. This core subproblem, which we call the rotten potato peeling problem, is a variant of the classic potato peeling problem. We solve it by finding maximum weighted paths in Directed Acyclic Graphs (DAGs) that correspond to visibility polygons, with the DAG construction carefully constrained to manage complexity. Our approach yields a substantial improvement in the overall running time and introduces techniques that may be of independent interest for other geometric covering problems.

In this paper, we study the many-to-many matching problem on planar point sets with integer coordinates: Given two disjoint sets $R,B \subset [\Delta]^2$ with $|R|+|B|=n$, the goal is to select a set of edges between $R$ and $B$ so that every point is incident to at least one edge and the total Euclidean length is minimized. In the general case that $R$ and $B$ are point sets in the plane, the best-known algorithm for the many-to-many matching problem takes $\tilde{O}(n^2)$ time. We present an exact $\tilde{O}(n^{1.5} \log \Delta)$ time algorithm for point sets in $[\Delta]^2$. To the best of our knowledge, this is the first subquadratic exact algorithm for planar many-to-many matching under bounded integer coordinates.

We consider a new treatment for making polyhedron nets referred to as ``apple peel unfolding'': drawing the nets as if we were peeling off appleskins. We define apple peel unfolding strictly and implement a program that derives the sequential selection of the polyhedral faces for a target polyhedron in accordance with the definition. Consequently, the program determines whether the polyhedron is peelable (can be peeled completely). We classify Archimedean solids and their duals (Catalan solids) as perfect (always peelable), possible (peelable for restricted cases), or impossible. The results show that three Archimedean and six Catalan solids are perfect, and three Archimedean and three Catalan ones are possible.

We study the patient zero problem in epidemic spreading processes in the independent cascade model and propose a geometric approach for source reconstruction. Using Johnson-Lindenstrauss projections, we embed the contact network into a low-dimensional Euclidean space and estimate the infection source as the node closest to the center of gravity of infected nodes. Simulations on Erd\H{o}s-R\'enyi graphs demonstrate that our estimator achieves meaningful reconstruction accuracy despite operating on compressed observations.

We study a variant of a polygon partition problem, introduced by Chung, Iwama, Liao, and Ahn [ISAAC'25]. Given orthogonal unit vectors $\mathbf{u},\mathbf{v}\in \mathbb{R}^2$ and a polygon $P$ with $n$ vertices, we partition $P$ into connected pieces by cuts parallel to $\mathbf{v}$ such that each resulting subpolygon has width at most one in direction $\mathbf{u}$. We consider the value version, which asks for the minimum number of strips, and the reporting version, which outputs a compact encoding of the cuts in an optimal strip partition. We give efficient algorithms and lower bounds for both versions on three classes of polygons of increasing generality: convex, simple, and self-overlapping. For convex polygons, we solve the value version in $O(\log n)$ time and the reporting version in $O\!\left(h \log\left(1 + \frac{n}{h}\right)\right)$ time, where $h$ is the width of $P$ in direction $\mathbf{u}$. We prove matching lower bounds in the decision-tree model, showing that the reporting algorithm is input-sensitive optimal with respect to $h$. For simple polygons, we present $O(n \log n)$-time, $O(n)$-space algorithms for both versions and prove an $\Omega(n)$ lower bound. For self-overlapping polygons, we extend the approach for simple polygons to obtain $O(n \log n)$-time, $O(n)$-space algorithms for both versions, and we prove a matching $\Omega(n \log n)$ lower bound in the algebraic computation-tree model via a reduction from the $\delta$-closeness problem. Our approach relies on a lattice-theoretic formulation of the problem. We represent strip partitions as antichains of intervals in the Clarke--Cormack--Burkowski lattice, originally developed for minimal-interval semantics in information retrieval. Within this lattice framework, we design a dynamic programming algorithm that uses the lattice operations of meet and join.

In the classical online model, the maximum independent set problem admits an $\Omega(n)$ lower bound on the competitive ratio even for interval graphs, motivating the study of the problem under additional assumptions. We first study the problem on graphs with a bounded independent kissing number $\zeta$, defined as the size of the largest induced star in the graph minus one. We show that a simple greedy algorithm, requiring no geometric representation, achieves a competitive ratio of $\zeta$. Moreover, this bound is optimal for deterministic online algorithms and asymptotically optimal for randomized ones. This extends previous results from specific geometric graph families to more general graph classes. Since this bound rules out further improvements through randomization alone, we investigate the power of randomization with access to geometric representation. When the geometric representation of the objects is known, we present randomized online algorithms with improved guarantees. For unit ball graphs in $\mathbb{R}^3$, we present an algorithm whose expected competitive ratio is strictly smaller than the deterministic lower bound implied by the independent kissing number. For $\alpha$-fat objects and for axis-aligned hyper-rectangles in $\mathbb{R}^d$ with bounded diameters, we obtain algorithms with expected competitive ratios that depend polylogarithmically on the ratio between the maximum and minimum object diameters. In both cases, the randomized lower bound implied by the independent kissing number grows polynomially with the ratio between the maximum and minimum object diameters, implying substantial performance guarantees for our algorithms.

Many algorithmic problems can be solved (almost) as efficiently in metric spaces of bounded doubling dimension as in Euclidean space. Unfortunately, the metric space defined by points in a simple polygon equipped with the geodesic distance does not necessarily have bounded doubling dimension. We therefore study the doubling dimension of fat polygons, for two well-known fatness definitions. We prove that locally-fat simple polygons do not always have bounded doubling dimension, while any $(\alpha,\beta)$-covered polygon does have bounded doubling dimension (even if it has holes). We also study the perimeter of geodesically convex sets in $(\alpha,\beta)$-covered polygons (possibly with holes), and show that this perimeter is at most a constant times the Euclidean diameter of the set. Using these two results, we obtain new results for several problems on $(\alpha,\beta)$-covered polygons, including an algorithm that computes the closest pair of a set of $m$ points in an $(\alpha,\beta)$-covered polygon with $n$ vertices that runs in $O(n + m\log{n})$ expected time.

Streamlines have been widely used to represent and analyze various steady vector fields. To sufficiently represent important features in complex vector fields (like flow), a large number of streamlines are required. Due to the lack of a rigorous definition of features or patterns in streamlines, user interaction and exploration are required to achieve effective interpretation. Existing approaches based on clustering or pattern search, while valuable for specific analysis tasks, often face challenges in supporting interactive and level-of-detail exploration of large-scale curve-based data, particularly when real-time parameter adjustment and iterative refinement are needed. To address this, we design and implement an interactive web-based system. Our system utilizes a Curve Segment Neighborhood Graph (CSNG) to encode the neighboring relationships between curve segments. CSNG enables us to adapt a fast community detection algorithm to identify coherent flow structures and spatial groupings in the streamlines interactively. CSNG also supports a multi-level exploration through an enhanced force-directed layout. Furthermore, our system integrates an adjacency matrix representation to reveal detailed inter-relations among segments. To achieve real-time performance within a web browser, our system employs matrix compression for memory-efficient CSNG storage and parallel processing. We have applied our system to analyze and interpret complex patterns in several streamline datasets. Our experiments show that we achieve real-time performance on datasets with hundreds of thousands of segments.

Cartograms are a technique for visually representing geographically distributed statistical data, where values of a numerical attribute are mapped to the size of geographic regions. Contiguous cartograms preserve the adjacencies of the original regions during the mapping. To be useful, contiguous cartograms also require approximate preservation of shapes and relative positions. Due to these desirable properties, contiguous cartograms are among the most popular ones. Most methods for constructing contiguous cartograms exploit a deformation of the original map. Aiming at the preservation of geographical properties, existing approaches are often algorithmically cumbersome and computationally intensive. We propose a novel deformation technique for computing time-varying contiguous cartograms based on integral images evaluated for a series of discrete density distributions. The density textures represent the given dynamic statistical data. The iterative application of the proposed mapping smoothly transforms the domain to gradually equalize the temporal density, i.e., region areas grow or shrink following their evolutionary statistical data. Global shape preservation at each time step is controlled by a single parameter that can be interactively adjusted by the user. Our efficient GPU implementation of the proposed algorithm is significantly faster than existing state-of-the-art methods while achieving comparable quality for cartographic accuracy, shape preservation, and topological error. We investigate strategies for transitioning between adjacent time steps and discuss the parameter choice. Our approach applies to comparative cartograms' morphing and interactive cartogram exploration.

In the realm of computer-aided design (CAD) software, the intersection of B-spline surfaces stands as a fundamental operation. Despite the extensive history of surface intersection algorithms, the challenge of handling complex intersection topologies persists. While subdivision algorithms have demonstrated strong robustness in computing surface/surface intersection and are capable of addressing singular cases, determining the topology of the intersection obtained through these methods is a key factor for calculating correct intersection, and remains a difficult issue. To address this challenge, we propose a Mapper-based method for determining the topology of the intersection between two B-spline surfaces. Our algorithm is designed to efficiently handle various common and complex intersection topologies. Experimental results verify the robustness and topological correctness of this method.

For a set $\mathcal{D}$ of disks in the plane, its disk graph $G(\mathcal{D})$ is the graph with vertex set $\mathcal{D}$, where two vertices are adjacent if and only if the corresponding disks intersect. Given a set $\mathcal{D}$ of $n$ weighted disks, computing a maximum independent set of $G(\mathcal{D})$ is NP-hard. In this paper, we present an $O(n^3\log n)$-time algorithm for this problem in a special setting in which the disks are in convex position, meaning that every disk appears on the convex hull of $\mathcal{D}$. This setting has been studied previously for disks of equal radius, for which an $O(n^{37/11})$-time algorithm was known. Our algorithm also works in the weighted case where disks have weights and the goal is to compute a maximum-weight independent set. As an application of our result, we obtain an $O(n^3\log^2 n)$-time algorithm for the dispersion problem on a set of $n$ disks in convex position: given an integer $k$, compute a subset of $k$ disks that maximizes the minimum pairwise distance among all disks in the subset.

We consider the parametric shortest paths problem in a linearly interpolated graph. Given two positively-weighted directed graphs $G_0=(V,E,\omega_0)$ and $G_1=(V,E,\omega_1),$ the linearly interpolated graph is the family of graphs $(1-\lambda)G_0+\lambda G_1$, parameterized by $\lambda\in [0,1]$. The problem is to compute all distinct parametric shortest paths. We compute a data structure in $\Theta(k|E|\log |V|)$ time, where~$k$ is the number of distinct parametric shortest paths over all~$\lambda\in [0,1]$ that exist for a nontrivial interval of parameters, each corresponding to a linear function in a maximal sub-interval of $[0,1]$. Using this data structure, a shortest path query takes~$\Theta(\log k)$ time.

This paper presents a novel space-filling polyhedron (SFPH), here named the Josehedron, derived from the extremal points of the Fischer-Koch S triply periodic minimal surface (TPMS). The Josehedron is a plesiohedron with 12 faces (4 isosceles triangles and 8 mirror-symmetric quadrilaterals), 12 vertices, and 22 edges. It tiles three-dimensional space with 12 instances per cubic unit cell in 6 distinct orientations. The generating point set exhibits a remarkable connection to the pentagonal Cairo tiling when projected onto any coordinate plane. Several additional geometric properties are described, including integer vertex coordinates, interwoven labyrinths, and chiral symmetry between the polyhedra obtained from the combined minima and maxima of the function. Finally, the paper presents a general method for finding novel SFPHs based on any periodic function, TPMS, or other functions. The described method is applied to a selection of TPMS, and 7 additional, previously undocumented SFPH are shown in the Appendix.

A linkage $\mathcal{L}$ consists of a graph $G=(V,E)$ and an edge-length function $\ell$. Deciding whether $\mathcal{L}$ can be realized as a planar straight-line embedding in $\mathbb{R}^2$ with edge length $\ell(e)$ for all $e \in E$ is $\exists\mathbb{R}$-complete [Abel et al., JoCG'25], even if $\ell \equiv 1$, but a considerable part of $\mathcal{L}$ is rigid. In this paper, we study the computational complexity of the realization question for structurally simpler, less rigid linkages inside an open polygonal domain $P$, where the placement of some vertices may be specified in the input. We show XP-membership and W[1]-hardness with respect to the size of $G$, even if $\ell \equiv 1$ and no vertex positions are prescribed. Furthermore, we consider the case where $G$ is a path with prescribed start and end position and $\ell \equiv 1$. Despite the absence of any rigid components, we obtain NP-hardness in general, and provide a linear-time algorithm for arbitrary $\ell$ if $G$ has only three edges and $P$ is convex.

The Solow--Polasky diversity indicator (or magnitude) is a classical measure of diversity based on pairwise distances. It has applications in ecology, conservation planning, and, more recently, in algorithmic subset selection and diversity optimization. In this note, we investigate the computational complexity of selecting a subset of fixed cardinality from a finite set so as to maximize the Solow--Polasky diversity value. We prove that this problem is NP-hard in general metric spaces. The reduction is from the classical Independent Set problem and uses a simple metric construction containing only two non-zero distance values. Importantly, the hardness result holds for every fixed kernel parameter $\theta_0>0$; equivalently, by rescaling the metric, one may fix the parameter to $1$ without loss of generality. A central point is that this is not a boilerplate reduction: because the Solow--Polasky objective is defined through matrix inversion, it is a nontrivial nonlinear function of the distances. Accordingly, the proof requires a dedicated strict-monotonicity argument for the specific family of distance matrices arising in the reduction; this strict monotonicity is established here for that family, but it is not assumed to hold in full generality. We also explain how the proof connects to continuity and monotonicity considerations for diversity indicators.

We provide a simple algorithm for computing a balanced separator for a set of segments that is $c$-packed, showing that the separator cuts only $O(c)$ segments. While the result was known before, arguably our proof is simpler.

We present a method for generating orthogonal quadrilateral meshes subject to user-defined feature alignment and sizing constraints. The approach relies on computing integrable orthogonal frame fields, whose symmetries are implicitly represented using orthogonally decomposable (odeco) tensors. We extend the existing 2D odeco integrability formulation to the 3D setting, and define the useful energies in a finite element approach. Our frame fields are shear-free (orthogonal) by construction, and we provide terms to minimize area and/or stretch distortion. The optimization naturally creates and places singularities to achieve integrability, obviating the need for user placement or greedy iterative methods. We validate the method on both smooth surfaces and feature-rich CAD models. Compared to previous works on integrable frame fields, we offer better performance in the presence of mesh sizing constraints and achieve lower distortion metrics.