Piecewise linear surface-to-plane maps, or parametrizations, are ubiquitous in computer graphics, geometry processing, mechanical engineering, and scientiic visualization. Depending on the applications, the maps are required to exhibit diferent properties, most commonly, low distortion, local injectivity, and global bijectivity.
The last two properties are challenging to guarantee for discrete maps. Most algorithms with guarantees use Tutte embedding as a component. Tutte embedding is a construction that is guaranteed to create bijective mappings under minimal assumptions, if both domains are simply connected and the target planar domain is convex. However, the guarantee only holds if the computation is performed in arbitrary precision rather than loating point arithmetic, as it is commonly done. Failure due to loating point approximation is not as uncommon as one would assume, as the algorithm is likely to create an extreme variation of scale and aspect ratios in complex mapping cases. To quantitatively evaluate this issue, we computed Tutte embeddings on 2718 models (all the genus 0 models from Thingi10k [Zhou and Jacobson 2016]) using double precision, and observed 80 failures. To the best of our knowledge, this problem has not been addressed before in the literature.
This rate of failure is problematic for batch processing large geometrical collections (for example for processing geometric deep learning datasets) or when the embedding has to be computed many times (for example in cross-parametrization [Kraevoy et al. 2003;Schreiner et al. 2004]). In these scenarios, a failure rate of 2.9% may not be tolerable, since it is not realistic to manually ix hundreds of problematic cases, and if failure happens on large meshes with millions of triangles it might not even be possible to ix them by hand.
A simple solution to this problem is the use of multi-precision (or rational) arithmetic [Granlund 2018]: if enough bits are used to represent the mantissa and exponent of the loating point representation, Tutte embedding will succeed, since the solution of a linear system can be computed exactly. However, the result in high precision is not directly usable by downstream applications, and requires to be rounded (or łsnappedž [Halperin and Packer 2002]) to loating point coordinates. This is a surprisingly challenging problem for which, to the best of our knowledge, no solution applicable to our setting exists (Section 3.1).
Instead, we propose a progressive algorithm to directly generate an embedding using loating point coordinates. We start from an initial, possibly invalid, loating point planar parametrization, and 32:2 • Hanxiao Shen, Zhongshi Jiang, Denis Zorin, and Daniele Panozzo we make it valid by collapsing all lipped and degenerate parts of the (possibly invalid) embedding produced, e.g., by a loatingpoint Tutte algorithm. We re-insert one vertex of the original mesh at a time, preserving the validity of the map at every step. This approach is inspired by [Schreiner et al. 2004], which proposes a progressive algorithm for computing cross-parametrizations based on progressive meshes [Hoppe 1996]. Our algorithm difers since we do not know a valid position for the inserted vertices, and we thus have to compute it as the vertices are added back. We provide a formal proof of correctness of our method in arbitrary precision (obtaining the same formal guarantees as Tutte embedding), and we practically demonstrate its superior robustness by parametrizing a large collection of 10k models.
Using our new embedding method and the matchmaker algorithm [Kraevoy and Shefer 2004;Kraevoy et al. 2003] as a foundation, we develop an algorithm for mapping between multiply-connected domains with arbitrary constraints, supporting fully general selfoverlapping domains as the target. We experimentally show that our algorithm is very robust, producing valid and distortion-optimized maps even for challenging cases where the original matchmaker algorithm fails due to numerical problems. We demonstrate the practical utility of our algorithm for UV mapping and quadrangulation applications.
To foster replicability of results and to maximize the practical impact of our algorithm, we also attach a reference implementation. https://github.com/hankstag/progressive_embedding
Fary's theorem [Fáry 1948] states that any planar graph can be embedded in the plane with straight edges. Tutte [Tutte 1963] extends this result to the case of ixed convex boundary with a spring analogue, and [Floater 1997] established its connection to the parameterization methods in the geometry processing community, and extend Tutte's uniform weight to arbitrary positive ones. In both cases, the problem is reduced to solving a linear system of equations and the resulting embeddings' minimum area might even be negative exponential with respect to the number of vertices. There has been active efort in the graph drawing community to address these issues, by bounding the total area when drawing on integer grids, or equivalently, controlling the minimum resolution [Chambers et al. 2011] under ixed diameter. Most notably, [Schnyder 1990] shows an algorithm to embed a planar graph onto integer grids inside a triangle region, and [Chambers et al. 2011] proves an upper polynomial bound on the area while keeping a speciied convex boundary shape: the proof is constructive and may (potentially) be used as a basis for a practical algorithm. In all cases, rounding problems will afect these algorithms as the size of the graph grows (Section 3.1).
Orbifold Tutte Embedding. Multiple extensions of Tutte's theorem to map surfaces to diferent co-domains have been proposed. In particular, the theorem has been extended to map surfaces to a Euclidean orbifold [Aigerman and Lipman 2015], to a hyperbolic orbifold [Aigerman and Lipman 2016], and to a spherical orbifold [Aigerman et al. 2017]. All three methods support hard positional constraints and ensure the generation of a bijective map between the surface and the orbifold in ininite precision arithmetic. These methods also sufer from similar numerical issue as Tutte's, and extending our algorithm to orbifold embeddings is an interesting direction for future work.
The well-known progressive meshes algorithm [Hoppe 1996;Sander et al. 2001] shows how a triangle mesh can be simpliied by collapsing one edge at a time, and reconstructed applying the inverse topological operations in the inverse order. This scheme has been introduced as an eicient way to store, transmit, and render large meshes, where the per-vertex properties of the removed vertex are stored together with the information required to insert them back. This work has been later applied to compute inter-surface mappings [Schreiner et al. 2004], by jointly simplifying two meshes into a common base mesh, then starting from optimizing their isometric distortion while reinserting the vertices in the base mesh.
We use the same idea to eliminate problematic regions of an existing embedding (either lipped, or with a high distortion), and then reinserting one vertex at a time, while preserving the quality of the triangulation. Diferently from progressive meshes, in our case we do not have geometrical information available that could help us decide where the vertex should be inserted to obtain a valid embedding.
In this section, we focus on the recent works closely related to generating distortion-minimizing discrete locally injective and globally bijective discrete maps, and we refer to [Floater and Hormann 2005;Hormann et al. 2007;Shefer et al. 2006] for a comprehensive treatment of earlier parametrization methods without these properties.
A discrete locally injective map requires that triangles maintain their orientation (i.e. they do not lip) and if the sum of (unsigned) triangle angles around each internal vertex is precisely 2π [Weber and Zorin 2014]. Three main families of methods have been proposed to deal with this challenging constraint: barrier, convexiication, and hybrid algorithms.
Barrier Algorithms. Barrier algorithms require a valid initial solution, and then optimize its quality without leaving the feasible space. The key idea is to adopt quality metrics diverging to ininity when triangles become degenerate, thus inhibiting lips. Popular choices strive to preserve angles [Degener et al. 2003;Hormann and Greiner 2000] or lengths [Aigerman et al. 2014;Poranne and Lipman 2014;Sander et al. 2001;Smith and Schaefer 2015;Sorkine et al. 2002]. Alternatively, a barrier functions can be added to existing energies to enforce local injectivity [Schüller et al. 2013]. These non-linear energies are diicult to minimize, stemming a series of methods speciically targeting this problem. They include coordinate descent [Hormann and Greiner 2000;Labsik et al. 2000], parallel gradient descent [Fu et al. 2015], Anderson Acceleration [Peng et al. 2018], as well as other quasi-newton approaches [Claici et al. 2017;Kovalsky et al. 2016;Liu et al. 2018;Rabinovich et al. 2017;Shtengel et al. 2017;Smith and Schaefer 2015;Zhu et al. 2018].
All these methods support hard-constraints if they are already satisied in the initial map, which is the key idea used in MatchMaker
• 32:3 [Kraevoy et al. 2003]. Our progressive embedding can be used to robustly generate the initial map, that can then be improved by any of the previous techniques (Section 5).
Projection Algorithms. An essential component of these methods is a convexiied form of the injectivity constraints [Kovalsky et al. 2015;Lipman 2012]. While these methods naturally support hard injectivity constraints, they might fail to ind a feasible solution, with no output generated. The only known way to guarantee that a feasible solution exists is to formulate the convexiied constraints using a reference frame derived from a valid (although potentially very high distortion) solution.
Hybrid Algorithms. Hybrid algorithms are an interesting mix between these two approaches [Fu and Liu 2016;Poranne et al. 2017]. The initial guess is produced by separating all triangles and isometrically rotating them into the UV space. A barrier method is then used to prevent them from lipping, while trying to seal the seams. This approach might fail to seal all the seams, not producing a valid map.
Globally Bijective Maps. For simply connected domains, bijective maps are locally injective maps whose boundary does not intersect. All embeddings described in Section 2.1 satisfy this property. These methods have been extended to non-convex, self-overlapping polygons [Weber and Zorin 2014] and polyhedrons [Campen et al. 2016], but they still require a ixed boundary. Few methods can produce bijective maps while letting the boundary free, relying on either collision detection [Smith and Schaefer 2015] or scafolding elements [Gotsman and Surazhsky 2001;Jiang et al. 2017;Müller et al. 2015;Zhang et al. 2005]. All free boundary methods require a starting point: our algorithm can be used to generate it, enabling these algorithms to create bijective maps with hard constraints (Section 4).
Hard Positional Constraints and Reinement. Matchmaker [Kraevoy et al. 2003] introduced hard positional constraints for texture mapping applications. The algorithm uses a two-step approach, irst generating a valid map, and then optimizing its geometrical quality. The method is one of the few using reinement to guarantee the existence of feasible solutions. The method has been extended by adding an intermediate warping stage to align the constraints in [Lee et al. 2008]. We show in section 4 how our embedding can be used withtin Matchmaker to increase its robustness, and we also show how to extend Matchmaker to support self-overlapping polygonal target domains.
Cross-Parametrization. Cross-parametrization, i.e. the computation of a map between two surfaces, is another problem that often relies on planar embeddings. [Schreiner et al. 2004] and [Kraevoy and Shefer 2004] proposed the irst provably guaranteed solutions to compute maps between surfaces, by reducing the problem to mapping both surfaces to a common subdomain by either using Tutte's embedding or a simpliication approach. A similar construction that cuts open the surface into a single topological disc has been proposed in [Aigerman et al. 2014], and extended to allow even the optimization of the seams positions in [Aigerman et al. 2015]. Floating point rounding errors have not been considered in any of these works, which are more prone to fail as the resolution of the mesh increase or whenever the user-provided constraints introduce a high distortion (Section 5).
Global Parametrization. Field-aligned parametrization methods [Bommes et al. 2009] strive to compute a locally injective map [Bommes et al. 2013] whose gradient is aligned with a user-provided directional ield. We refer an interested reader to [Bommes et al. 2012] for a comprehensive overview of these techniques. Our embedding algorithm can be used to compute parametrizations to a target self-overlapping polygon, enabling to robustly generate these parametrization if a valid boundary polygon is provided (Section 5).
We discuss in detail when Tutte embedding implemented in loating point may fail, and also show that straightforward solutions with of-the-shelves geometry processing tools do not solve these issues.
Tutte Embedding implemented in Floating Points. We use the implementation of Tutte embedding in libigl [Jacobson et al. 2016], and apply it to all the 2718 genus 0 models of the Thingi10k dataset [Zhou and Jacobson 2016], after cleaning them up and improving their quality using TetWild [Hu et al. 2018], to ensure that no degenerate triangles are present. We also ensure that the meshes are 3-connected by reining them locally. For every model, we randomly pick and delete a triangle, and map the resulting boundary to an equilateral triangle. We compute the Tutte embedding, and check for lips using CGAL's exact loating point predicates [Brönnimann et al. 2018]. The check fails for 80 models, due to the numerical errors introduced in the mapping. In retrospect, this is not surprising since it is well-known that Tutte planar drawing may admit exponential area when drawing on integer grids. Two problematic cases are shown in Figure 2, where the embedding introduces a large variation of scale, and the lip occurs on triangles with small areas.
Multi-Precision Tutte Embedding with Snap Rounding. A straightforward way to address this problem is to increase the number of bits used in the loating point representation. We double the number of bits using the library MPFR [Fousse et al. 2007], which is directly integrated into Eigen, and can thus be used with the Tutte embedding in libigl with minimal code changes. With this setup, all the problematic cases are solved. However, the runtime is increased by around one order of magnitude, and, most importantly, the results generated cannot be rounded back to loating point since trivial rounding introduces lips. Snap rounding [Packer 2018] could be used to avoid them, but it will collapse possibly large regions of the mesh: 6.3% of the vertices of the model shown in the bottom left of Figure 8 are collapsed when using a snap rounding resolution of 10 -16 times the diagonal of the bounding box of the embedding.
Multi-Precision Tutte Embedding with Quality Optimization. The problem with rounding to loats is induced by the small triangles (and correspondingly small edges) which leads to lips after snapping. A possible way to address this issue is to use a mesh optimization algorithm, using multi-precision representation, before rounding to loats. We tested two approaches: (1) SLIM [Rabinovich et al. 2017] adapted to run in multiprecision, and (2) minimizing the symmetric Dirichlet energy by moving one vertex at a time using coordinate descent [Hormann and Greiner 2000]. The irst approach is prohibitively slow, due to the linear solve in high precision and the very small steps due to the elements with almost zero area. The second one succeeds on 57 models, but still fails on 23, even after 24 hours of running time.
Our approach draws on the ideas of progressive meshes [Hoppe 1996] and inter-surface mappings [Schreiner et al. 2004], which are, in turn, closely related to theoretical ideas from PL topology (e.g., [Hudson and Shaneson 1969]).
Our algorithm does not require Tutte embedding and can be used to construct an embedding from scratch or from a random initial mapping (Figure 3). It can be accelerated by using an existing, possibly invalid, embedding as a starting point. A triangle is invalid if its signed area is negative, or if its quality measure is below a threshold (we use the symmetric Dirichlet energy [Smith and Schaefer 2015] with respect to a canonical equilateral triangle whose area is the area of the target boundary polygon divided by the number of triangles and mark invalid if it is above τ = 1e20). Using a quality measure in addition to signed area is important, since triangles with small, positive areas might cause numerical problems during the vertex insertion (Phase 2 below).
Starting from an invalid embedding, our algorithm (1) performs edge collapses until the simpliied mesh has no invalid triangles (Algorithm 1), and (2) progressively inserts back each vertex in the same order (Algorithm 2), with feasibility of insertion ensured at each step.
Stage 1: Simpliication. At this stage, we iteratively ind an interior edge that can be collapsed until all invalid elements are removed from the initial embedding, or a single interior vertex is left (Algorithm 1). A theorem in [Mijatović 2003] and our Theorem A.8 guarantees that a sequence of collapses reducing the mesh to a mesh Fig. 5. The two admissible insertion positions from Lemma A.11. The dark region on the let shows the valid positions for v m while fixing v 0 . The right case is the opposite. Our algorithm opts for the let case for stability, since the calculation of the valid sector in the right case involves intersection of the prolonged edge (dashed lines) and the 1-ring neighbors. We pick the valid sector as the one that has an inner angle sum smaller than π . with a single interior vertex can always be found; as the boundary embedding is convex, we can also always ind a position for this vertex to create a valid embedding for the fully simpliied mesh. The simpliication algorithm starts by tagging all invalid triangles (Line 1), and attempts to collapse all their edges. If this procedure is successful in eliminating all invalid triangles, the algorithm terminates, otherwise all the triangles adjacent to tagged triangles are tagged (Line 12), and their edges collapsed. Note that we only allow edge collapses on internal edges to avoid changes to the boundary. This algorithm is guaranteed to terminate, since in the worst case it will tag the entire mesh and Theorem A.8 ensures that at least one edge will be collapsible. If only one internal vertex is left (Line 4), we move it to the barycenter of the boundary vertices (which is inside the convex boundary by construction, but might fail in degenerate cases, as discussed in Section 6). Once the algorithm terminates, the resulting simpliied mesh has no inverted triangles and by Proposition A.9, also has sums of triangle angles at each vertex equal to 2π .
Stage 2: Insertion. Starting from the valid embedding computed after Stage 1, we perform a sequence of splits reverting the collapses, while maintaining embedding validity at every step (Algorithm 2). Lemma A.11 ensures that this is always possible (in ininite precision), as the area of each triangle after insertion is always positive, whenever the newly inserted points lie in the valid sector (Figure 5). The algorithm irst computes candidate directions in the valid sector, then performs a lip-avoiding line search [Smith and Schaefer 2015] (Line 3) to ind candidate positions P(v m ) for the newly inserted vertex v m , such that the 1-ring neighborhoods are valid. In our experiments, we use step length α = 0.8, and cap the number of line search iterations to 75. If at least one candidate is found, the split is performed using the candidate resulting in a mesh minimizing the error measured as the maximum of the 1-ring energy. Since a candidate position always exist in ininite precision (Lemma A.11), the only possible cause for not inding it is a lack of representation power in the loating point representation. We thus improve the quality of the mesh (Line 11) until a candidate is found.
This algorithm may still fail to ind a candidate in degenerate conigurations. However, we experimentally found that the constrained mesh smoothing is very efective at ameliorating this issue, keeping the mesh quality suiciently high during the insertion to allow split operations to succeed (Figure 6). In all our experiments we only found one failure case, where the prescribed target boundary is a numerically degenerate triangle (Section 6). All our other experiments, even on a large data set and with complex boundary conditions (Section 5) were successful.
Local Smoothing. To improve the quality of the map in the insertion step (Algorithm 2, line 8 and 11), we minimize the symmetric Dirichlet energy [Smith and Schaefer 2015], optimizing one vertex position at a time using Newton iterations, similarly to [Fu et al. 2015;Hormann and Greiner 2000;Labsik et al. 2000]. We favor this local approach since it is more robust to low quality elements, which would otherwise badly afect both the numerical stability and the step size of global optimization methods. Since our goal is to improve the minimal quality of the mesh, we minimize the symmetric Dirichlet energy only for the invalid triangles (using reference shape as an equilateral triangle with area equal to the average of triangles in the 2D domain), and use the scafold energy [Jiang et al. 2017] (i.e. we use the element itself as the reference triangle for the symmetric Dirichlet energy) for the valid ones, which allows them to move more freely. In our experiments, the local smoothing is performed for 10 iterations after every insertion step. If a valid insertion candidate cannot be found, we keep improving the quality with batches of 50 smoothing iterations until a candidate is found (Algorithm 2 line 11). Furthermore, at each smoothing phase, we perform a greedy coloring of the edge graph [Kucera 1991], and the vertices inside each color are optimized in parallel.
The computation of locally injective maps is important in geometry processing (Section 2), with the majority of the methods focusing on eicient and scalable quality optimization. However, few method guarantees positional constraints: the notable MatchMaker algorithm [Kraevoy et al. 2003] reduces the problem to a number of convex planar embeddings, which are computed with Tutte's algorithm. However, as we observe in some cases (Section 5.2), such embeddings can be numerically challenging. Replacing Tutte embedding with progressive embedding enables matchmaker to robustly compute maps with very challenging conigurations of constraints. In this Section, we describe an extension of [Kraevoy et al. 2003] that (1) makes use of progressive embedding to increase robustness, and (2) supports weakly self-overlapping polygons as co-domains [Weber and Zorin 2014].
Overview. Combining our progressive embedding algorithm and the matchmaker algorithm, we describe an algorithm for solving the following problem: Given a simply-connected 3d mesh, equipped with with a set of user-deined hard positional constraints at vertices, compute a valid piecewise-linear parametrization, such that (1) the map is valid in the following sense: there are no lipped triangles, and for each vertex, the map restricted to the one ring of triangle of that vertex is bijective, unless it is a singular boundary vertex, as deined below and (2) the parametrization bitwise exactly satisies the user-deined positional constraints. We tackle this in three steps: we decompose the target domain into convex polygonal subdomains, match these domains to the subdomains of the source domain, compute an initial bijective map by stitching progressive embeddings for each subdomain, and inal globally optimize the mapping distortion.
User input. We distinguish between two cases, chosen by the user (1) the required map is a global embedding, (2) the map is an immersion. For the irst case, the constraint speciication is more lexible: the user only has to provide a set of point or line constraints. For the second case, the target domain is not a subset of the plane, but rather, an everywhere lat surface with overlaps. We require the user to prescribe constraints for the whole boundary of the polygon to deine the target domain unambiguously (some parts may be marked as movable, but an initial position is needed) and to provide a path connecting each point or line constraint to the boundary, which allows to deine its location on the target surface implied by the boundary speciication. In this case, target boundary polygon has to be weakly self-overlapping [Weber and Zorin 2014], otherwise, the map does not exist.
Phase 1: Subdivision of the Target Domain. In the global embedding case, the target domain is generated by triangulating the bounding box of the input with Triangle [Shewchuk 1996]. In the second case, the self-overlapping domain is triangulated using a modiication of the Shor-Van Wyck algorithm [Shor and Van Wyk 1992], described in [Weber and Zorin 2014]. In both cases we ensure that the hard positional constraints are vertices of the triangulation. We then merge triangles into convex polygons in a greedy manner, by dropping edges if the resulting subdomains are convex. While this merging step is optional (the algorithm works also without it) it reduces the number of subdomains and vertices, making the next steps more eicient. In the case when an immersion is computed, by construction, it will be an embedding on subdomains computed starting from the Shor-Van Wyck triangulation.
Phase 2: Path Tracing on the Original Mesh. After the target domain is subdivided into convex polygons, we match this decomposition to the input mesh. The goal is to ind non-intersecting paths Fig. 7. Starting from a triangulation generated from only boundary segments and internal constraint points (let). Instead of treating triangles as sub-domains as in [Kraevoy et al. 2003], we merge triangles to convex polygons (middle). Then we find paths (bold, right) connecting constraint points to the boundary without new cycles, and prioritize their tracing.
connecting each of these pairs in 3D mesh, subdividing the whole mesh into same number of patches.
At the tracing stage, we perform a reordering of the paths to make sure that no previously traced path will block the future ones (Figure 7). In the case of embedding, the algorithm inds paths that connect all positional constraints to the boundary without creating additional loops (than the existing boundary). In practice, these paths are found by dropping one segment on the boundary, and then grow the minimum spanning tree (over the edges of the polygonal mesh) from the incomplete boundary loop. We irst trace the paths on the minimum spanning tree and then the remaining ones, connecting the boundary to the constraints. The correctness of this procedure can be found in [Praun et al. 2001].
For the tracing of each path on the surface, we follow [Kraevoy et al. 2003] to ind the shortest path connecting two endpoints, and add Steiner points on the edges if no path, not intersecting other paths, can be found.
Phase 3: Bijective Mapping. After establishing a correspondence between each patch of the 3D mesh and a convex polygon in the target domain, we can irst subdivide the 2D paths in the target domain to match the number of vertices on the corresponding 3D path to obtain the one-to-one correspondence between them. We observe that up to this point, the algorithm is largely combinatorial (while some vertices are inserted on edges, their geometric position is trivially determined and is very unlikely to result in numerical problems; none were observed in our experiments). At this point, the map is deined for boundaries of the subdomains corresponding to the convex subdomains in the target. Next, we extend the map to the interior of each region using our progressive embedding algorithm (Section 3.2). Notice that when the mesh patch is not 3-connected, we need to split the edges with two endpoints on the boundary.
Phase 4: Quality Optimization. The map obtained in the previous steps is valid, according to the deinition of the weakly selfoverlapping map [Weber and Zorin 2014]. Therefore, its quality can be optimized using any locally injective map improvement algorithm (Section 2). We opt for [Rabinovich et al. 2017], since it is eicient for large models and the implementation is readily available [Jacobson et al. 2016]. The implementation is modiied to support hard positional constraints, by eliminating the corresponding variables.
We implemented our algorithm in C++, using Eigen [Guennebaud et al. 2010] for linear algebra, and libigl [Jacobson et al. 2016] 1. Statistics of the input and output meshes in the planar embedding test (Section 5.1). From let to right: Name of the dataset, number of vertices, number of faces, number of invalid elements (positive area, but with energy above 1e20) ater Tute embedding, number of flipped elements ater Tute embedding, progressive embedding (Section 3) running time in seconds. geometry processing and visualization. The reference source code, the data used, and the scripts to reproduce the results are attached in the additional material. The timings and statistics for the datasets shown in the paper are summarized in Table 1 andTable 2. We irst present results computed using only our progressive embedding (Section 5.1), and then demonstrate the generation of low distortion, locally bijective maps created with our extension of MatchMaker (Section 5.2).
Planar Embedding for the Thingi10k Dataset [Zhou and Jacobson 2016]. By computing Tutte's embedding for the genus-zero models in 2718 surface mesh models on a triangle boundary, we observed there are 80 cases where the generated parametrization has lipped elements due to loating point rounding errors. Using our progressive strategy, we are able to ix all failed cases. A selection of the parametrization results are shown in Figure 2.
Integration with OptCuts [Li et al. 2018]. OptCuts is a joint optimization method to create UV seam from a 3D model, balancing seam length and parameterization quality. It relies on a valid initialization, which for genus 0 model, is compute through randomly cutting two adjacent edges as seams, then latten it on the plane with Tutte embedding. In Figure 8, we show two examples where this initialization fail. Both models can be processed if progressive Fig. 8. Three UV maps generated by OptCuts [Li et al. 2018] using an initial embedding created by our algorithm. OptCuts fails to process both models if Tute embedding is used instead.
embedding is used instead of Tutte embedding, allowing OptCuts to proceed and optimize the UV map.
Mapping an Hele-Shaw Polygon to a Square. Hele-Shaw low is a two-dimensional Stokes low of mixing liquids between two parallel lat surfaces separated by a small gap. In Figure 1, we show an example mesh generated using the Hele-Shaw simulation proposed in [Segall et al. 2016]. One way to compute a bijective map of the interior of the polygon between diferent frames is a crossparameterization using a square as the common domain, with no internal constraints. Tutte embedding fails in this case, introducing 46 lipped faces (Figure 1, left), while progressive embedding produces a valid map with lower distortion (Figure 1, right). [Kraevoy et al. 2003] Self-Overlapping Locally-Injective Maps. By introducing Shor Van Wyck algorithm into the matchmaker pipeline, we are able to mapping a surface mesh with disk topology to self-overlaping boundaries as in [Weber and Zorin 2014]. Similarly to [Weber and Zorin 2014] our algorithm can generate locally-injective, self-overlapping parametrizations (Figure 9), which are commonly used by quadrangulation algorithms [Bommes et al. 2012].
Comparison with [Kovalsky et al. 2015]. We parametrized the global parametrization benchmark introduced in [Myles et al. 2014], using the seams in the obj iles, and ixing in random positions 3 random points of each mesh. This is a challenging task, since the random constraints introduce a large distortion. Our method succeeded on all 102 models: a selection of the most challenging ones is shown in Figure 10. We also run the same experiment using the most recent projection method [Kovalsky et al. 2015] (which is one of the few methods that supports similar constraints without requiring a fully speciied target domain), using LSCM [Lévy et al. 2002] as an initial guess. The method failed on 28 models over 102 (27%). We show three failed cases using their method with lipped elements in the output, and the quality is considerably lower than our approach, as shown in Figure 11. Note that this is a comparison that favours our method, since we are allowed to remesh the map, while [Kovalsky et al. 2015] preserves the original connectivity.
Stress Test. To further evaluate the robustness and applicability of our algorithm, we performed an additional stress test, by parametrizing the 102 models of [Myles et al. 2014] into a planar space illing curve, and adding 3 random positional constraints. These experiments push the algorithm to the limit: MatchMaker fails on 5 if Tutte embedding is used, while it succeeds in all cases, producing bijective maps exactly satisfying the hard positional constraints, with progressive embedding (Figure 12).
We introduced a robust algorithm to compute planar embeddings, and demonstrated its practical utility in common geometry processing tasks. Our algorithm is provably correct in ininite precision and is designed to work robustly with loating point coordinates: unfortunately we cannot guarantee that an output is produced in the latter case since a solution of the local point placement problem might not exist. Consider the example in Figure 13: the bounding box of the triangle has short sides (the diference between the loating point coordinate representation is only in the least signiicant bit of the mantissa). Assume that our algorithm needs to split of a vertex from the vertex with numerically lat angle A, placing the resulting point in the interior. In this situation, our algorithm fails, since the average of the coordinates (in loating point) of the boundary triangle does not lie inside the triangle due to numerical rounding.
Except for this extreme case, we have not observed any other failure cases for our algorithm, which produced robustly thousands of embeddings, and, when paired with matchmaker, enables the robust generation of constrained locally injective maps. /2, 1), and(b/2, 1) resp., where h = 2 -53 (The illustration is not to scale.)
The proof of the existence of the collapse sequence for two-dimensional manifold meshes can be found, e.g., in [Mijatović 2003], where it is derived from the shellability of two-dimensional manifold meshes homeomorphic to a disk (i.e., the possibility of removing triangles one-by-one, keeping the topology of the remaining part of the mesh unchanged), and make use of a speciic composition of Pachner moves equivalent to edge collapse. We present a diferent proof, based on proving the existence of a collapsible edge, which is aligned with the structure of our algorithm and helps us to show the existence of the inverse vertex split sequence.
We assume that the input mesh connectivity (V 0 , F 0 ) is manifold, i.e., each edge is shared by no more than two triangles, and the triangles incident at a vertex can be arranged in a sequence so that two sequential triangles share an edge. For interior vertices, the sequence is circular, i.e. the irst and last triangles also share an edge. With the topology of a 2D disc, the graphs of edges of such meshes are planar i.e., can be embedded in the plane, with positions P 0 = {p i ∈ R 2 } assigned to vertices v i . By the Fáry's theorem, [Fáry 1948] there is a straight-edge embedding of this graph in the plane with non-intersecting edges (Fáry embedding). In subsequent lemmas, we use geometric images of vertices and edges under this embedding. Only the existence of this embedding, but not the speciic construction, is used to prove the existence of the collapse sequence.
The following sequence of lemmas focuses on the one-ring neighborhood of an interior vertex, and shows that at least one of the adjacent edges satisies the link condition. This observation further leads to Theorem A.8: a valid sequence of edge collapses can be used to reduce the mesh to a mesh with a single interior vertex. A vertex of the mesh is interior if it does not lie on the boundary, and an edge is interior if its two endpoints are interior vertices.
Deinition A.1. An interior edge v i v j satisies the link condition if |N i ∩ N j | = 2, where N i is the set of the adjacent vertices of v i .
Lemma A.2. Let v 0 be an interior vertex of degree d (Figure 14). We enumerate its neighbors counterclockwise around the vertex (using Fáry embedding), denoting them v 1 , v 2 , . . . , v d . Assume v 0 v 1 violates the link condition, i.e.,
Proof. Without the loss of generality, consider ∆v 0 v 1 v k orients counterclockwise. Consider the segment v 0 v i , for 1 < i < k. The half-line starting at v 0 and containing this segment is between halflines containing v 1 and v k , because the vertices were numbered counterclockwise. Therefore the half-line contains a point in the interior of ∆v 0 v 1 v k , by continuity. If v i is outside or on the boundary of ∆v 0 v 1 v k then the half-line connects an interior and non-interior point diferent from v 0 , and intersects v 1 v k . which contradicts the assumption on the embedding. Thus, v i is in the interior of ∆v 0 v 1 v k .
32:12 • Hanxiao Shen, Zhongshi Jiang, Denis Zorin, and Daniele Panozzo interior vertex forms a complete boundary loop. As we assume the mesh to be simply connected, then there is only one boundary loop. So the whole boundary has to coincide with the link of any interior vertex, from which it follows that there is only one. □
We remark that when there is only one interior vertex left, if the boundary vertices v i , i ≥ 1, are assigned positions p i , so that they form a star-shaped simple polygon, there is a position (within the interior of the kernel of the boundary) p 0 for the remaining interior vertex v 0 that results in a valid straight-edge embedding.
As a result of sequentially collapsing edges, we obtain a sequence (V i , F i , C i ), i = 1, 2, . . . , k with the following properties:
) is a valid triangulation, boundary vertices are the same for all V i , and C i is a valid collapse.
Proposition A.9. Suppose the vertices v i of the disk-topology manifold mesh (V , F ) are assigned parametric positions p i in the plane, and the map is bijective on the boundary so that the triangles all have positive orientation. Then the sum of the angles of triangles incident at an interior vertex is 2π .
Proof. Assign, e.g., unit length to all edges of the mesh; this associates a a surface M with the mesh, with each combinatorial triangle corresponding to an equilateral triangle. Then the positions p i deine a PL map from M to the plane. By Theorem 1 from [Lipman 2014], this map is globally bijective; the statement of the proposition immediately follows. □
Deinition A.10. Let (V , F ) be a mesh with a valid straight-edge embedding in the plane given by vertex positions P. Consider a closed fan of triangles F (v 0 ; v 1 . . . v d+1 ) centered at an interior vertex v 0 , with v d +1 = v 1 , and an open sub-fan F (v 0 ; v 1 . . . v k ). Vertex split introduces a new vertex v m with a position p m , F (v 0 ; v 1 . . . v k ), replaces it with a fan F (v m ; v 1 . . . v k ), and adds triangles ∆v 0 v 1 v m and ∆v 0 v m v k . We denote such a split S by (v m , p m , F (v 0 ; v 1 . . . v k )).
the connectivity of the mesh obtained by applying the split is identical to the mesh that the collapse was applied to.
The following lemma establishes that we can perform a split reversing any collapse while maintaining the validity of the embedding, if the initial embedding is valid.
Lemma A.11. Consider a fan of triangles F (v 0 ; v 1 . . . v k ) (Figure 5), with positive signed areas {A(∆v 0 v i-1 v i )} 1<i ≤k and with angles of triangles incident at v 0 summing up to 2π . The kernel of the fan has a non-empty interior. Then a new vertex position p m corresponding to a new vertex v m located in the interior of the fan, can be split of p 0 , so that min 1<i ≤k A(∆v
and angles of triangles in both resulting fans at v 0 and v 1 sum up to 2π .
Proof. Deine a function f (p x ) = min i A(∆v x v i-1 v i ). This is a continuous function of the coordinates p x of the point v x . Because f (p 0 ) > 0, there is a disk B(p 0 , ε), of radius ε > 0, such that f (p x ) > 0 for any p x ∈ B(p 0 , ε), i.e. for all i, A(∆v m v i-1 v i ) > 0, if we pick p m inside B(v 0 , ε). Suppose we initially place p m at p 0 , with new triangles added as a result of the split having zero angles at v 1 and v k . We note that the for each of v 0 and v 1 in this degenerate coniguration the angles of incident triangles sum up to 2π . The angles of triangles also change continuously as functions of vertex position p x , so does their sum. On the other hand, if each triangle remains positively oriented (A(∆v m v i-1 v i ) > 0), then the sum of the angles can only change discretely, and has to be of the form 2πn, n ∈ Z (n-fold cover). We conclude that n has to remain one, as it is one for the initial position.
Consider the intersection C of the half-planes bounded by lines containing p 0 p 1 and p 0 p k (for each segment, we choose the half-line on the side of the interior of the fan). If
The following theorem is a straightforward application of Lemma A.11.
Theorem A.12. Suppose we have a sequence of valid collapses (V i , F i , C i ), i = 0 . . . N -2, where N = |V 0 |, the number of interior vertices in the initial mesh, all (V i , F i ) i = 0 . . . N -1 are 3-connected planar, and the last mesh (V N -1 , F N -1 ) with a single interior vertex has a valid straight-edge embedding in the plane with vertex positions P N -1 . Suppose C i = (v i m , F (v i+1 0 , v i+1 1 . . . v i+1 k )). Then the sequence of vertex splits S i that are inverses of C i , results in a valid straight-edge embedding of (V 0 , F 0 ), given by vertex positions P 0 .
Proof. V N -1 , F N -1 with positions P N -1 is valid by assumption. Each step of vertex split with S i results in a straight-edge embedding of (V i , F i ) by Lemma A.11. By induction, the embedding of (V 0 , F 0 ) obtained by the sequence of splits reverting the sequence of collapses is a straight-edge embedding. MatchMaker [Kraevoy et al. 2003] tackles the problem of planar parametrizaton of a surface with hard positional constraints. The general idea is to decompose the mesh into patches and map them to convex domains using Tutte embedding to generate a bijective parametrization. We briely summarize the pipeline of the original MatchMaker method.
Triangulate Virtual Boundary. Use a conventional unconstrained parametrization method to get a planar mesh M 0 (the dark region in Figure 16(b)). Embed M 0 in a rectangular bounding box, and use the constrained Delaunay triangulation to triangulate the region between the planar mesh and the virtual boundary, to obtain M ′ 1 (Figure 16(a)). Similarly, create a constrained Delaunay triangluation of the same bounding box with hard constraints at prescribed position using [Shewchuk 1996], to obtain the mesh in Figure 16 (a). We call this mesh M 1 .
Matching Patches. For each internal edge of M 1 , trace an edge path between the two corresponding points on M ′ 0 using Dijsktra's shortest path algorithm. This process is performed sequentially, one edge each time. Every new paths traced must meet the following criteria: (1) it does not intersect with any other paths; (2) it does not block necessary future paths. This is achieved using the minimum spanning tree as described in Section 4. Steiner points may be needed to make sure a valid path can always be found (red dot in Figure 17 (c)). As the result of this step, extended planar mesh M ′ 0 is subdivided into patches, and every patch is matched with a triangle in M 1 .
Embedding. At the previous step, the problem is reduced to mapping a patch of a planar mesh to a convex boundary, i.e., a problem solved by Tutte embedding [Tutte 1963]. Additional edge splitting operation is needed to match the boundary vertices number of triangles in M 0 to vertices number on corresponding paths in M ′ 1 , as shown in Figure 16 (b).
Smoothing and Post-processing. The parametrization generated at the embedding stage can be optimized using any technique preventing triangles from changing their orientations.
ACM Trans. Graph., Vol. 38, No. 4, Article 32. Publication date: July 2019.