We introduce a Monte Carlo method for evaluating derivatives of the solution to a partial differential equation (PDE) with respect to problem parameters (such as domain geometry or boundary conditions). Derivatives can be evaluated at arbitrary points without performing a global solve, or constructing a volumetric grid or mesh. The method is hence well-suited to inverse problems with complex geometry, such as PDE-constrained shape optimization. Like other walk on spheres (WoS) algorithms, our method is trivial to parallelize, and is agnostic to boundary representation (meshes, splines, implicit surfaces, etc.), supporting large topological changes. We focus in particular on screened Poisson equations, which model diverse problems from scientific and geometric computing. As in differentiable rendering, we jointly estimate derivatives with respect to all parameters—hence, cost does not grow significantly with parameter count. In practice, even noisy derivative estimates exhibit fast, stable convergence for stochastic gradient-based optimization—as we show via examples from thermal design, shape from diffusion, and computer graphics.

The sphere tracing algorithm provides a fast and high-quality strategy for visualizing surfaces encoded by signed distance functions (SDFs), which have become a centerpiece in a wide range of visual computing algorithms. In this paper we introduce a sphere tracing algorithm for a completely different class of functions, harmonic functions, opening up a whole new set of possibilities. Harmonic functions are found throughout geometric and visual computing, where they arise naturally as the solution to interpolation problems, and in the physical sciences, where they appear as solutions to the Laplace equation. Our algorithm for harmonic functions is similar in spirit to the sphere tracing algorithm for SDFs: by using a conservative lower bound on the distance to the level set, we can take much larger steps than with naïve ray marching. Our key observation is that for harmonic functions such a bound is given by Harnack's inequality. Unlike Lipschitz bounds used in traditional sphere tracing, this Harnack bound requires only the value of the function at a point—we use this bound to develop a sphere tracing algorithm that can also handle jump discontinuities arising in angle-based harmonic functions. We show how this algorithm can be used to directly visualize smooth surfaces reconstructed from point clouds (via Poisson surface reconstruction) or polygon soup (via generalized winding numbers) without performing linear solves or mesh extraction. We also show how it can be used to render nonplanar polygons (including those with holes), and to visualize key objects from mathematics, including knots, links, spherical harmonics, and Riemann surfaces.

This paper develops a shape space framework for collision-aware geometric modeling, where basic geometric operations automatically avoid interpenetration. Shape spaces are a powerful tool for surface modeling, shape analysis, nonrigid motion planning, and animation, but past formulations permit nonphysical intersections. Our framework augments an existing shape space using a repulsive energy such that collision avoidance becomes a first-class property, encoded in the Riemannian metric itself. In turn, tasks like intersection-free shape interpolation or motion extrapolation amount to simply computing geodesic paths via standard numerical algorithms. To make optimization practical, we develop an adaptive collision penalty that prevents mesh self-intersection, and converges to a meaningful limit energy under refinement. The final algorithms apply to any category of shape, and do not require a dataset of examples, training, rigging, nor any other prior information. For instance, to interpolate between two shapes we need only a single pair of meshes with the same connectivity. We evaluate our method on a variety of challenging examples from modeling and animation.

Numerous scientific and engineering applications require solutions to boundary value problems (BVPs) involving elliptic partial differential equations, such as the Laplace or Poisson equations, on geometrically intricate domains. We develop a Monte Carlo method for solving such BVPs with arbitrary first-order linear boundary conditions—Dirichlet, Neumann, and Robin. Our method directly generalizes the walk on stars (WoSt) algorithm, which previously tackled only the first two types of boundary conditions, with a few simple modifications. Unlike conventional numerical methods, WoSt does not need finite element meshing or global solves. Similar to Monte Carlo rendering, it instead computes pointwise solution estimates by simulating random walks along star-shaped regions inside the BVP domain, using efficient ray-intersection and distance queries. To ensure WoSt produces bounded-variance estimates in the presence of Robin boundary conditions, we show that it is sufficient to modify how WoSt selects the size of these star-shaped regions. Our generalized WoSt algorithm reduces estimation error by orders of magnitude relative to alternative grid-free methods such as the walk on boundary algorithm. We also develop bidirectional and boundary value caching strategies to further reduce estimation error. Our algorithm is trivial to parallelize, scales sublinearly with increasing geometric detail, and enables progressive and view-dependent evaluation.

We introduce a method for approximating the signed distance function (SDF) of geometry corrupted by holes, noise, or self-intersections. The method implicitly defines a completed version of the shape, rather than explicitly repairing the given input. Our starting point is a modified version of the heat method for geodesic distance, which diffuses normal vectors rather than a scalar distribution. This formulation provides robustness akin to generalized winding numbers (GWN), but provides distance function rather than just an inside/outside classification. Our formulation also offers several features not common to classic distance algorithms, such as the ability to simultaneously fit multiple level sets, a notion of distance for geometry that does not topologically bound any region, and the ability to mix and match signed and unsigned distance. The method can be applied in any dimension and to any spatial discretization, including triangle meshes, tet meshes, point clouds, polygonal meshes, voxelized surfaces, and regular grids. We evaluate the method on several challenging examples, implementing normal offsets and other morphological operations directly on imperfect curve and surface data. In many cases we also obtain an inside/outside classification dramatically more robust than the one obtained provided by GWN.

This paper describes an optimization-based framework for finding arrangements of 2D shapes subject to pairwise constraints. Such arrangements naturally arise in tasks such as vector illustration and diagram generation, but enforcing these criteria robustly is surprisingly challenging. We approach this problem through the minimization of novel energetic penalties, derived from the signed distance function of the Minkowski difference between interacting shapes. This formulation provides useful gradients even when initialized from a wildly infeasible state, and, unlike many common collision penalties, can handle open curves that do not have a well-defined inside and outside. Moreover, it supports rich features beyond the basic no-overlap condition, such as tangency, containment, and precise padding, which are especially valuable in the vector illustration context. We develop closed-form expressions and efficient approximations of our penalty for standard vector graphics primitives, yielding efficient evaluation and easy implementation within existing automatic differentiation pipelines. The method has already been “battle-tested” as a component of public-facing open source software; we demonstrate the utility of the framework via examples from illustration, data visualization, diagram generation, and geometry processing.

Mesh processing algorithms are often communicated via concise mathematical notation (e.g., summation over mesh neighborhoods). However, conversion of notation into working code remains a time consuming and error-prone process which requires arcane knowledge of low-level data structures and libraries—impeding rapid exploration of high-level algorithms. We address this problem by introducing a domain-specific language (DSL) for mesh processing called I❤️️️MESH, which resembles notation commonly used in visual and geometric computing, and automates the process of converting notation into code. The centerpiece of our language is a flexible notation for specifying and manipulating neighborhoods of a cell complex, internally represented via standard operations on sparse boundary matrices. This layered design enables natural expression of algorithms while minimizing demands on a code generation back-end. In particular, by integrating I❤️️️MESH with the linear algebra features of the I❤️️️LA DSL, and adding support for automatic differentiation, we can rapidly implement a rich variety of algorithms on point clouds, surface meshes, and volume meshes.

Grid-free Monte Carlo methods based on the walk on spheres (WoS) algorithm solve fundamental partial differential equations (PDEs) like the Poisson equation without discretizing the problem domain or approximating functions in a finite basis. Such methods hence avoid aliasing in the solution, and evade the many challenges of mesh generation. Yet for problems with complex geometry, practical grid-free methods have been largely limited to basic Dirichlet boundary conditions. This paper introduces the walk on stars (WoSt) method, which solves linear elliptic PDEs with arbitrary mixed Neumann and Dirichlet boundary conditions. The key insight is that one can efficiently simulate reflecting Brownian motion (which models Neumann conditions) by replacing the balls used by WoS with star-shaped domains; we identify such domains by locating the closest visible point on the geometric silhouette. Overall, WoSt retains many attractive features of other grid-free Monte Carlo methods, such as progressive and view-dependent evaluation, trivial parallelization, and sublinear scaling to increasing geometric detail.

In the plane, the winding number is the number of times a curve wraps around a given point. Winding numbers are a basic component of geometric algorithms such as point-in-polygon tests, and their generalization to data with noise or topological errors has proven valuable for geometry processing tasks ranging from surface reconstruction to mesh booleans. However, standard definitions do not immediately apply on surfaces, where not all curves bound regions. We develop a meaningful generalization, starting with the well-known relationship between winding numbers and harmonic functions. By processing the derivatives of such functions, we can robustly filter out components of the input that do not bound any region. Ultimately, our algorithm yields (i) a closed, completed version of the input curves, (ii) integer labels for regions that are meaningfully bounded by these curves, and (iii) the complementary curves that do not bound any region. The main computational cost is solving a standard Poisson equation, or for surfaces with nontrivial topology, a sparse linear program. We also introduce special basis functions to represent singularities that naturally occur at endpoints of open curves.

Grid-free Monte Carlo methods such as walk on spheres can be used to solve elliptic partial differential equations without mesh generation or global solves. However, such methods independently estimate the solution at every point, and hence do not take advantage of the high spatial regularity of solutions to elliptic problems. We propose a fast caching strategy which first estimates solution values and derivatives at randomly sampled points along the boundary of the domain (or a local region of interest). These cached values then provide cheap, output-sensitive evaluation of the solution (or its gradient) at interior points, via a boundary integral formulation. Unlike classic boundary integral methods, our caching scheme introduces zero statistical bias and does not require a dense global solve. Moreover we can handle imperfect geometry (e.g., with self-intersections) and detailed boundary/source terms without repairing or resampling the boundary representation. Overall, our scheme is similar in spirit to virtual point light methods from photorealistic rendering: it suppresses the typical salt-and-pepper noise characteristic of independent Monte Carlo estimates, while still retaining the many advantages of Monte Carlo solvers: progressive evaluation, trivial parallelization, geometric robustness, etc.  We validate our approach using test problems from visual and geometric computing.

This paper describes a method for fast simplification of surface meshes. Whereas past methods focus on visual appearance, our goal is to solve equations on the surface. Hence, rather than approximate the extrinsic geometry, we construct a coarse intrinsic triangulation of the input domain. In the spirit of the quadric error metric (QEM), we perform greedy decimation while agglomerating global information about approximation error. In lieu of extrinsic quadrics, however, we store intrinsic tangent vectors that track how far curvature “drifts” during simplification. This process also yields a bijective map between the fine and coarse mesh, and prolongation operators for both scalar- and vector-valued data. Moreover, we obtain hard guarantees on element quality via intrinsic retriangulation—a feature unique to the intrinsic setting. The overall payoff is a “black box” approach to geometry processing, which decouples mesh resolution from the size of matrices used to solve equations. We show how our method benefits several fundamental tasks, including geometric multigrid, all-pairs geodesic distance, mean curvature flow, geodesic Voronoi diagrams, and the discrete exponential map.

Partial differential equations (PDEs) with spatially varying coefficients arise throughout science and engineering, modeling rich heterogeneous material behavior. Yet conventional PDE solvers struggle with the immense complexity found in nature, since they must first discretize the problem—leading to spatial aliasing, and global meshing/sampling that is costly and error-prone. We describe a method that approximates neither the domain geometry, the problem data, nor the solution space, providing the exact solution (in expectation) even for problems with extremely detailed geometry and intricate coefficients. Our main contribution is to extend the walk on spheres (WoS) algorithm from constant- to variable-coefficient problems, by drawing on techniques from volumetric rendering. In particular, an approach inspired by null-scattering yields unbiased Monte Carlo estimators for a large class of 2nd order elliptic PDEs, which share many attractive features with Monte Carlo rendering: no meshing, trivial parallelism, and the ability to evaluate the solution at any point without solving a global system of equations.

We introduce a new approach to deep learning on 3D surfaces, based on the insight that a simple diffusion layer is highly effective for spatial communication. The resulting networks automatically generalize across different samplings and resolutions of a surface—a basic property which is crucial for practical applications. Our networks can be discretized on various geometric representations such as triangle meshes or point clouds, and can even be trained on one representation then applied to another. We optimize the spatial support of diffusion as a continuous network parameter ranging from purely local to totally global, removing the burden of manually choosing neighborhood sizes. The only other ingredients in the method are a multi-layer perceptron applied independently at each point, and spatial gradient features to support directional filters. The resulting networks are simple, robust, and efficient. Here, we focus primarily on triangle mesh surfaces, and demonstrate state-of-the-art results for a variety of tasks including surface classification, segmentation, and non-rigid correspondence.

This paper describes a numerically robust data structure for encoding triangulations of polyhedral surfaces. Existing data structures either rely on floating point values to encode connectivity, or do not support remeshing operations beyond basic edge flips. We instead provide an integer-based data structure that guarantees valid connectivity, even for meshes with near-degenerate elements. Our starting point is the framework of normal coordinates from geometric topology, which we extend to the broader set of operations needed for mesh processing (vertex insertion, edge splits, etc.). The resulting data structure can be used as a drop-in replacement for earlier schemes, automatically improving reliability across a wide variety of applications. As a stress test, we successfully compute an intrinsic Delaunay refinement and associated subdivision for all manifold meshes in the Thingi10k dataset. In turn, we can compute reliable and highly accurate solutions to partial differential equations even on extremely low-quality meshes.

Functionals that penalize bending or stretching of a surface play a key role in geometric and scientific computing, but to date have ignored a very basic requirement: in many situations, surfaces must not pass through themselves or each other. This paper develops a numerical framework for optimization of surface geometry while avoiding (self-)collision. The starting point is the tangent-point energy, which effectively pushes apart pairs of points that are close in space but distant along the surface. We develop a discretization of this energy for triangle meshes, and introduce a novel acceleration scheme based on a fractional Sobolev inner product. In contrast to similar schemes developed for curves, we avoid the complexity of building a multiresolution mesh hierarchy by decomposing our preconditioner into two ordinary Poisson equations, plus forward application of a fractional differential operator. We further accelerate this scheme via hierarchical approximation, and describe how to incorporate a variety of constraints (on area, volume, etc.). Finally, we explore how this machinery might be applied to problems in mathematical visualization, geometric modeling, and geometry processing.

This paper describes a numerical method for surface parameterization, yielding maps that are locally injective and discretely conformal in an exact sense. Unlike previous methods for discrete conformal parameterization, the method is guaranteed to work for any manifold triangle mesh, with no restrictions on triangulation quality or cone singularities. In particular we consider maps from surfaces of any genus (with or without boundary) to the plane, or globally bijective maps from genus zero surfaces to the sphere. Recent theoretical developments show that each task can be formulated as a convex problem where the triangulation is allowed to change—we complete the picture by introducing the machinery needed to actually construct a discrete conformal map. In particular, we introduce a new scheme for tracking correspondence between triangulations based on normal coordinates, and a new interpolation procedure based on layout in the light cone. Stress tests involving difficult cone configurations and near-degenerate triangulations indicate that the method is extremely robust in practice, and provides high-quality interpolation even on meshes with poor elements.

Curves play a fundamental role across computer graphics, physical simulation, and mathematical visualization, yet most tools for curve design do nothing to prevent crossings or self-intersections. This paper develops efficient algorithms for (self-)repulsion of plane and space curves that are well-suited to problems in computational design. Our starting point is the so-called tangent-point energy, which provides an infinite barrier to self-intersection. In contrast to local collision detection strategies used in, e.g., physical simulation, this energy considers interactions between all pairs of points, and is hence useful for global shape optimization: local minima tend to be aesthetically pleasing, physically valid, and nicely distributed in space. A reformulation of gradient descent, based on a Sobolev-Slobodeckij inner product enables us to make rapid progress toward local minima—independent of curve resolution. We also develop a hierarchical multigrid scheme that significantly reduces the per-step cost of optimization. The energy is easily integrated with a variety of constraints and penalties (e.g., inextensibility, or obstacle avoidance), which we use for applications including curve packing, knot untangling, graph embedding, non-crossing spline interpolation, flow visualization, and robotic path planning.

This 3-hour course provides a first introduction to intrinsic triangulations and their use in mesh processing algorithms. As geometric data becomes more ubiquitous, e.g., in applications such as augmented reality or machine learning, there is a pressing need to develop algorithms that work reliably on low-quality data. Intrinsic triangulations provide a powerful framework for these problems, by de-coupling the mesh used to encode geometry from the one used for computation. The basic shift in perspective is to encode the geometry of a mesh not in terms of ordinary vertex positions, but instead only in terms of edge lengths. Intrinsic triangulations have a long history in mathematics, but only in recent years have been applied to practical geometric computing. The course begins by giving motivation for intrinsic triangulations in terms of recent problems in computer graphics, followed by an interactive coding session where participants can make first contact with the idea of intrinsic meshes. We then give some mathematical background, and describe key data structures (overlay, signpost, normal coordinates). Using this machinery, we translate algorithms from computational geometry and scientific computing into cutting-edge algorithms for curved surfaces. For instance, we look at mesh parameterization, vector field processing, finding geodesics, solving partial differential equations (PDEs), and more. We also discuss processing of nonmanifold meshes and point clouds; participants can explore these algorithms via interactive demos. We conclude with a discussion of open questions and opportunities for future work.

Discrete Differential Geometry (DDG) is an emerging discipline at the boundary between mathematics and computer science. It aims to translate concepts from classical differential geometry into a language that is purely finite and discrete, and can hence be used by algorithms to reason about geometric data. In contrast to standard numerical approximation, the central philosophy of DDG is to faithfully and exactly preserve key invariants of geometric objects at the discrete level. This process of translation from smooth to discrete helps to both illuminate the fundamental meaning behind geometric ideas and provide useful algorithmic guarantees. This volume illustrates the principles of DDG via several recent topics: discrete nets, discrete differential operators, discrete mappings, discrete conformal geometry, and discrete optimal transport.

This paper introduces a new approach to computing geodesics on polyhedral surfaces—the basic idea is to iteratively perform edge flips, in the same spirit as the classic Delaunay flip algorithm. This process also produces a triangulation conforming to the output geodesics, which is immediately useful for tasks in geometry processing and numerical simulation. More precisely, our FlipOut algorithm transforms a given sequence of edges into a locally shortest geodesic while avoiding self-crossings (formally: it finds a geodesic in the same isotopy class). The algorithm is guaranteed to terminate in a finite number of operations; practical runtimes are on the order of a few milliseconds, even for meshes with millions of triangles. The same approach is easily applied to curves beyond simple paths, including closed loops, curve networks, and multiply-covered curves. We explore how the method facilitates tasks such as straightening cuts and segmentation boundaries, computing geodesic Bézier curves, extending the notion of constrained Delaunay triangulations (CDT) to curved surfaces, and providing accurate boundary conditions for partial differential equations (PDEs). Evaluation on challenging datasets such as Thingi10k indicates that the method is both robust and efficient, even for low-quality triangulations.

We describe a discrete Laplacian suitable for any triangle mesh, including those that are nonmanifold or nonorientable (with or without boundary). Our Laplacian is a robust drop-in replacement for the usual cotan matrix, and is guaranteed to have nonnegative edge weights on both interior and boundary edges, even for extremely poor-quality meshes. The key idea is to build what we call a “tufted cover” over the input domain, which has nonmanifold vertices but manifold edges. Since all edges are manifold, we can flip to an intrinsic Delaunay triangulation; our Laplacian is then the cotan Laplacian of this new triangulation. This construction also provides a high-quality point cloud Laplacian, via a nonmanifold triangulation of the point set. We validate our Laplacian on a variety of challenging examples (including all models from Thingi10k), and a variety of standard tasks including geodesic distance computation, surface deformation, parameterization, and computing minimal surfaces.

This paper explores how core problems in PDE-based geometry processing can be efficiently and reliably solved via grid-free Monte Carlo methods. Modern geometric algorithms often need to solve Poisson-like equations on geometrically intricate domains. Conventional methods most often mesh the domain, which is both challenging and expensive for geometry with fine details or imperfections (holes, self-intersections, etc.). In contrast, grid-free Monte Carlo methods avoid mesh generation entirely, and instead just evaluate closest point queries. They hence do not discretize space, time, nor even function spaces, and provide the exact solution (in expectation) even on extremely challenging models. More broadly, they share many benefits with Monte Carlo methods from photorealistic rendering: excellent scaling, trivial parallel implementation, view-dependent evaluation, and the ability to work with any kind of geometry (including implicit or procedural descriptions). We develop a complete “black box” solver that encompasses integration, variance reduction, and visualization, and explore how it can be used for various geometry processing tasks. In particular, we consider several fundamental linear elliptic PDEs with constant coefficients on solid regions of Rⁿ. Overall we find that Monte Carlo methods significantly broaden the horizons of geometry processing, since they easily handle problems of size and complexity that are essentially hopeless for conventional methods.

We introduce a system called Penrose for creating mathematical diagrams. Its basic functionality is to translate abstract statements written in familiar math-like notation into one or more possible visual representations. Rather than rely on a fixed library of visualization tools, the visual representation is user-defined in a constraint-based specification language; diagrams are then generated automatically via constrained numerical optimization. The system is user-extensible to many domains of mathematics, and is fast enough for iterative design exploration. In contrast to tools that specify diagrams via direct manipulation or low-level graphics programming, Penrose enables rapid creation and exploration of diagrams that faithfully preserve the underlying mathematical meaning. We demonstrate the effectiveness and generality of the system by showing how it can be used to illustrate a diverse set of concepts from mathematics and computer graphics.

We present a data structure that makes it easy to run a large class of algorithms from computational geometry and scientific computing on extremely poor-quality surface meshes. Rather than changing the geometry, as in traditional remeshing, we consider intrinsic triangulations which connect vertices by straight paths along the exact geometry of the input mesh. Our key insight is that such a triangulation can be encoded implicitly by storing the direction and distance to neighboring vertices. The resulting signpost data structure then allows geometric and topological queries to be made on-demand by tracing paths across the surface. Existing algorithms can be easily translated into the intrinsic setting, since this data structure supports the same basic operations as an ordinary triangle mesh (vertex insertions, edge splits, etc.). The output of intrinsic algorithms can then be stored on an ordinary mesh for subsequent use; unlike previous data structures, we use a constant amount of memory and do not need to explicitly construct an overlay mesh unless it is specifically requested. Working in the intrinsic setting incurs little computational overhead, yet we can run algorithms on extremely degenerate inputs, including all manifold meshes from the Thingi10k data set. To evaluate our data structure we implement several fundamental geometric algorithms including intrinsic versions of Delaunay refinement and optimal Delaunay triangulation, approximation of Steiner trees, adaptive mesh refinement for PDEs, and computation of Poisson equations, geodesic distance, and flip-free tangent vector fields.

A basic challenge in field-guided hexahedral meshing is to find a spatially-varying frame that is adapted to the domain geometry and is continuous up to symmetries of the cube. We introduce a fundamentally new representation of such 3D cross fields based on Cartan's method of moving frames. Our key observation is that cross fields and ordinary frame fields are locally characterized by identical conditions on their Darboux derivative. Hence, by using derivatives as the principal representation (and only later recovering the field itself), one avoids the need to explicitly account for symmetry during optimization. At the discrete level, derivatives are encoded by skew-symmetric matrices associated with the edges of a tetrahedral mesh; these matrices encode arbitrarily large rotations along each edge, and can robustly capture singular behavior even on coarse meshes. We apply this representation to compute 3D cross fields that are as smooth as possible everywhere but on a prescribed network of singular curves—since these fields are adapted to curve tangents, they can be directly used as input for field-guided mesh generation algorithms. Optimization amounts to an easy nonlinear least squares problem that behaves like a convex program in the sense that it always appears to produce the same result, independent of initialization. We study the numerical behavior of this procedure, and perform some preliminary experiments with mesh generation.

This paper describes a method for efficiently computing parallel transport of tangent vectors on curved surfaces, or more generally, any vector-valued data on a curved manifold. More precisely, it extends a vector field defined over any region to the rest of the domain via parallel transport along shortest geodesics. This basic operation enables fast, robust algorithms for extrapolating level set velocities, inverting the exponential map, computing geometric medians and Karcher/Fréchet means of arbitrary distributions, constructing centroidal Voronoi diagrams, and finding consistently ordered landmarks. Rather than evaluate parallel transport by explicitly tracing geodesics, we show that it can be computed via a short-time heat flow involving the connection Laplacian. As a result, transport can be achieved by solving three prefactored linear systems, each akin to a standard Poisson problem. To implement the method we need only a discrete connection Laplacian, which we describe for a variety of geometric data structures (point clouds, polygon meshes, \etc). We also study the numerical behavior of our method, showing empirically that it converges under refinement, and augment the construction of intrinsic Delaunay triangulations (iDT) so that they can be used in the context of tangent vector field processing.

Conformal geometry studies geometric properties that are invariant with respect to angle-preserving transformations. In the discrete setting of polyhedral surfaces, the idea of naively preserving angles leads to a notion of conformal geometry that is far too rigid, i.e., it does not capture the rich behavior found in the smooth setting. These notes explore several alternative notions of discrete conformal structure, culminating with a recent uniformization theorem for general simplicial surfaces. Topics covered include circle packings, circle patterns, inversive distance, discrete Yamabe flow, and connections to variational principles for ideal hyperbolic polyhedra. (This chapter is an extended version of course notes developed for the 2018 AMS Short Course on Discrete Differential Geometry.)

Developable surfaces are those that can be made by smoothly bending flat pieces without stretching or shearing. We introduce a definition of developability for triangle meshes which exactly captures two key properties of smooth developable surfaces, namely flattenability and presence of straight ruling lines. This definition provides a starting point for algorithms in developable surface modeling—we consider a variational approach that drives a given mesh toward developable pieces separated by regular seam curves. Computation amounts to gradient descent on an energy with support in the vertex star, without the need to explicitly cluster patches or identify seams. We briefly explore applications to developable design and manufacturing.

This paper develops a global variational approach to cutting curved surfaces so that they can be flattened into the plane with low metric distortion. Such cuts are a critical component in a variety of algorithms that seek to parameterize surfaces over flat domains, or fabricate structures from flat materials. Rather than evaluate the quality of a cut solely based on properties of the curve itself (e.g., its length or curvature), we formulate a flow that directly optimizes the distortion induced by cutting and flattening. Notably, we do not have to explicitly parameterize the surface in order to evaluate the cost of a cut, but can instead integrate a simple evolution equation defined on the cut curve itself. We arrive at this flow via a novel application of shape derivatives to the Yamabe equation from conformal geometry. We then develop an Eulerian numerical integrator on triangulated surfaces, which does not restrict cuts to mesh edges and can incorporate user-defined data such as importance or occlusion. The resulting cut curves can be used to drive distortion to arbitrarily low levels, and have a very different character from cuts obtained via purely discrete formulations. We briefly explore potential applications to computational design, as well as connections to space filling curves and the problem of uniform heat distribution.

Angle-preserving or conformal surface parameterization has proven to be a powerful tool across applications ranging from geometry processing, to digital manufacturing, to machine learning, yet conformal maps can still suffer from severe area distortion. Cone singularities provide a way to mitigate this distortion, but finding the best configuration of cones is notoriously difficult. This paper develops a strategy that is globally optimal in the sense that it minimizes total area distortion among all possible cone configurations (number, placement, and size) that have no more than a fixed total cone angle. A key insight is that, for the purpose of optimization, one should not work directly with curvature measures (which naturally represent cone configurations), but can instead apply Fenchel-Rockafellar duality to obtain a formulation involving only ordinary functions. The result is a convex optimization problem, which can be solved via a sequence of sparse linear systems easily built from the usual cotangent Laplacian. The method supports user-defined notions of importance, constraints on cone angles (e.g., positive, or within a given range), and sophisticated boundary conditions (e.g., convex, or polygonal). We compare our approach to previous techniques on a variety of challenging models, often achieving dramatically lower distortion, and demonstrating that global optimality leads to extreme robustness in the presence of noise or poor discretization.

Deployable structures are physical mechanisms that can easily transition between two or more geometric configurations; such structures enable industrial, scientific, and consumer applications at a wide variety of scales. This paper develops novel deployable structures that can approximate a large class of doubly-curved surfaces and are easily actuated from a flat initial state via inflation or gravitational loading. The structures are based on two-dimensional rigid mechanical linkages that implicitly encode the curvature of the target shape via a user-programmable pattern that permits locally isotropic scaling under load. We explicitly characterize the shapes that can be realized by such structures—in particular, we show that they can approximate target surfaces of positive mean curvature and bounded scale distortion relative to a given reference domain. Based on this observation, we develop efficient computational design algorithms for approximating a given input geometry. The resulting designs can be rapidly manufactured via digital fabrication technologies such as laser cutting, CNC milling, or 3D printing. We validate our approach through a series of physical prototypes and present several application case studies, ranging from surgical implants to large-scale deployable architecture.

Conformal parameterizations over the sphere provide high-quality maps between genus zero surfaces, and are essential for applications such as data transfer and comparative shape analysis. However, such maps are not unique: to define correspondence between two surfaces, one must find the Möbius transformation that best aligns two parameterizations—akin to picking a translation and rotation in rigid registration problems. We describe a simple procedure that canonically centers and rotationally aligns two spherical maps. Centering is implemented via elementary operations on triangle meshes in R³, and minimizes area distortion. Alignment is achieved using the FFT over the group of rotations. We examine this procedure in the context of spherical conformal parameterization, orbifold maps, non-rigid symmetry detection, and dense point-to-point surface correspondence.

A conformal flattening maps a curved surface to the plane without distorting angles—such maps have become a fundamental building block for problems in geometry processing, numerical simulation, and computational design. Yet existing methods provide little direct control over the shape of the flattened domain, or else demand expensive nonlinear optimization. Boundary first flattening (BFF) is a linear method for conformal parameterization which is faster than traditional linear methods, yet provides control and quality comparable to sophisticated nonlinear schemes. The key insight is that the boundary data for many conformal mapping problems can be efficiently constructed via the Cherrier formula together with a pair of Poincaré-Steklov operators; once the boundary is known, the map can be easily extended over the rest of the domain. Since computation demands only a single factorization of the real Laplace matrix, the amortized cost is about 50x less than any previously published technique for boundary-controlled conformal flattening. As a result, BFF opens the door to real-time editing or fast optimization of high-resolution maps, with direct control over boundary length or angle. We show how this method can be used to construct maps with sharp corners, cone singularities, minimal area distortion, and uniformization over the unit disk; we also demonstrate for the first time how a surface can be conformally flattened directly onto any given target shape.

We introduce the heat method for solving the single- or multiple-source shortest path problem on both flat and curved domains. A key insight is that distance computation can be split into two stages: first find the direction along which distance is increasing, then compute the distance itself. The heat method is robust, efficient, and simple to implement since it is based on solving a pair of standard sparse linear systems. These systems can be factored once and subsequently solved in near-linear time, substantially reducing amortized cost. Real-world performance is an order of magnitude faster than state-of-the-art methods, while maintaining a comparable level of accuracy. The method can be applied in any dimension, and on any domain that admits a gradient and inner product—including regular grids, triangle meshes, and point clouds. Numerical evidence indicates that the method converges to the exact distance in the limit of refinement; we also explore smoothed approximations of distance suitable for applications where greater regularity is desired.

The emerging field of discrete differential geometry (DDG) studies discrete analogues of smooth geometric objects, providing an essential link between analytical descriptions and computation. In recent years it has unearthed a rich variety of new perspectives on applied problems in computational anatomy/biology, computational mechanics, industrial design, computational architecture, and digital geometry processing at large. The basic philosophy of discrete differential geometry is that a discrete object like a polyhedron is not merely an approximation of a smooth one, but rather a differential geometric object in its own right. In contrast to traditional numerical analysis which focuses on eliminating approximation error in the limit of refinement (e.g., by taking smaller and smaller finite differences), DDG places an emphasis on the so-called “mimetic” viewpoint, where key properties of a system are preserved exactly, independent of how large or small the elements of a mesh might be. Just as algorithms for simulating mechanical systems might seek to exactly preserve physical invariants such as total energy or momentum, structure-preserving models of discrete geometry seek to exactly preserve global geometric invariants such as total curvature. More broadly, DDG focuses on the discretization of objects that do not naturally fall under the umbrella of traditional numerical analysis. This article provides an overview of some of the themes in DDG.

Creating mathematical diagrams is essential for both developing one’s intuition and conveying it to others. However, formalizing diagrams in most general-purpose tools requires painstaking low-level manipulation of shapes and positions. We report on early work on PENROSE, a system we are building to automatically visualize mathematics from notation. PENROSE comprises two languages: SUBSTANCE, a domain-specific language that mimics the declarativeness of mathematical notation, and STYLE, a styling language that concisely specifies the visual semantics of the notation. Our system can automatically visualize set theory expressions with user-defined styles, and it can visualize abstract definitions of functions by producing concrete examples. We plan to extend the system to more domains of math.

The eigenfunctions and eigenvalues of the Laplace-Beltrami operator have proven to be a powerful tool for digital geometry processing, providing a description of geometry that is essentially independent of coordinates or the choice of discretization. However, since Laplace-Beltrami is purely intrinsic it struggles to capture important phenomena such as extrinsic bending, sharp edges, and fine surface texture. We introduce a new extrinsic differential operator called the relative Dirac operator, leading to a family of operators with a continuous trade-off between intrinsic and extrinsic features. Previous operators are either fully or partially intrinsic. In contrast, the proposed family spans the entire spectrum: from completely intrinsic (depending only on the metric) to completely extrinsic (depending only on the Gauss map). By adding an infinite potential well to this (or any) operator we can also robustly handle surface patches with irregular boundary. We explore use of these operators for a variety of shape analysis tasks, and study their performance relative to operators previously found in the geometry processing literature.

In science, a well-chosen illustration can turn bafflement into enlightenment. Yet technical exposition remains largely textual, due to the tremendous expertise required to produce high-quality figures. We propose PENROSE, a system to automatically generate professional-quality mathematical illustrations from high-level, purely semantic descriptions of mathematical objects. Unlike low-level tools where diagrams are specified via graphical primitives, a mathematically-inclined user should not require any graphic design skill to create beautiful diagrams. PENROSE comprises two extensible domain-specific languages (DSLs): SUBSTANCE, which diagrammers use to specify mathematical objects and relationships, and STYLE, which implementers use to encode various ways of realizing these relationships visually, akin to the separation of content and style in modern HTML/CSS. To compile diagrams, we are developing a sophisticated constraint solver incorporating techniques from optimization and computer graphics.

Telescoping structures are valuable for a variety of applications where mechanisms must be compact in size and yet easily deployed. So far, however, there has been no systematic study of the types of shapes that can be modeled by telescoping structures, nor practical tools for telescopic design. We present a novel geometric characterization of telescoping curves, and explore how free-form surfaces can be approximated by networks of such curves. In particular we consider piecewise helical space curves with torsional impulses, which significantly generalize the linear telescopes found in typical engineering designs. Based on this principle we develop a system for computational design and fabrication which allows users to explore the space of telescoping structures; inputs to our system include user sketches or arbitrary meshes, which are then converted to a curve skeleton. We prototype applications in animation, fabrication, and robotics, using our system to design a variety of both simulated and fabricated examples.

This paper presents a method for constructing smooth unit n-direction fields (line fields, cross fields, etc.) on surfaces that is an order of magnitude faster than state-of-the-art methods, while still producing fields of equal or better quality. The method is based on a simple quadratic energy whose minimizers are globally optimal in the sense that they produce the smoothest fields over all possible configurations of singularities (number, location, and index). The method is fully automatic and can optionally produce fields aligned with a given guidance field, for example, principal curvature directions. Computationally the smoothest field is found via a sparse eigenvalue problem involving a matrix similar to the cotan-Laplacian. When a guidance field is present, finding the optimal field amounts to solving a single linear Poisson problem.

This course provides an introduction to geometry processing using discrete exterior calculus (DEC). DEC provides a simple, flexible, and efficient framework within which one can build a unified platform for geometry processing. The course provides essential mathematical background as well as a large array of real-world examples. It also provides a short survey of the most relevant recent developments in digital geometry processing and discrete differential geometry.

This paper introduces a new method for computing conformal transformations of triangle meshes in R³. Conformal maps are desirable in digital geometry processing because they do not exhibit shear, and therefore preserve texture fidelity as well as the quality of the mesh itself. Traditional discretizations consider maps into the complex plane, which are useful only for problems such as surface parameterization and planar shape deformation where the target surface is flat. We instead consider maps into the quaternions, which allows us to work directly with surfaces sitting in R³. In particular, we introduce a quaternionic Dirac operator and use it to develop a novel integrability condition on conformal deformations. Our discretization of this condition results in a sparse linear system that is simple to build and can be used to efficiently edit surfaces by manipulating curvature and boundary data, as demonstrated via several mesh processing applications.

Connections provide a way to compare local quantities defined at different points of a geometric space. This thesis develops a discrete theory of connections that naturally leads to practical, efficient numerical algorithms for geometry processing, including texture synthesis and quadrilateral remeshing. Our formulation is motivated by real-world applications where meshes may be noisy or coarsely discretized. Further, because our discrete framework closely parallels the smooth theory, we can draw upon a huge wealth of existing knowledge to develop and interpret mesh processing algorithms. The solutions we produce are globally optimal in the sense that they describe the trivial connection closest to Levi-Civita among all solutions with the prescribed set of singularities. Relative to previous methods our algorithm is surprisingly simple, and can be implemented using standard operations from mesh processing and linear algebra.

This paper presents a straightforward algorithm for constructing connections on surfaces that are as smooth as possible everywhere but on a prescribed set of isolated singularities with given index. Such connections can be used to design rotationally symmetric direction fields, which are essential in applications such as quadrilateral remeshing and texture synthesis. We compute connections by solving a single sparse linear system built from standard operators. Our solutions are globally optimal in the sense that they describe the trivial connection closest to Levi-Civita among all connections with the prescribed set of singularities. Relative to previous methods our algorithm is surprisingly simple, and can be implemented using standard operations from mesh processing and linear algebra.

This paper is concerned with the animation and control of vehicles with complex dynamics such as helicopters, boats, and cars. Motivated by recent developments in discrete geometric mechanics we develop a general framework for integrating the dynamics of holonomic and nonholonomic vehicles by preserving their state-space geometry and motion invariants. We demonstrate that the resulting integration schemes are superior to standard methods in numerical robustness and efficiency, and can be applied to many types of vehicles. In addition, we show how to use this framework in an optimal control setting to automatically compute accurate and realistic motions for arbitrary user-specified constraints.

Numerical viscosity has long been a problem in fluid animation. Existing methods suffer from intrinsic artificial dissipation and often apply complicated computational mechanisms to combat such effects. Consequently, dissipative behavior cannot be controlled or modeled explicitly in a manner independent of time step size, complicating the use of coarse previews and adaptive time-stepping methods. This paper proposes simple, unconditionally stable, fully Eulerian integration schemes with no numerical viscosity that are capable of maintaining the liveliness of fluid motion without recourse to corrective devices. Pressure and fluxes are solved efficiently and simultaneously in a time-reversible manner on simplicial grids, and the energy is preserved exactly over long time scales in the case of inviscid fluids. These integrators can be viewed as an extension of the classical energy-preserving Harlow-Welch / Crank-Nicolson scheme to simplicial grids.

Reference grids are commonly used in design software to help users to judge distances and to understand the orientation of the virtual workspace. However, despite their ubiquity in 3D graphics authoring applications, little research has gone into important design considerations of the 3D reference grids themselves, which directly impact their usefulness. In an effort to resolve some of these outstanding issues, we have developed two new techniques; the multiscale reference grid and position pegs that form a consistent foundation for presenting relative scale and position information to the user.
Our design of a multiscale reference grid consistently subdivides and coalesces gridlines, based on the computation of a closeness metric, while ensuring that there are neither too many nor too few subdivisions. Position pegs extend the grid so that objects that are lying above or below the ground plane can be brought into a common environmental frame of reference without interfering with the grid or object data. We introduce an analytical system, similar to MIP mapping, to provide a stable viewpoint-determined result. Our design solves several depth cue problems in a way that is independent of viewing projection.

We capture the shape of moving cloth using a custom set of color markers printed on the surface of the cloth. The output is a sequence of triangle meshes with static connectivity and with detail at the scale of individual markers in both smooth and folded regions. We compute markers' coordinates in space using marker correspondence across multiple synchronized video cameras. Correspondence is determined from color information in small neighborhoods and refined using a novel strain pruning process. Final correspondence does not require neighborhood information. We use a novel data driven hole-filling technique to fill occluded regions. Our results include several challenging examples: a wrinkled shirt sleeve, a dancing pair of pants, and a rag tossed onto a cup. Finally, we demonstrate that cloth capture is reusable by animating a pair of pants using human motion capture data.

We present a new method for cloth animation based on data driven synthesis. In contrast to approaches that focus on physical simulation, we animate cloth by manipulating short sequences of existing cloth animation. While our source of data is cloth animation captured using video cameras, the method is equally applicable to simulation data. The approach has benefits in both cases: current cloth capture is limited because small tweaks to the data require filming an entirely new sequence. Likewise, simulation suffers from long computation times and complications such as tangling. In this sketch we create new animations by fitting cloth animation to human motion capture data, i.e., we drive the cloth with a skeleton.

Physically based animation of fluids such as smoke, water, and fire provides some of the most stunning visuals in computer graphics, but has historically been the domain of high-quality offline rendering due to great computational cost. In this chapter we not only show how these effects can be simulated and rendered in real time, but also how they can be seamlessly integrated into real time applications.

Many mesh parameterization algorithms have focused on minimizing distortion and utilizing texture area, but few have addressed issues related to processing a signal on the mesh surface. We present an algorithm which partitions a mesh into rectangular charts while preserving a one-to-one texel correspondence across chart boundaries. This mapping permits any computation on the mesh surface which is typically carried out on a regular grid, and prevents seams by ensuring resolution continuity along the boundary. These features are also useful for traditional texture applications such as surface painting where continuity is important. Distortion is comparable to other parameterization schemes, and the rectangular charts yield efficient packing into a texture atlas. We apply this parameterization to texture synthesis, fluid simulation, mesh processing and storage, and locating geodesics.

Using the GPU to accelerate ray tracing may seem like a natural choice due to the highly parallel nature of the problem. However, determining the most versatile GPU data structure for scene storage and traversal is a challenge. In this paper, we introduce a new method for quick intersection of triangular meshes on the GPU. The method uses a threaded bounding volume hierarchy built from a geometry image, which can be efficiently traversed and constructed entirely on the GPU. This acceleration scheme is highly competitive with other GPU ray tracing methods, while allowing for both dynamic geometry and an efficient level of detail scheme at no extra cost.

Pluto and its satellite Charon regularly occulted or transited each other's disks from 1985 through 1990. The light curves resulting from these events (collectively called "mutual events") have been used to determine albedo maps of Pluto's sub-Charon hemisphere. We now use a data set of four light curves that were obtained in both B and V Johnson filters to construct a two-color map of Pluto's surface. We are able to resolve the central part of Pluto's sub-Charon hemisphere. We find that the dark albedo feature that forms a band below Pluto's equator is comprised of several distinct color units. We detect ratios of V-filter/B-filter normal reflectances ranging from 1.15 to 1.39 on Pluto's sub-Charon hemisphere.