Algebraic Representations for Volumetric Frame Fields
AAlgebraic Representations for Volumetric Frame Fields
DAVID PALMER,
Massachusetts Institute of Technology, USA
DAVID BOMMES,
University of Bern, Switzerland
JUSTIN SOLOMON,
Massachusetts Institute of Technology, USA
Fig. 1. Octahedral fields generated with our methods on various models.
Field-guided parametrization methods have proven effective for quad mesh-ing of surfaces; these methods compute smooth cross fields to guide themeshing process and then integrate the fields to construct a discrete mesh.A key challenge in extending these methods to three dimensions, however,is representation of field values. Whereas cross fields can be represented bytangent vector fields that form a linear space, the 3D analog—an octahedralframe field—takes values in a nonlinear manifold. In this work, we describethe space of octahedral frames in the language of differential and algebraicgeometry. With this understanding, we develop geometry-aware tools foroptimization of octahedral fields, namely geodesic stepping and exact projec-tion via semidefinite relaxation. Our algebraic approach not only providesan elegant and mathematically-sound description of the space of octahedralframes but also suggests a generalization to frames whose three axes scaleindependently, better capturing the singular behavior we expect to see involumetric frame fields. These new odeco frames , so-called as they are rep-resented by orthogonally decomposable tensors, also admit a semidefiniteprogram–based projection operator. Our description of the spaces of octahe-dral and odeco frames suggests computing frame fields via manifold-basedoptimization algorithms; we show that these algorithms efficiently producehigh-quality fields while maintaining stability and smoothness.CCS Concepts: •
Computing methodologies → Computer graphics ; Mesh geometry models ; Volumetric models .Additional Key Words and Phrases: Hexahedral meshing, octahedral framefields, convex relaxations, convex algebraic geometry
Authors’ addresses: David Palmer, Massachusetts Institute of Technology, 77 Mas-sachusetts Avenue, Cambridge, MA, 02139, USA, [email protected]; David Bommes, Univer-sity of Bern, Hochschulstrasse 6, Bern, 3012, Switzerland, [email protected];Justin Solomon, Massachusetts Institute of Technology, 77 Massachusetts Avenue,Cambridge, MA, 02139, USA, [email protected] to make digital or hard copies of part or all of this work for personal orclassroom use is granted without fee provided that copies are not made or distributedfor profit or commercial advantage and that copies bear this notice and the full citationon the first page. Copyrights for third-party components of this work must be honored.For all other uses, contact the owner/author(s).© 2020 Copyright held by the owner/author(s).0730-0301/2020/3-ART16https://doi.org/10.1145/3366786
Fig. 2. We compute frame fields as maps into the octahedral and odecovarieties , a projected slice of which is depicted here.
ACM Reference Format:
David Palmer, David Bommes, and Justin Solomon. 2020. Algebraic Repre-sentations for Volumetric Frame Fields.
ACM Trans. Graph.
39, 2, Article 16(March 2020), 17 pages. https://doi.org/10.1145/3366786
Inspired by the success of field–based approaches to quadrilateralmeshing on surfaces ( cf. [Vaxman et al. 2016]), recent research inapplied geometry has focused on developing an analogous approachto hexahedral meshing. Motivated by applications in finite-elementmodeling, hexahedral meshing is the problem of dividing a givenvolume into hexahedral elements (deformed cubes) with minimaldistortion and such that mesh boundary faces are aligned to theboundary of the volume.Hexahedral meshing couples a geometry problem—minimizingdistortion of mesh elements—to a combinatorial problem—placingmesh elements to achieve a desired connectivity structure. As in 2D,field-based approaches first ignore combinatorial constraints andsolve for a frame field , which represents the local alignment andsingular structure of a mesh (see Figure 1). Then, they integrate thatfield to guide the placement of hex elements [Nieser et al. 2011]. Soas not to impose unnatural constraints, the space of frame fields
ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020. a r X i v : . [ c s . G R ] A p r must be expressive enough to represent the range of possible singu-larities that may appear in hexahedral meshes. These singularitiesare described by gluing relations restricted to the symmetries of acube, i.e., the octahedral group.One might hope that 2D field-based methods would extend easilyto 3D. There are at least two obstructions to transferring ideas fromthe 2D case. First, the singular structure of a 3D frame field can bemuch more complicated than that of a cross field, comprising anembedded graph rather than a set of isolated points. Second, wemust understand the geometry of the space of frames to measureand optimize smoothness of a frame field.Complementing the recent results of Liu et al. [2018] on framefield singularities and necessary conditions for hex-meshability, inthe present work we study the second challenge, namely charac-terizing the geometry of the space of frames. A 3D frame field is,intuitively, an assignment of three mutually orthogonal directionsto each point in a volume. An orthonormal basis of three vectors—comprising an orthogonal matrix—is sufficient to specify a frame,but thanks to octahedral symmetry, multiple orthogonal matricescan specify the same frame. Formally, the space of octahedral framescan be viewed as the quotient of the group of 3D rotations SO ( ) by the right action of the octahedral group O comprising the rota-tional symmetries of a cube (see §3). The non-commutativity of 3Drotations makes the geometry of this quotient more complicatedthan that of its 2D counterpart. Recent attempts to lift ideas andtechniques from 2D have either ignored the geometry of the framespace entirely or treated it as a black box for nonlinear optimization.Consider perhaps the simplest form of optimization on a space,projection—finding the closest point in the space to a given pointin an ambient space. Previous work on frame fields has treated thisprojection problem as a nonlinear, nonconvex optimization problemover frames parametrized by Euler angles, with no guarantees onconvergence or global optimality. Our description of the octahedralframe space as an algebraic variety suggests a different approachto projection based on semidefinite programming, which yields acertificate of global optimality in polynomial time. Our semidefiniterelaxation of projection is exact in a neighborhood of the octahedralvariety, and we conjecture—with strong empirical evidence—that itis so universally.Even when conducting local optimization on the space of oc-tahedral frames, parametrization by Euler angles may not be thebest approach. We show that the map from SO ( ) to the octahedralvariety is a local isometry, enabling us to compute geodesics onthe octahedral variety in closed form. Manifold-based optimizationthat moves along geodesics can then be used to accelerate localoptimization dramatically.Beyond precisely characterizing the space of octahedral frames,our algebraic approach admits a generalization to frames whose axesscale independently. This larger space better captures frame fieldgeometry, e.g. allowing for a nonzero direction aligned to singulararcs even if the directions orthogonal to the arcs must vanish. Wecall these new objects odeco frames , thanks to their constructionusing orthogonally decomposable tensors, and we derive relevantprojection operators.Our experiments show how the theoretical objects we study en-able volumetric frame field design in practice. In particular, we apply standard manifold-based optimization algorithms to field design,built on our differential and algebraic descriptions of octahedral andodeco frames (see Figure 2). The end result is an efficient suite oftechniques for producing smooth fields that obey typical constraintsfor our target application. Outline.
In §3.1 we reintroduce the spherical harmonic represen-tation of octahedral frames and demonstrate how this amounts toan equivariant isometric embedding of the quotient SO ( )/ O into R . In §3.2 we introduce the odeco frames, whose axes can scaleindependently, and we exhibit the spaces of octahedral and odecoframes as nested varieties cut out by quadratic equations. In §4we describe essential primitives for optimization over octahedraland odeco frames, namely geodesics and projection via semidefi-nite programming. In §5 we formulate the frame field optimizationproblem. In §6 we describe two algorithms for optimizing fields.§7 describes experiments using these algorithms. A discussion andconclusion follow in §8. Further experimental results are includedin the supplemental material. In summary, our contributions are • a proof of isometric embedding of SO ( )/ O in R ; • descriptions of the spaces of octahedral and more general odecoframes as nested algebraic varieties; and • new state-of-the-art optimization techniques for volumetric framefields valued in both varieties, featuring geodesics and SDP-basedprojection as primitives. Cross fields, and their application to quad meshing, have been stud-ied extensively in geometry processing ( cf. [Vaxman et al. 2016]). Auseful insight from cross field research is that it is advantageous toreplace field values defined up to some symmetry with a representa-tion vector —that is, some function of the field value invariant underthe relevant symmetry. In two dimensions, four vectors forming aright-angled cross can be represented in a unified manner by theircommon complex fourth power. This amounts to an embedding ofa quotient manifold into a Euclidean space; we show that the octa-hedral variety generalizes this idea to 3D in a natural and isometricmanner.Recently, an effort has been made to unify and formalize thevarious algorithmic approaches to cross fields, borrowing from theGinzburg–Landau theory from physics [Beaufort et al. 2017; Ostingand Wang 2017; Viertel and Osting 2019]. This amounts to replacingan ill-posed unit norm constraint with a penalty term and taking thelimit as the penalty parameter goes to infinity. Viertel and Osting[2019] show that local optima of such a procedure have isolatedsingularities with indices ± /
4, as appropriate for quad meshing.They propose a diffusion-generated algorithm to compute such localoptima. Inspired by the work of Merriman, Bence, and Osher (MBO)[1992] on mean curvature flow, this algorithm alternates betweenfinite-time diffusion and pointwise normalization. Osting and Wang[2017] study an analogous method applied to orthogonal matrix–valued fields. Our projection operators enable us to develop similarMBO diffusion-generated methods for optimization of octahedraland odeco fields. In §7 we demonstrate that a modified strategy,
ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020. lgebraic Representations for Volumetric Frame Fields • 16:3 where the diffusion parameter is adjusted on-the-fly, frequentlyhelps to avoid local minima.
Huang et al. [2011] introduce a representation of frames as functionsover the sphere that exhibit octahedral symmetry, parametrized bycoefficients in the spherical harmonic basis. As an initializationstep, they solve a Laplace equation, resulting in coefficients in theinterior that do not correspond to octahedral frames. These mustbe projected via nonconvex optimization over frames parametrizedby Euler angles. They further optimize the results by minimizing adiscrete Dirichlet energy over Euler angles. Ray et al. [2016] refinethe three stages of this approach, with improved boundary con-straints in the initialization, an efficient projection technique, andan L-BFGS optimization algorithm.Solomon et al. [2017] reformulate the initialization step of Rayet al. [2016] using the boundary element method (BEM). This pro-vides a way to harmonically interpolate the boundary conditionsexactly, ignoring the constraints in the interior. Finally, sampledinterior values are approximately projected onto the constraintsas in the previous methods. As the constraints are ignored at theDirichlet energy–minimization stage, there is no sense in which thefinal frame fields are optimal.A related but distinct problem is the computation of fields ofsymmetric second order tensors (i.e., symmetric matrices) [Palacioset al. 2017]. Every symmetric matrix has an orthonormal basis ofeigenvectors, which corresponds to an octahedral frame. One mightthus think symmetric matrix fields can be used to parametrize octa-hedral frame fields. While a symmetric matrix field corresponds toat least one frame field, the correspondence is not one-to-one—forexample, a field of identity matrices corresponds to all frame fields.Moreover, symmetric matrix fields can only represent singularitieswhose indices are multiples of ± /
2, while indices ± / de novo framefield design. However, the connection associated to a field is a natu-ral object for the study of such properties as integrability. We leaveto future work the study of connections associated to fields valuedin the octahedral and odeco varieties and the extraction of suchconnections directly from the field coefficients. Complementing the spherical harmonic–based representation inthe papers above on octahedral frame fields, Chemin et al. [2018]propose an equivalent representation of octahedral frames as certainsymmetric tensors of order four. They introduce algebraic equationscharacterizing the octahedral frames among symmetric fourth-ordertensors. These equations are equivalent under a change of basis toour defining equations for the octahedral variety. However, by usinga basis suggested by the structure of SO ( ) , we are able to not onlypresent the defining equations of the octahedral variety, but alsoilluminate how it is an isometric embedding of the quotient spaceSO ( )/ O, enabling us to compute geodesics. Additionally, we useour algebraic description of the octahedral variety to introduce asemidefinite relaxation of projection and to place it in the contextof the more general odeco variety.In concurrent work, Golovaty et al. [2019] also represent framesby fourth-order symmetric tensors subject to some algebraic condi-tions, and they propose gradient flow on a Ginzburg-Landau–typeenergy to optimize for smooth fields.
ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020.
An additional alternative representation is proposed by Gao et al.[2017], who use quaternions to encode octahedral frames. Theirrepresentation uses relatively few values, but a matching procedurehas to be embedded in their optimization objective to account forthe fact that 48 quaternions correspond to the same frame. Theconcurrent work [Beaufort et al. 2019] uses quaternions to deriveanother parametrization of octahedral frames by points of a varietyin three complex dimensions.All previous work on 3D frame fields has only considered octa-hedral frames, which do not capture the unidirectional behaviorof frame fields near singular curves ( cf.
Figure 3). In the algebraicgeometry community, Robeva [2016] and Boralevi et al. [2017] havecompletely characterized a family of algebraic varieties of orthogo-nally decomposable (odeco) tensors. We show how the octahedralvariety is embedded in one of the odeco varieties, and we intro-duce a technique for optimization over the odeco frames. For ourpurposes, odeco frames generalize octahedral frames by allowingindependent scaling of the “axes” of a frame, including degenerationto unidirectionality at singular curves and to zero at singular nodes.Finally, Shen et al. [2016] consider fields having different localdiscrete symmetry groups, such as the tetrahedral and icosahedralgroups. They extend the methods of [Huang et al. 2011; Ray et al.2016] to compute such fields.
Relaxation of algebraic optimization problems to semidefinite pro-grams has been studied extensively in the field of real algebraicgeometry and optimization. Blekherman et al. [2012] provide anintroduction to this discipline. The efficacy of semidefinite relax-ation in computer science was demonstrated dramatically in theseminal paper of Goemans and Williamson [1995] on the maximumcut problem. Since then, semidefinite relaxation has continued to beemployed to solve both combinatorial and continuous optimizationproblems, such as angular synchronization [Singer 2011]. In geom-etry processing and graphics, this machinery has been applied tosuch problems as correspondence [Kezurer et al. 2015], consistentmapping [Huang and Guibas 2013], registration [Maron et al. 2016],camera motion estimation/calibration [Agrawal and Davis 2003;Ozyesil et al. 2015], and deformation [Kovalsky et al. 2014].Most relevant to the present work are relaxations of the Euclideanprojection problem onto a variety defined by quadratic equations, anexample of a quadratically-constrained quadratic program (QCQP).Cifuentes et al. [2017] have recently shown a stability result imply-ing that the semidefinite relaxation of Euclidean projection onto asmooth, quadratically-defined variety is exact in a neighborhood ofthe variety. Cifuentes et al. [2018] have additionally shown that theregion in which the relaxation is exact is a semialgebraic set , and theyhave provided a formula for the degree of its algebraic boundary.Unfortunately, computing this boundary is generally intractible forinteresting varieties. A deeper theoretical understanding of whensemidefinite relaxations of Euclidean projection are globally exactis still lacking.
As previously studied in [Chemin et al. 2018; Huang et al. 2011;Ray et al. 2016; Solomon et al. 2017], the basic unknown in thevolumetric frame field problem is a tuple of three mutually orthogo-nal directions at each point in a volume. These directions may berepresented by an orthonormal basis of vectors, but the signs andorder of the vectors are irrelevant thanks to octahedral symmetry.This redundancy makes detecting smoothness of a field of tuplesdifficult. Hence, it pays to use a unified representation invariant tothe symmetries of the frame.Below, we apply machinery from differential and algebraic geom-etry to derive a succinct description of this basic octahedral frameand show how it is related to rotations of the function x + x + x onthe unit sphere expressed in the spherical harmonic basis (as used torepresent frames in [Huang et al. 2011; Ray et al. 2016; Solomon et al.2017]) and to tensorial representations (as used in [Chemin et al.2018]). Algebraic language not only provides a succinct descriptionof previous representations but also suggests a means of generaliz-ing to frames whose three axes scale independently (e.g. rotationsof the function (cid:205) i λ i x i for varying λ ∈ R ), whereas previous workrequires them to have the same length ( λ = λ = λ ). This broaderset better aligns with the realities of the frame field problem, sincesingular edges have nontrivial directionality that cannot be capturedby existing representations. Relevant background material may befound in Appendix A. Intuitively, an octahedral frame field is a smooth assignment of(unoriented) orthonormal coordinate axes to each point x in a region Ω ⊂ R . Such frames are called octahedral because they exhibitsymmetries described by the octahedral group O. The space ofsuch octahedral frames can be identified with the quotient space F = SO ( )/ O — that is, the quotient of the group of orientedrotations by the right action of the octahedral group. Since O isnot a normal subgroup of SO ( ) (indeed SO ( ) is simple) F is not agroup. However, as O is a finite group acting freely on SO ( ) , F isa manifold, and the surjective map SO ( ) → F is a covering map.The universal cover of F is that of SO ( ) , namely SU ( ) , and soits fundamental group π (F ) is the lift of O to SU ( ) , the binaryoctahedral group BO. In particular, BO classifies singular curves ofoctahedral fields [Mermin 1979].SO ( ) admits a bi-invariant Riemannian metric (23), which meansthat the action of O by right-multiplication is isometric. Thus theRiemannian metric on SO ( ) descends to F , making it a Riemannianmanifold. In particular, F has geodesics that lift to geodesics ofSO ( ) . We will employ these geodesics to compute optimizationsub-steps on F ( cf. §4.1).We have introduced F as an abstract smooth manifold, but an em-bedding of F in a Euclidean space is necessary to make it amenableto computation. Previous works have introduced this equivariantembedding via the action of SO ( ) on polynomials and sphericalharmonic coefficients. For completeness, we reintroduce it here, butin a slightly different way that highlights the role of representationtheory and the orbit-stabilizer theorem. Later, we will show that Technically, their stabilizers are conjugate to O .ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020. lgebraic Representations for Volumetric Frame Fields • 16:5 the embedding of F in R is an equivariant isometric embedding ( cf. [Moore 1976]).Consider the irreducible representation ρ : SO ( ) → SO ( ) cor-responding to the fourth band of spherical harmonics. The 9 × ρ ( r ) for r ∈ SO ( ) are sometimes referred to asWigner D -matrices. Form a linear operator H ∈ R × as H = | O | (cid:213) o ∈ O ρ ( o ) . This H is a projection operator onto the subspace of R invariantunder all octahedral rotations:Lemma 3.1. For any o ′ ∈ O , ρ ( o ′ ) H = H . Proof. ρ ( o ′ ) H = | O | (cid:213) o ∈ O ρ ( o ′ ) ρ ( o ) = | O | (cid:213) o ∈ O ρ ( o ′ o ) = | O | (cid:213) o ∈ O ρ ( o ) = H . □ Corollary 3.2.
For q ∈ R , ρ ( O ) · q = q if and only if q ∈ Im ( H ) . Proof. The forward implication follows by writing q = Hq . Thereverse implication follows directly from Lemma 3.1. □ Because O is a maximal subgroup of SO ( ) , we have the follow-ing corollary. Here stab q denotes the stabilizer of q —that is, thesubgroup of SO ( ) leaving q invariant (see Appendix A).Corollary 3.3. Let q ∈ Im ( H ) and nonzero. Then stab ( q ) = O . It happens that H ∈ R × has rank one, i.e., its image is one-dimensional, motivating the following definitions. Definition 3.4.
The canonical octahedral frame is the normal-ized vector q = (cid:32) , , , , (cid:114) , , , , (cid:114) (cid:33) ⊤ ∈ R . such that H = q q ⊤ .Our canonical frame is the same as the normalized projection ofthe polynomial (cid:205) i x i into the fourth band of spherical harmonics,as in, e.g., [Solomon et al. 2017]. Definition 3.5.
The octahedral variety is the orbit of q underthe action of SO ( ) , F = ρ ( SO ( )) q = { ρ ( r ) q | r ∈ SO ( )} . To summarize, F is the orbit of q , whose stabilizer is O. A smoothversion of the orbit-stabilizer theorem [Lee 2012, Theorem 21.18]shows that there is an equivariant diffeomorphism ϕ : F → F :SO ( ) SO ( )F F ρ · q ϕ Next, we will characterize the Riemannian geometry of F byshowing that ϕ is an isometry up to a uniform scale factor. To ourknowledge, this observation has not appeared in previous work. Proposition 3.6. Let α = (cid:112) / and F α (cid:66) αF = { αq : q ∈ F } .Let π α : R × → R denote matrix multiplication by the scaledcanonical octahedral frame q α (cid:66) αq . That is, π α ( A ) = Aq α . Thenthe map π α ◦ ρ : SO ( ) → F α is a local isometry, making the induced diffeomorphism ϕ α : F → F α an isometry. Proof. Taking the differential of ρ at the identity yields the as-sociated Lie algebra representation ( Dρ ) I : so ( ) = T I SO ( ) → so ( ) ⊂ R × . ( Dρ ) I is characterized by L i (cid:66) ( Dρ ) I ( l i ) , the images of the Liealgebra generators l i (see Appendix A for definitions and supple-mental material for explicit expressions). π α is a linear map, so itsdifferential is also multiplication by q α .To see that π α ◦ ρ is a local isometry, first recall that the metric onSO ( ) is bi-invariant. In particular, for each д ∈ SO ( ) , an orthonor-mal basis for T д SO ( ) is given by the right-translated Lie algebragenerators {( D R д ) I l i } i = . So it suffices to prove that their images under π α ◦ ρ form an or-thonormal basis in T ρ ( д ) q α F α . But D ( π α ◦ ρ ) д ( D R д ) I l i = ( Dπ α )( Dρ ) д ( D R д ) I l i = ( Dπ α )( D R ρ ( д ) ) I ( Dρ ) I l i = ( Dπ α ) L i ρ ( д ) = L i ρ ( д ) q α , (1)where the second equality follows by differentiating the representa-tion property ρ ◦ R д = R ρ ( д ) ◦ ρ at the identity. So the equation weneed to prove is (cid:10) L i q , L j q (cid:11) = δ ij ∀ q ∈ F α . (2)This isometry condition can be checked explicitly at q α : (cid:10) L i q α , L j q α (cid:11) = δ ij = (cid:10) l i , l j (cid:11) . (3)For a general q = ρ ( д ) q α ∈ F α , we compute (cid:10) L i ρ ( д ) q α , L j ρ ( д ) q α (cid:11) = (cid:68) ρ ( д ) ⊤ L i ρ ( д ) q α , ρ ( д ) ⊤ L j ρ ( д ) q α (cid:69) = (cid:68) Dρ I ( д ⊤ l i д ) q α , Dρ I ( д ⊤ l j д ) q α (cid:69) . (4)The first equality follows because ρ ( д ) ∈ SO ( ) , and the second fol-lows by Lemma A.1. Let a ik = a ik ( д ) be the coefficients of theadjoint representation of SO ( ) (see Appendix A), i.e., д ⊤ l i д (cid:67) (cid:205) k a ik l k ; these coefficients form an orthogonal matrix. Now us-ing linearity of Dρ along with (3) and (4), we obtain (cid:10) L i ρ ( д ) q α , L j ρ ( д ) q α (cid:11) = (cid:42)(cid:213) r a ir L r q α , (cid:213) s a js L s q α (cid:43) = (cid:213) r , s a ir a js ⟨ L r q α , L s q α ⟩ = (cid:213) k a ik a jk = δ ij (5)as required. □ ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020.
Fig. 3. An odeco field on a triangular prism. The norm of the second bandcoefficients indicates degree of anisotropy (left). Odeco frames close to theapproximate position of the singular curve, indicated by a dot (right), scaletoward zero in the directions normal to the curve, but remain nondegeneratein the direction along the curve.
In summary, we have described the octahedral variety as an em-bedded submanifold in R isometric to F = SO ( )/ O. This isometrymeans that we can do manifold optimization over frames by work-ing in the embedding ( cf. §6.1). To show that F is really an algebraicvariety, we will need to exhibit equations that cut it out. We willdelay this until §3.2.1, when we can give a unified description of theoctahedral and odeco varieties. The previous section provides more insight into the octahedralframes used in all previous work. Octahedral frames are suitablefor representing singularity-free frame fields. However, frame fieldscommonly encountered in applications have singularities compris-ing an embedded graph. Consider the triangular prism shown inFigure 3. Near the singular curve, a smooth octahedral field wouldrotate infinitely quickly, and it would not have a well-defined valuealong the curve. This issue is analogous to the case of unit crossfields—note that for cross fields, the hairy ball theorem requires singularities on simply-connected surfaces.One solution to this is to replace the hard constraint that the fieldvalues lie in the octahedral variety with a penalty term, as motivatesthe MBO method for octahedral fields detailed below (§6.2). An-other solution is homogenization—i.e., allowing field values to scale,replacing singularities by zeros, as is considered for cross fieldsin [Knöppel et al. 2013]. But consider the triangular prism again:a scaled octahedral field would vanish completely at the singularcurve since all three orthogonal axes must scale uniformly. Thismakes the octahedral representation unable to capture the align-ment of the field to the singular curve. To cope with this problemand to show the value of our algebraic approach, we now describe asuperset of the octahedral frames. This set allows the axes to scaleindependently; for instance, as shown in Figure 3, the frame axesorthogonal to the singular curve scale toward zero while the axistangent to the singular curve remains nonzero.The symmetric orthogonally decomposable ( odeco ) tensor vari-eties, introduced in [Robeva 2016], parametrize symmetric tensors T ∈ Sym d R n that can be written T = n (cid:213) i = λ i ( v i ) ⊗ d for some set of n orthonormal vectors v i ∈ R n , where v ⊗ d denotesthe d -wise tensor power of the vector v . That is, an odeco tensorencodes a set of orthogonal vectors up to permutation. If d is even, T is also invariant under sign changes to the v i . Moreover, an odecotensor assigns weights λ i to the vectors v i . This property allows usto encode frames whose axes scale independently.There is a one-to-one correspondence between symmetric tensorsof order d over R n and homogeneous polynomials of degree d in n variables ( cf [Robeva 2016, §1.2]). This correspondence is given inone direction by taking a polynomial p ∈ R [ x , . . . , x n ] to its tensorof d th derivatives, and in the other direction by symbolic evaluationon the vector of formal variables x = ( x , . . . , x n ) : T ∈ Sym d R n (cid:55)→ p T ( x ) = d ! T ( x , . . . , x ) ∈ R [ x ] . This is a generalization of the correspondence between symmetricbilinear forms (i.e. symmetric 2-tensors) and quadratic forms (i.e.,homogeneous quadratic polynomials). Note that T is odeco if andonly if p T can be written as a sum of d th powers of linear forms p T ( x ) = (cid:213) i λ i ( v ⊤ i x ) d where the v i are orthonormal as above. In this case, we also refer to p T as an odeco polynomial.The defining equations of the odeco varieties are quadrics—homogeneousquadratic equations—in the tensor coefficients [Boralevi et al. 2017,Theorem 4] or equivalently in the coefficients of the associatedpolynomial p T . That is, a homogeneous polynomial p ( x ) = (cid:213) d + ··· + d n = d u d ,..., d n x d . . . x d n n is odeco if and only if the coefficients u satisfy u ⊤ A i u = A i . In the case relevant to us,where n = d =
4, there are exactly 27 such defining equations.While dimension counting might suggest otherwise, these equationsare not redundant—as can be seen by computing a Gröbner basis ofthe ideal they generate. The matrices A i are listed explicitly in thesupplemental document. We will henceforth refer to this particularodeco variety simply as the odeco variety ˜ F . Figure 4 plots severalfourth-order odeco polynomials over the unit sphere. Fig. 4. Examples of z -aligned odeco polynomials, plotted over the sphere. ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020. lgebraic Representations for Volumetric Frame Fields • 16:7
Recallthat our octahedral frames were represented by coefficients in anirreducible representation of SO ( ) , while the odeco variety was de-fined using monomial or tensor coefficients. To see the relationshipbetween the two varieties, it is beneficial to recast the odeco varietyin the irreducible representation basis. This corresponds to lookingat the coefficients of odeco polynomials in the basis of sphericalharmonics.The quartic polynomials comprising the odeco variety decomposeas linear combinations of the spherical harmonics of bands 0, 2, and4. Consider an odeco polynomial (cid:205) i = λ i ( v ⊤ i x i ) , and let q = ( q , q , q ) ∈ V × V × V be its coefficients in this basis of even spherical harmonics. These 15coefficients give us a different representation of odeco frames in R ,where each band has a clear meaning. In particular, q = C (cid:205) i λ i ,where C is a constant, and similarly ∥ q ∥ = C (cid:169)(cid:173)(cid:171)(cid:213) i λ i − (cid:213) i < j λ i λ j (cid:170)(cid:174)(cid:172) . The parenthesized expression is the squared distance between ( λ , λ , λ ) and the line λ = λ = λ . Thus q = q is a scalarmultiple of an octahedral frame. The set { q ∈ ˜ F | q = } consists of scalar multiples of the octahedral variety indexed by q . Fix q = q be a constant such that ∥ q ∥ =
1. Theoctahedral variety is the intersection of this affine subspace withthe odeco variety. Reducing the equations (6) of the odeco varietywith respect to this subspace yields 15 inhomogeneous quadraticequations (cid:18) q (cid:19) ⊤ P i (cid:18) q (cid:19) = , i = , . . . ,
15 (7)cutting out the octahedral variety F . As was the case for the odecovariety, these equations are not redundant. The symmetric matrices P i are listed explicitly in the supplemental document. We have introduced two spaces, the octahedral and odeco varieties,that can serve as target sets for frame fields. To compute such fields,we will have to solve optimization problems over products of manycopies of these varieties. Naïvely, one might plug the quadratic con-straints in (6) and (7) directly into a generic quadratic optimizationsolver. However, the A i or P i are not positive semidefinite, nor cantheir span be rewritten as the span of positive semidefinite matrices.So the constraints are nonconvex and challenging to enforce.As an alternative, we use optimization algorithms that are tailoredto the manifold-valued variables in our problem. These algorithmsemploy sub-steps such as geodesic traversal and projection. Wederive these operations for a single frame below. By Proposition 3.6, the scaled octahedral variety F α is locally isomet-ric to SO ( ) via the map r (cid:55)→ ρ ( r ) · q α . It follows that geodesics ofSO ( ) push forward to geodesics of F α . The relation (17) in Appendix A then allows us to compute geodesics on F α in closed form withoutevaluating the representation map ρ explicitly. Let v ∈ T q F α . v canbe written in a basis induced by the SO ( ) action, v = (cid:205) i = v i L i q ,Here the coefficient vector v = ( v , v , v ) ∈ R is the “axis-angle”representation of a rotation, and the SO ( ) exponential maps it tothe corresponding rotation matrix. This can be computed by con-jugation with a rotation r taking the unit vector v /∥ v ∥ to ( , , ) .Composing with the representation map, we haveexp q ( v ) = ρ ( r ⊤ exp (∥ v ∥ l ) r ) q = ρ ( r ) ⊤ exp (∥ v ∥ L ) ρ ( r ) q , (8)where (unsubscripted) exp denotes the ordinary matrix exponential.To compute ρ ( r ) , define spherical coordinates θ , φ such that v = ∥ v ∥( cos θ sin φ , sin θ sin φ , cos φ ) . Then one possible choice for r is r = exp (− ϕl ) exp (− θl ) = r ⊤ exp (− ϕl ) r exp (− θl ) , where r = exp (( π / ) l ) . So ρ ( r ) = R ⊤ exp (− ϕL ) R exp (− θL ) , (9)where R = exp (( π / ) L ) . Combining (9) with (8), we can computegeodesics by products of two simple ingredients: R and exp ( tL ) for t ∈ R . Closed-form expressions for both appear in the supplemen-tal document, §2. Note that a formula similar to (9) is common in thegraphics literature, e.g., for rotating BRDFs expressed in sphericalharmonic coefficients ( cf. [Kautz et al. 2002], Appendix). F and ˜ F are both varieties defined by quadratic equations. F ⊂ R iscut out by 15 inhomogeneous quadratic equations (7), while ˜ F ⊂ R is cut out by 27 homogeneous quadratic equations (6). Consider theproblem of finding the closest point in F to a given point y ∈ R : Π F ( y ) = arg min q ∈ F ∥ q − y ∥ . (P F )This problem is called Euclidean projection onto a quadratic variety.It is an example of a QCQP, and the general recipe for semidefiniterelaxation detailed in §A.2 automatically applies. The SDP will havethe form arg min Q ∈ R × ⟨ Y , Q ⟩ subject to Q = ⟨ P i , Q ⟩ = , i = , . . . , Q ⪰ , (SDP F )where P i are the symmetric matrices from (7), and Y = (cid:18) ∥ y ∥ − y ⊤ − y I × (cid:19) . Let Q ∗ be an optimal solution to (SDP F ). Since (SDP F ) is a relax-ation of (P F ), ⟨ Y , Q ∗ ⟩ is a lower bound on the objective value of(P F ). If it happens that rank ( Q ∗ ) =
1, then Q ∗ can be written in theform Q ∗ = (cid:18) q ∗ (cid:19) (cid:18) q ∗ (cid:19) ⊤ , ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020. from which it follows that (cid:18) q ∗ (cid:19) ⊤ P i (cid:18) q ∗ (cid:19) = , i = , . . . , , i.e., q ∗ ∈ F . In this case, q ∗ is the globally optimal solution to (P F ).This situation is called exact recovery. The foregoing discussion alsoapplies, mutatis mutandis , to ˜ F .In our MBO algorithm presented below ( cf. §6.2), we alternateprojection with stepping off the variety. We refer to the followingtheorem, which suggests that when taking small enough steps fromsmooth points of a variety, projection will be generically exact.Theorem 4.1 (Cifuentes et al. [2017], Theorem 1.2). Considerthe Euclidean projection problem arg min q ∈V ∥ q − y ∥ , y ∈ R n , where V ⊂ R n is a real variety defined by quadratic equations f ( q ) = · · · = f m ( q ) = . Let ¯ y ∈ V be a point at which the rank of theJacobian ∇ f ( ¯ y ) is equal to the codimension of V . Then the semidefiniterelaxation of projection is exact for y ∈ R n sufficiently close to ¯ y . In addition to this theoretical motivation, we have ample em-pirical evidence that our relaxations are exact under much moregeneral conditions. The blue histogram in Figure 5 shows the resultsof attempting to project 10 random points onto F via semidefiniterelaxation. We use the ratio of the second largest eigenvalue λ ( Q ∗ ) to the largest eigenvalue λ ( Q ∗ ) as a proxy for exactness, as it in-dicates whether Q ∗ was rank-one up to machine precision. For alloctahedral projections, λ ( Q ∗ )/ λ ( Q ∗ ) was ≈ − or less. Based onthese results, we make the following conjecture.Conjecture 4.2. For a generic point q ∈ R , the solution to (SDP F ) has rank , and therefore yields the exact projection Π F ( q ) . We also tried projecting 10 random quartic polynomials ontothe odeco variety, as shown in the red histogram in Figure 5. Thevast majority of odeco projections were also exact up to numericalprecision: out of 10 solutions, only 60 had an eigenvalue ratiogreater than 10 − . Theorem 4.1 gives us some intuition as to whythis might happen. The stability result only holds near smoothpoints of the variety, but whereas the octahedral variety is a smoothmanifold, the odeco variety has a singularity at the origin, separatingpolynomials of different signs.Conjecture 4.3. For a generic point q ∈ R representing a positive polynomial, the SDP relaxation yields the exact projection Π ˜ F ( q ) . Indeed, the green histogram in Figure 5 shows the results of odecoprojection on 10 random positive quartic polynomials, generatedas sums of squares of random quadratic polynomials. For all suchpositive initial points, the SDP solution had λ ( Q ∗ )/ λ ( Q ∗ ) < − ,supporting Conjecture 4.3.For octahedral frames, we also compare our SDP-based projectionto the previous state of the art method proposed by Ray et al. [2016].Because Ray et al.’s method is based on nonconvex optimization,we would expect it to get stuck in local minima. Indeed, out of100 000 trials, we found 600 cases for which the result of Ray et al.’smethod was at least 10 − further from the initial point than our −
16 10 −
13 10 −
10 10 − − − . . . . λ / λ P r o b a b i l i t y D e n s i t y OdecoOctahedralOdeco (positive initial polynomial)
Fig. 5. Histogram of eigenvalue ratio λ ( Q ∗ )/ λ ( Q ∗ ) for solutions to theSDP relaxations of Euclidean projection onto F and ˜ F . Projections of random points were tested for each. The maximum ratio for octahedralprojection was . × − . The maximum ratio for odeco projection ofpositive quartic polynomials was . × − . See §4.2. d Ray et al. = .
963 1.07 1.08 1.07 d Ours = .
496 0.850 0.857 0.871
Fig. 6. For some query polynomials, our globally optimal SDP-based projec-tion onto the octahedral variety (blue) yields dramatically different resultsfrom Ray et al.’s approximate projection (red). Our projections are closerto the query points y (distances shown in blue) compared to Ray et al.’s(distances in red). The query polynomials are plotted on the sphere, withcolor and distance from the center proportional to magnitude. projection—a nearly 1% error rate. Moreover, the difference betweenthe computed projections can be substantial, as illustrated in Figure6. The overall aspiration of this work is to compute smooth frame fieldsin a region Ω ⊆ R aligned to the normal n on its boundary ∂ Ω .We assume Ω to be compact with ∂ Ω a union of smooth, embeddedsurfaces. We think of a frame field as a map ϕ : Ω → V into somespace of frames V satisfying alignment boundary conditions.We’ve examined the geometry of two candidates for the fiber V —the octahedral variety F and the odeco variety ˜ F . It remains todescribe the boundary conditions and to define an energy we want ϕ to minimize. The term “fiber” is intended to be suggestive. It may be fruitful to consider a bundle π : P → Ω with fiber V , whose restriction ∂ P → ∂ Ω is nontrivial and encodes theboundary alignment condition. The map ϕ would be replaced with a section of P .ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020. lgebraic Representations for Volumetric Frame Fields • 16:9 In either F and ˜ F , the frames aligned to a given direction n ∈ R formthe intersection of an affine subspace with the respective variety. Ineach case, this intersection is a lower-dimensional variety definedby a different set of quadratic equations. We can impose alignmentboundary conditions by working over this lower-dimensional vari-ety. For boundaries with sharp creases, we can optionally excludecreased vertices, where the normal direction is ill-defined, fromthe alignment constraint. Other constraints can also be imposed atcreases—for example, alignment to the crease tangent. Let ρ : SO ( ) → SO ( ) bethe irreducible representation as in §3. First consider the octahedralframes aligned to the vertical z = ( , , ) ⊤ . The space of vertical-aligned frames F z (cid:66) ρ ( A z ) q , where A z is the coaxial subgroupconsisting of all rotations about z , i.e., A z = { exp ( tl ) | t ∈ R } , and q is the canonical octahedral frame. As derived in [Huang et al.2011; Ray et al. 2016], such frames have the formexp ( tL ) q = (cid:32)(cid:114)
512 sin ( t ) , , , , (cid:114) , , , , (cid:114)
512 cos ( t ) (cid:33) ⊤ = q z + Bs z , (10)where q z = ( , , , , (cid:112) / , , , , ) ⊤ , s z = ( cos ( t ) , sin ( t )) ⊤ , and B z is the obvious 9 × q ∈ R , the closest vertical-aligned frame is Π F z ( q ) = q z + B z B ⊤ z q ∥ B ⊤ z q ∥ . For a general unit normal n , the aligned frames can be written ρ ( r n )( q z + B z s z ) , (11)where r n ∈ SO ( ) is any rotation taking z to n . The projection of q ∈ R onto the n -aligned frames is then Π F n = ρ ( r n ) Π F z ( ρ ( r n ) ⊤ q ) (12)During optimization, this projection is used for boundary-alignedframes. Consider the standard odecoframe (cid:205) i λ i x i rotated by some r ∈ A z , and fix λ =
1. As describedin §3.2.1, it has coefficients in the even-numbered irreducible repre-sentations (bands of spherical harmonics) q = ( q , q , q ) ∈ V × V × V = R × R × R . Just as in the octahedral case, q can be decomposed into a fixed partand a rotating part parametrized by a lower-dimensional vector q = q z + B z s z , (13)where now s z ∈ R and B z ∈ R × . The equations (6) reduceto three quadratic equations in the coefficients of s z , which canbe used to construct a semidefinite program to project onto thevertical-aligned odeco frames, as in §4.2. The details are given inthe supplementary material. For frames aligned to a direction n , q = ρ ( r n )( q z + B z s z ) , where ρ is now the direct sum representation on V × V × V .Projection onto the n -aligned odeco frames can be constructed fromprojection onto the z -aligned odeco frames as in (12). As we are searching for smoothest fields, a natural choice for theenergy is Dirichlet energy ∫ Ω ∥∇ ϕ ∥ dx , where the norm is takenwith respect to a metric on V . There are multiple problems thisformulation will need to confront. As we have seen, F is a smoothmanifold. But being a quotient of the 3-sphere, it has positive cur-vature, making mere existence of harmonic maps Ω → F a hardproblem. Even more fundamentally, it is not clear how to representsingularities of ϕ —the map cannot be consistently defined alongsingular curves. ˜ F attempts to solve this problem by representingsingular frames with only one direction, while allowing the othertwo to vanish. Along a singular curve, we would expect an odecofield to take rank-one values of the form λv ⊗ , where v is tangentto the singular curve. Let T = ( V , T ) be a tetrahedral mesh of Ω with vertices V andtetrahedra T . The set of boundary vertices will be denoted by ∂ V . Adiscrete octahedral frame field on T is specified by a map q : V →V = F or ˜ F , giving a frame q i for each vertex i . As the octahedralvariety F ⊂ R and the odeco variety ˜ F ⊂ R , we can think of q asa d × n matrix, where n = | V | and d = V .This is equivalent to measuring L distance between the correspond-ing polynomials over the sphere, as employed in [Huang et al. 2011;Ray et al. 2016; Solomon et al. 2017]. Note that this would not bethe case if we compared odeco frames in the monomial basis. Thediscrete energy is E ( q ) =
12 tr (cid:16) qLq ⊤ (cid:17) , where L is the linear finite element Laplacian on T , with the usualcotangent weights. We now have two variety-constrained optimization problems of theform min 12 tr (cid:16) qLq ⊤ (cid:17) s.t. q i ∈ V , ∀ i ∈ Vq i ∈ V n i , ∀ i ∈ ∂ V , (14)where the variety V is either F or ˜ F , q i are the columns of q , ∂ V denotes the set of boundary vertices, and V n i is the variety of framesaligned to the normal at boundary vertex i . In the case V = F , (14) becomes a manifold-constrained optimization problem. The oc-tahedral variety is a Riemannian manifold,and we can compute geodesics on it as de-scribed in §4.1. So we can apply the standard ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020.
Riemannian trust-region (RTR) algorithm [Absil et al. 2007] as im-plemented in the Manopt toolbox for MATLAB [Boumal et al. 2014].We consider the field to be a point in the product manifold q ∈ (cid:214) i ∈ V \ ∂ V F × (cid:214) i ∈ ∂ V F n i . In addition to the geodesics computed in §4.1, we also need a wayto project vectors (e.g. gradients) in the ambient space R into thetangent space at a point q , proj q : R → T q F . The vectors { L i q } i = form an orthogonal basis for this tangent space ( cf. Proposition 3.6),so projection is just given by taking the inner product with each ofthese vectors.The odeco variety ˜ F is not a manifold, but it is smooth almosteverywhere. We know of no way to compute exact geodesics on˜ F , but we can implement a version of RTR by replacing the exactexponential map with a retraction , i.e., a first-order approximation( cf. [Absil et al. 2007, Definition 2.1]). In doing so, we give up thesuperlinear local convergence guarantees of RTR ( cf. [Absil et al.2007, Theorem 4.12]) but retain a global convergence guarantee. Inpractice, we find that odeco RTR converges at a reasonable rate.At a smooth point of ˜ F , it is easy to project vectors onto thetangent space using the fact that ˜ F is defined by quadratic equations.Let q ∈ ˜ F , and assume the coefficients of q are expressed in sphericalharmonic coefficients. Then differentiating the equations q ⊤ ˜ A i q = q is N q ˜ F = span { ˜ A i q } , where ˜ A i are expressed in the spherical harmonic basis. Then T q ˜ F = ( N q ˜ F ) ⊥ , where the orthogonal complement is taken with respect tothe standard metric under which the spherical harmonic functionsare orthonormal.We compute retractions as follows. The tangent space to theodeco variety decomposes into a rotational part and a scaling part : T rq ˜ F = span { ˜ L i q } i = , T sq ˜ F = ( T rq ˜ F ) ⊥ , where the orthogonal complement is taken within T q ˜ F . We then setretr q ( v ) = e v r · ˜ L ( q + v s ) , where v s is the component of v in T sq ˜ F , v r is the component of v in T rq ˜ F , and v r · ˜ L (cid:66) (cid:213) j ( v r ) j ˜ L j , where v r = (cid:205) j ( v r ) j ˜ L j q . It is simple to verify that retr q ( v ) ∈ ˜ F andthat ∂∂ t retr q ( tv ) (cid:12)(cid:12)(cid:12)(cid:12) t = = v . We compare the convergence behavior of octahedral RTR to thatof the method of Ray et al. [2016] under two sets of conditions—Rayet al.’s initialization and combinatorial Laplacian (Figure 7); and ourrandom initialization and linear finite element Laplacian (Figure 8).The quadratic local convergence of RTR stands in stark contrast tothe slower linear convergence behavior of Ray et al.’s method. RTRconverges faster but finds solutions of similar Dirichlet energy toprevious work (see Table 1 in the supplemental document). − − − − Time (s) G r a d i e n t N o r m Ray et al. (Bone)RTR (Bone)Ray et al. (Cup)RTR (Cup)Ray et al. (Rockerarm)RTR (Rockerarm)
Fig. 7. Convergence behavior of octahedral RTR and the method of Rayet al. [2016] on various models, starting from Ray et al.’s initialization andusing the combinatorial Laplacian.
10 20 30 40 50 60 7010 − − − − Time (s) G r a d i e n t N o r m Ray et al. (Bone)RTR (Bone)Ray et al. (Cup)RTR (Cup)Ray et al. (Rockerarm)RTR (Rockerarm)
Fig. 8. Convergence behavior of octahedral RTR and the method of Rayet al. [2016] on various models, starting from random initialization andusing the geometric Laplacian.
As we have observed, it is possible to do optimization on the octa-hedral and odeco varieties by moving along curves that stay on thevarieties exactly. However, this hard constraint sometimes causespure manifold optimization to get stuck in local minima. To avoidsuch local minima, an approach that allows “tunneling” is required.The method of Merriman, Bence, and Os-her [1992] is a diffusion-generated methodfor computing (possibly singular) maps intoa target space embedded in a Euclideanspace. Following [Osting and Wang 2017;
ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020. lgebraic Representations for Volumetric Frame Fields • 16:11
Viertel and Osting 2019], we first replace the hard constraint thatthe field values lie in the variety with a penalty term, yielding anenergy of Ginzburg-Landau type, E ϵ ( q ) = ∫ Ω ∥∇ q ∥ dx + ϵ ∫ Ω dist ( q , V) dx , (15)where V is our variety. Consider this energy in the limit ϵ → E ϵ . Gradient flow on thefirst term—Dirichlet energy—yields componentwise heat diffusion, ∂ t q ( x , t ) = ∆ q ( x , t ) q ( x , ) = q k − . (16)with the boundary constraints given in §5.1. In practice, we doone step of implicit Euler integration per overall algorithm step.Gradient flow on the second term of (15) in the limit ϵ → Algorithm 1:
MBO over a variety
V ⊂ R d input : initial d × n field q , diffusion parameter τ , optimization schedule β ( k ) result : q k k ← repeat τ k ← β ( k ) τ Diffusion step:
Solve the linear system ( M − τ k L ) q ⊤ k = Mq ⊤ k − withcolumns ( q k ) i constrained to be in the affine span (11) or (13) foreach i ∈ ∂ V . Projection step:
Project q k into the variety: ( q k ) i ← (cid:40) Π V (( q k ) i ) i (cid:60) ∂ V Π V n i (( q k ) i ) i ∈ ∂ V . ∆ k ← tr (( q k − q k − ) M ( q k − q k − ) ⊤ ) tr ( q k Mq ⊤ k ) E k ← tr ( q k Lq ⊤ k ) ∆ E k ← E k − E k − k ← k + until ∆ E k / E k < δ or ∆ k < δ We follow the suggestion of [van Gennip et al. 2014; Vierteland Osting 2019] to set τ relative to the inverse of the smallestnonzero eigenvalue of the Laplacian. In ordinary MBO, β ( k ) = τ does not change over the course of the algorithm. In Figure9, we test the effects of different values of τ . As we decrease τ ,the field is able to relax more and reduce the Dirichlet energy. Onthe other hand, larger values of τ allow more tunneling so thatthe field can escape local minima. This suggests a modified MBOscheme (mMBO) in which we start with a large τ and progressivelyreduce it over the course of optimization—analogously to whathappens to the temperature in simulated annealing. Our mMBOuses an optimization schedule β ( k ) = k − . This optimizationschedule starts with a large diffusion time for robustness to randominitialization, then sweeps τ through various orders of magnitudequickly. We have found that this heuristic produces a good balancebetween quick convergence and field quality. τ ≈ . τ ≈ . τ ≈ . + RTR E = . E = . E = . Fig. 9. Results of octahedral MBO on a torus. With a relatively large dif-fusion time τ , MBO produces a field with tightly packed singular regions(left). At a smaller value of τ , singular curves relax toward the boundary,reducing the Dirichlet energy (center). RTR reduces the energy even further(right), pushing the singular curves to the boundary.
20 40 60 80 1003 , , , , Iteration D i r i c h l e t E n e r g y E MBO ( τ ≈ τ ≈ τ ≈ . Fig. 10. Energy convergence on the rockerarm_91k model. MBO’s conver-gence to lower energies is limited by the diffusion time τ . Decreasing τ according to the optimization schedule used in mMBO achieves the bestresults overall. Figure 10 depicts the convergence behavior of RTR, MBO withvarious τ values, and mMBO starting from a projection of an ap-proximately harmonic field. While the energy value achieved byordinary MBO is limited for each τ , mMBO achieves a lower valueby sweeping through various values of τ . Note that the iterationcounts shown do not reflect wall-clock time; in particular, RTR runsat least an order of magnitude faster than the other algorithms. We initialize our algorithms with vertexwise random octahedralfields, generated by—separately for each vertex—starting with thecanonical frame, rotating it by a random angle between 0 and 2 π about the positive z -axis, and then rotating the positive z -axis to arandom direction. We do this even for optimization of odeco fields ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020. to avoid encountering odeco frames that have negative weights λ i ( cf. Conjecture 4.3). When starting from octahedral initialization,the odeco frames we compute always have nonnegative weights.
Fig. 11. Robustness to random initialization. Odeco MBO transforms ver-texwise random initial fields into qualitatively identical results up to thesymmetry of the sphere. All results were computed on a sphere with vertices. Integral curves are depicted in regions of rapid rotation, i.e., inproximity to singularities.
In Figure 11, we demonstrate the robustness of MBO to randominitialization. On the sphere, global rotations of a field do not changethe objective value. Initializing odeco MBO with different fields ofvertexwise random frames, we find that it converges to randomlyrotated copies of a qualitatively identical, symmetric field.Figure 12 illustrates the behavior of octahedral and odeco fields asthe density of the underlying tetrahedral mesh increases. Due to theunit-norm constraint, the gradient of an octahedral field becomesunbounded close to its singularities. This leads to logarithmic diver-gence of the total energy as the tet mesh becomes finer. In contrast,the additional scaling degrees of freedom available to odeco fieldsallow renormalization of singularities, as shown in Figure 13. Thus,odeco field energy plateaus as mesh density increases. Note alsothe much smaller variance in energy between runs for odeco fields,quantitatively illustrating robustness to random initialization.Table 1 (supplemental document) compares timings and energyvalues for our methods and the method of Ray et al. [2016], com-prising 18 types of experiments on 15 different models, for 270 totalruns. Direct comparison of energy values with Ray et al. [2016]is not possible, as their method uses the graph Laplacian and ini-tializes by solving a linear system, whereas our method uses thegeometric (finite element) Laplacian and random initialization. Toattempt a fair comparison, we report results of their method alone,their method followed by RTR with the geometric Laplacian, and their method substituted with the geometric Laplacian and randominitialization. All experiments in the table were conducted on anUbuntu workstation with a four-core, 3.6 GHz Intel Core i7-7700and 32 GB of RAM. A MATLAB implementation of our algorithmsaccompanies this paper and also includes our implementation of[Ray et al. 2016].The results in Table 1 show that RTR converges very quickly but,like previous work, easily gets stuck in local minima. In contrast,mMBO takes longer to converge but produces higher-quality fieldsthat reflect the symmetries of the volume. Running RTR to polishthe results of mMBO produces the best energies.Despite sometimes getting stuck in local minima (see Figure 6), weobserve that Ray et al.’s approximation of projection can be usefulin practice. Table 1 (supplemental document) includes experimentswith MBO and mMBO (see §6.2) substituting Ray et al.’s projection − − Tet Volume (Geometric Mean) E n e r g y OctahedralOdeco
Fig. 12. Energy density diverges for octahedral fields, but plateaus for odecofields, as mesh density increases. Also notice the higher variance for oc-tahedral fields. Ten fields of each type were computed by RTR on each ofthirteen tetrahedral meshes of the unit sphere of various densities.Fig. 13. The energy density at singularities of an octahedral field dominatesthe total energy (left). In contrast, an odeco field’s energy is distributedmore uniformly (right). Results computed using MBO + RTR.
ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020. lgebraic Representations for Volumetric Frame Fields • 16:13 [ R a y e t a l . ] m M B O + R T R Fig. 14. Twisted singular curves on the handle of the cup in Ray et al.’sresult lead to twisting in the resulting mesh (top). Our field yields a meshwith a more regular structure throughout (bottom). [Ray et al. 2016] mMBO + RTR
Fig. 15. On the bone model, twisted singular curves in a field computedby the method of Ray et al. [2016] lead to a degenerate integer grid map,producing highly-distorted hexahedra (left). Our field yields a mesh withoutthis degeneracy (right). for the true projection. In most cases, the resulting energy is verysimilar, suggesting that the diffusion iterations in MBO smooth outany errors resulting from incorrect projection. This hybrid MBOcan be a useful tradeoff between correctness and speed.In Figures 14–17, we compare fields computed by octahedralmMBO + RTR to those computed by the method of Ray et al. [2016].Our fields not only have lower Dirichlet energy, but also bettersingular structures. To visualize singular structures, we use the visu-alization technique of Liu et al. [2018]. To illustrate the importanceof singular structure, we have generated hexahedral meshes fromboth sets of fields, following the methods laid out in [Nieser et al.2011] and [Lyon et al. 2016]. The meshes are computed from theraw field data, with no preprocessing other than tetrahedral meshrefinement to resolve and localize singular curves. In particular,we do not “correct” singularities—thus, both sets of meshes showsome degeneracies resulting from collapsed or flipped tetrahedraduring the parametrization step. However, our fields yield mesheswith fewer, smaller degeneracies, less distortion, and more regularstructure. [ R a y e t a l . ] m M B O + R T R Fig. 16. Due to the twisting of singular curves in the field generated byRay et al.’s method (top), the mesh is highly distorted on the right arm. Incontrast, our field (bottom) has a regular cubic singular structure, leadingto fewer mesh degeneracies. [ R a y e t a l . ] m M B O + R T R Fig. 17. Note the simpler, more regular singular curves in our result ascompared to the result of Ray et al. [2016]—especially on the bunny’s nose.The degeneracy in Ray et al.’s result leads to a collapse of the integer gridmap on the head.
Figures 18 and 19 compare fields produced by mMBO + RTR tofields produced by the method of Gao et al. [2017]. Our fields betterrespect the symmetries of the underlying models. We note thatmMBO + RTR yields a higher-quality result without requiring amultiscale method like the one Gao et al. [2017] employ.
ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020. E = E = . Fig. 18. Comparison of octahedral field algorithms on the gear model. Afield produced with mMBO + RTR (right) exhibits lower energy and greatersymmetry than one produced by the method of Gao et al. [2017] (left). E = . E = . Fig. 19. Comparison of octahedral field algorithms on the space-filling torusmodel. A field produced with mMBO + RTR (right) exhibits lower energyand greater symmetry than one produced by the method of Gao et al. [2017](left).
A stronger understanding of the unknowns in the volumetric framefield problem enables both theoretical and practical developments.On the theoretical side, our reexamination of the typical representa-tion of octahedral frames not only yields useful geodesic formulasand projection operators but also suggests generalization to odecoframes. For both frame representations, optimization algorithmsdesigned explicitly to traverse the manifold of unknowns yield gainsin efficiency and in the quality of computed results.Our work is intended not only to improve on existing frame fieldcomputation algorithms but also to inspire additional research intothe structure of spaces of frames and to highlight the significance ofthis structure to practical aspects of field computation and meshing.Most immediately, a few open questions are left from our discus-sion. On the theoretical side, conjectures 4.2 and 4.3 remain to beproven. We anticipate that both can be embedded in a more generalframework explaining when relaxations of projection problems areexact. On the algorithmic side, RTR seems to get stuck in local min-ima much more often than MBO, especially on denser meshes. We [ R a y e t a l . ] m M B O + R T R Fig. 20. A field of lower Dirichlet energy (bottom) may still result in moremesh degeneracies than a field of higher Dirichlet energy (top) due to topo-logical impediments to meshability, namely the presence of an additionalvalence – junction. hypothesize that this is because RTR strictly moves along the man-ifold, while MBO is free to tunnel through the ambient space. Onthe other hand, RTR converges much faster and yields high-qualityfields on coarse meshes. Given these observations, we anticipatethat it may be possible to incorporate RTR into a multiscale methodthat leverages its efficiency while avoiding local minima that appearat the finest scales.As with most existing frame field computation algorithms, evenwhen mMBO+RTR converges to a smooth field, the topology of thefield is not always hex-meshable. Figure 20 shows a case where themethod of Ray et al. [2016] yields a field of higher Dirichlet energy,but which leads to fewer degeneracies than our field. In particular,the presence of a singular curve whose type changes from valence3 to 5 leads to a degeneracy in the integer-grid map. While ourmethod succeeds at finding fields of lower Dirichlet energy—theobjective function of our method and that of Ray et al. [2016]—minimizing this energy is an incomplete proxy for our ultimategoal, namely to obtain smooth, meshable fields. We would like toinvestigate additional metrics, e.g., integrability of frame fields, thatmight make it possible to express meshability constraints rigorously.To define such additional metrics, it might be fruitful to considerfurther alternative frame field representations. For example, givenour use of Lie algebra representations, a logical next step mightbe to introduce an SO ( ) –principal bundle and work with connec-tions on that bundle as variables. In this theory, quantities such ascurvature, torsion, or the Chern-Simons functional might encodeimportant features of frame fields. The discretization of connectionsin [Corman and Crane 2019] represents a promising first step.Finally, while our algorithms do not explicitly take account ofsymmetry, we find that fields computed by MBO consistently re-produce the symmetries of the volume. It would be interesting todevelop a better theoretical understanding of this behavior and todevelop machinery for explicitly promoting conformation to near-symmetry in deformed domains.These and other challenges appear when extending well-knownmachinery from geometry processing on surfaces to volumetricdomains, as demanded by applications in engineering, simulation, ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020. lgebraic Representations for Volumetric Frame Fields • 16:15 and other disciplines. Over the longer term, careful considerationof problems like the ones we consider in this paper will lead to aversatile and generalizable approach to geometry processing.
A MATHEMATICAL BACKGROUND
Representation theory and algebraic geometry provide a conciselanguage for the discussion of the octahedral and odeco varietiesabove. We begin with definitions of a few terms used throughoutthe paper. While a detailed introduction to either of these subjects isoutside the scope of our discussion, we refer the reader to [Stillwell2008] for an introduction to matrix Lie groups and to [Blekhermanet al. 2012] for a primer on real algebraic geometry and optimization.
A.1 Lie Groups and Representations A Lie group is a group G that is also a smooth manifold, and suchthat the multiplication · : G × G → G and inversion − : G → G aresmooth maps. A matrix Lie group is a subgroup of the invertiblematrices GL ( n ) ⊂ R n × n that is simultaneously a submanifold. Exam-ples include the orthogonal groups O ( n ) and the special orthogonalgroups SO ( n ) . The primary Lie group of interest to us in this paperwill be SO ( ) .For a matrix Lie group G ⊂ R n × n , the Lie algebra is a linearsubspace of R n × n identified with the tangent space to G at theidentity I , denoted g = T I G . The inner automorphism Ad д : G → G taking h (cid:55)→ дhд − preserves the identity, and its derivative at theidentity induces the Lie bracket [· , ·] : g × g → g . In the case of amatrix Lie group, this is just the matrix commutator.As left translation L д : h (cid:55)→ дh is a diffeomorphism, every el-ement of g is also associated to a left-invariant vector field on G .The exponential map exp : g → G is given by starting at I andintegrating this vector field for one unit of time. For a matrix Liegroup, exp is the ordinary matrix exponential.A representation of a matrix Lie group G is a smooth grouphomomorphism ρ : G → GL ( n ) for some n . A representationof G induces a representation of each subgroup H ⊂ G . A Liegroup representation comes with a Lie algebra representation ( Dρ ) I : g → R n × n , which preserves the Lie bracket. Moreover,representations commute with exponentials: ρ ◦ exp = exp ◦( Dρ ) I . (17)An irreducible representation is one that cannot be decomposed asa direct sum of sub-representations. The irreducible representationsof compact Lie groups such as SO ( ) are finite-dimensional. The adjoint representation is the irreducible representation of G onits own Lie algebra given by conjugation:Ad : G → Aut ( g ) Ad ( д ) a (cid:66) д − aд . (18)Lemma A.1. Ad commutes with representations of G : ( Dρ ) I ( д − aд ) = ρ ( д ) − (( Dρ ) I a ) ρ ( д ) . (19)Proof. Using (17) and conjugating by д , we obtain ρ ( д − exp ( ta ) д ) = ρ ( д ) − exp ( t ( Dρ ) I a ) ρ ( д ) . (20)Differentiating in t at t = □ The stabilizer stab ( v ) of a vector v is the subgroup consisting ofall elements д ∈ G that preserve v —that is, such that ρ ( д ) v = v .The Lie algebra so ( ) associated to SO ( ) consists of the skew-symmetric 3 × so ( ) has a basis consisting of infinitesimalrotations about the three coordinate axes: l = (cid:169)(cid:173)(cid:171) − (cid:170)(cid:174)(cid:172) l = (cid:169)(cid:173)(cid:171) − (cid:170)(cid:174)(cid:172) l = (cid:169)(cid:173)(cid:171) − (cid:170)(cid:174)(cid:172) . Their Lie brackets are [ l i , l i ] = [ l , l ] = l [ l , l ] = l [ l , l ] = l . (21)These generators might be familiar as the angular momentum op-erators from quantum mechanics. Indeed, SO ( ) acts on functionson the sphere by rotating them, and its Lie algebra elements act asrotational derivatives. This action decomposes into irreducible rep-resentations. A basis for each irreducible representation of SO ( ) isprovided by spherical harmonics ; each irreducible representationis referred to as a band of spherical harmonics. The real represen-tations are odd-dimensional vector spaces, and the representationmatrices ρ (·) are also called Wigner D -Matrices.SO ( ) admits a bi-invariant Riemannian metric whose valueat the identity is given by ⟨ u , v ⟩ = (cid:213) i , j u ij v ij . (22)Note that the generators l i form an orthonormal basis under thismetric. Bi-invariance means that ⟨( D L д ) I v , ( D L д ) I w ⟩ = ⟨ v , w ⟩ = ⟨( D R д ) I v , ( D R д ) I w ⟩ , (23)where v , w ∈ so ( ) and R д is right translation by д ∈ SO ( ) , definedanalogously to left translation. (23) completely describes the Rie-mannian metric on SO ( ) . Under this metric, the matrix exponentialis also the Riemannian exponential map at the identity, giving riseto geodesics.The octahedral group O is a discrete subgroup of SO ( ) gener-ated by right-angle rotations about the three axes, (cid:110) r i = exp (cid:16) π l i (cid:17)(cid:111) i = . O has 24 elements comprising all the symmetries of the cube, orequivalently, of its dual, the octahedron.
A.2 Semidefinite Relaxations A (real) algebraic variety is the set of solutions of a system ofpolynomial equations over n real variables: V ( p , . . . , p k ) = { x ∈ R n | p i ( x ) = , i = , . . . , k } , where p , . . . , p k are polynomials.If p , . . . , p k are quadratic in x , they may be written p i ( x ) = (cid:16) x ⊤ (cid:17) A i (cid:18) x (cid:19) = ⟨ A i , X ⟩ , where A i are symmetric matrices, ⟨ , ⟩ denotes the usual inner prod-uct on matrices inducing the Frobenius norm, and X = (cid:18) x (cid:19) (cid:16) x ⊤ (cid:17) = (cid:18) x ⊤ x xx ⊤ (cid:19) . ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020.
A quadratically-constrained quadratic program (
QCQP ) is of theform arg min { f ( x ) | x ∈ V ( p , . . . , p k )} , where f , p , . . . , p k are all quadratic. A QCQP may be rewritten inthe more illuminating formarg min X ⟨ C , X ⟩ subject to X = ⟨ A i , X ⟩ = , i = , . . . , kX ⪰ ( X ) = . (24)Absent the rank constraint, (24) would be a semidefinite program( SDP ). Semidefinite programs may be solved in polynomial timeby interior point methods, implemented in common optimizationsoftware packages like Mosek [2020]. Thus it is natural to ignorethe rank constraint and solve the associated SDP; this technique iscalled semidefinite relaxation . The optimal objective value of therelaxed problem is a lower bound for the globally optimal value ofthe QCQP (24), and if the recovered solution to the SDP has rankone, the solution is said to be exact , as it is the globally optimalsolution of (24).
ACKNOWLEDGMENTS
The authors would like to thank Elina Robeva for helping us un-derstand her work on the odeco varieties, Pablo Parrilo and DiegoCifuentes for sharing their insights into semidefinite relaxation,and Nicolas Boumal for helping us get started with manifold opti-mization. We thank Heng Liu for providing tetrahedral refinement,parametrization, and hexahedral mesh extraction code. We are grate-ful to Xifeng Gao and Daniele Panozzo for guiding us in using theircode from [Gao et al. 2017]. Finally, we thank Mirela Ben-Chen,Amit Bermano, Edward Chien, Etienne Corman, Keenan Crane, Fer-nando de Goes, Steven Gortler, Leif Kobbelt, Braxton Osting, NirSharon, Nir Sochen, Amir Vaxman, Josh Vekhter, and Paul Zhangfor valuable discussions.David Palmer appreciates the generous support of the Fannieand John Hertz Foundation through the Hertz Graduate Fellowship.David Bommes appreciates the generous support of the EuropeanResearch Council (ERC) under the European UnionâĂŹs Horizon2020 research and innovation program (AlgoHex, grant agreementno. 853343). Justin Solomon acknowledges the generous support ofArmy Research Office grant W911NF-12-R-0011, National ScienceFoundation grant IIS-1838071, Air Force Office of Scientific Researchaward FA9550-19-1-0319, and a gift from Adobe Systems. Any opin-ions, findings, and conclusions or recommendations expressed inthis material are those of the authors and do not necessarily reflectthe views of these organizations.
REFERENCES
P.-A. Absil, C.G. Baker, and K.A. Gallivan. 2007. Trust-Region Methods on RiemannianManifolds.
Foundations of Computational Mathematics
7, 3 (01 Jul 2007), 303–330.https://doi.org/10.1007/s10208-005-0179-9Motilal Agrawal and Larry S Davis. 2003. Camera calibration using spheres: A semi-definite programming approach. In
Proc. ICCV . IEEE, 782.MOSEK ApS. 2020.
The MOSEK Fusion API for C++ manual. Version 9.1. https://docs.mosek.com/9.1/cxxfusion/index.html Cecil G. Armstrong, Harold J. Fogg, Christopher M. Tierney, and Trevor T. Robinson.2015. Common Themes in Multi-block Structured Quad/Hex Mesh Generation.
Pro-cedia Engineering arXiv preprint arXiv:1910.06240 (2019).Pierre-Alexandre Beaufort, Jonathan Lambrechts, François Henrotte, ChristopheGeuzaine, and Jean-François Remacle. 2017. Computing cross fields: A PDE ap-proach based on the Ginzburg–Landau theory.
Procedia Engineering
203 (2017),219–231.Grigoriy Blekherman, Pablo A Parrilo, and Rekha R Thomas. 2012.
Semidefinite opti-mization and convex algebraic geometry . SIAM.Ada Boralevi, Jan Draisma, Emil Horobeţ, and Elina Robeva. 2017. Orthogonal andunitary tensor decomposition from an algebraic perspective.
Israel Journal ofMathematics
Journal of Machine Learning Research
International Meshing Roundtable (Aug. 2018).Diego Cifuentes, Sameer Agarwal, Pablo A Parrilo, and Rekha R Thomas. 2017. On thelocal stability of semidefinite relaxations. arXiv preprint arXiv:1710.04287 (2017).Diego Cifuentes, Corey Harris, and Bernd Sturmfels. 2018. The Geometry of SDP-Exactness in Quadratic Optimization. arXiv preprint arXiv:1804.01796 (2018).Etienne Corman and Keenan Crane. 2019. Symmetric Moving Frames.
ACM Trans.Graph.
38, 4 (2019).Xianzhong Fang, Weiwei Xu, Hujun Bao, and Jin Huang. 2016. All-hex Meshing UsingClosed-form Induced Polycube.
ACM Trans. Graph.
35, 4, Article 124 (July 2016),9 pages. https://doi.org/10.1145/2897824.2925957Xifeng Gao, Wenzel Jakob, Marco Tarini, and Daniele Panozzo. 2017. Robust Hex-dominant Mesh Generation Using Field-guided Polyhedral Agglomeration.
ACMTrans. Graph.
36, 4, Article 114 (July 2017), 13 pages.Michel X Goemans and David P Williamson. 1995. Improved approximation algorithmsfor maximum cut and satisfiability problems using semidefinite programming.
Jour-nal of the ACM (JACM)
42, 6 (1995), 1115–1145.Dmitry Golovaty, Jose Alberto Montero, and Daniel Spirn. 2019. A variationalmethod for generating n -cross fields using higher-order Q -tensors. arXiv preprintarXiv:1909.00922 (2019).Jin Huang, Yiying Tong, Hongyu Wei, and Hujun Bao. 2011. Boundary Aligned Smooth3D Cross-frame Field. ACM Trans. Graph.
30, 6, Article 143 (Dec. 2011), 8 pages.https://doi.org/10.1145/2070781.2024177Qi-Xing Huang and Leonidas Guibas. 2013. Consistent shape maps via semidefiniteprogramming. In
Proc. Symposium on Geometry Processing . Eurographics Association,177–186.Tengfei Jiang, Jin Huang, Yuanzhen Wang, Yiying Tong, and Hujun Bao. 2014. FrameField Singularity Correction for Automatic Hexahedralization.
IEEE Trans. on Visu-alization & Computer Graphics
20, 8 (Aug. 2014), 1189–1199.Jan Kautz, John Snyder, and Peter-Pike J Sloan. 2002. Fast Arbitrary BRDF Shadingfor Low-Frequency Lighting Using Spherical Harmonics.
Rendering Techniques
Computer Graphics Forum , Vol. 34. Wiley Online Library,115–128.Felix Knöppel, Keenan Crane, Ulrich Pinkall, and Peter Schröder. 2013. Globally optimaldirection fields.
ACM Trans. Graph.
32, 4 (2013), 59.Shahar Z Kovalsky, Noam Aigerman, Ronen Basri, and Yaron Lipman. 2014. Controllingsingular values with semidefinite programming.
ACM Trans. Graph.
33, 4 (2014),68–1.J. Lee. 2012.
Introduction to Smooth Manifolds (second ed.). Springer New York.Yufei Li, Yang Liu, Weiwei Xu, Wenping Wang, and Baining Guo. 2012. All-hex meshingusing singularity-restricted field.
ACM Trans. Graph.
31, 6 (2012), 177.Heng Liu, Paul Zhang, Edward Chien, Justin Solomon, and David Bommes. 2018.Singularity-constrained Octahedral Fields for Hexahedral Meshing.
ACM Trans.Graph.
37, 4, Article 93 (July 2018), 17 pages. https://doi.org/10.1145/3197517.3201344Max Lyon, David Bommes, and Leif Kobbelt. 2016. HexEx: Robust Hexahedral MeshExtraction.
ACM Trans. Graph.
35, 4, Article 123 (July 2016), 11 pages.Haggai Maron, Nadav Dym, Itay Kezurer, Shahar Kovalsky, and Yaron Lipman. 2016.Point registration via efficient convex relaxation.
ACM Trans. Graph.
35, 4 (2016),73.N. D. Mermin. 1979. The topological theory of defects in ordered media.
Reviews ofModern Physics
51, 3 (July 1979), 591–648. https://doi.org/10.1103/RevModPhys.51.591ACM Trans. Graph., Vol. 39, No. 2, Article 16. Publication date: March 2020. lgebraic Representations for Volumetric Frame Fields • 16:17
Barry Merriman, James Kenyard Bence, and Stanley Osher. 1992.
Diffusion generatedmotion by mean curvature . Department of Mathematics, University of California,Los Angeles.John Douglas Moore. 1976. Equivariant Embeddings of Riemannian HomogeneousSpaces.
Indiana University Mathematics Journal
Computer Graphics Forum , Vol. 30. Wiley OnlineLibrary, 1397–1406.Braxton Osting and Dong Wang. 2017. A generalized MBO diffusion generated motionfor orthogonal matrix-valued fields. arXiv preprint arXiv:1711.01365 (2017).Onur Ozyesil, Amit Singer, and Ronen Basri. 2015. Stable camera motion estimationusing convex programming.
SIAM Journal on Imaging Sciences
8, 2 (2015), 1220–1262.Jonathan Palacios, Lawrence Roy, Prashant Kumar, Chen-Yuan Hsu, Weikai Chen,Chongyang Ma, Li-Yi Wei, and Eugene Zhang. 2017. Tensor Field Design in Volumes.
ACM Trans. Graph.
36, 6, Article 188 (Nov. 2017), 15 pages. https://doi.org/10.1145/3130800.3130844Nicolas Ray, Dmitry Sokolov, and Bruno Lévy. 2016. Practical 3D Frame Field Generation.
ACM Trans. Graph.
35, 6, Article 233 (Nov. 2016), 9 pages. https://doi.org/10.1145/2980179.2982408E. Robeva. 2016. Orthogonal Decomposition of Symmetric Tensors.
SIAM J.Matrix Anal. Appl.
37, 1 (2016), 86–102. https://doi.org/10.1137/140989340arXiv:https://doi.org/10.1137/140989340Zhongwei Shen, Xianzhong Fang, Xinguo Liu, Hujun Bao, and Jin Huang. 2016. Har-monic functions for rotational symmetry vector fields. In
Computer Graphics Forum , Vol. 35. Wiley Online Library, 507–516.Amit Singer. 2011. Angular synchronization by eigenvectors and semidefinite program-ming.
Applied and Computational Harmonic Analysis
30, 1 (2011), 20.Dmitry Sokolov, Nicolas Ray, Lionel Untereiner, and Bruno Lévy. 2016. Hexahedral-Dominant Meshing.
ACM Trans. Graph.
35, 5, Article 157 (June 2016), 23 pages.Justin Solomon, Amir Vaxman, and David Bommes. 2017. Boundary Element OctahedralFields in Volumes.
ACM Trans. Graph.
36, 3, Article 114b (May 2017). https://doi.org/10.1145/3065254John Stillwell. 2008.
Naive Lie Theory . Springer Science & Business Media.Yves van Gennip, Nestor Guillen, Braxton Osting, and Andrea L. Bertozzi. 2014. MeanCurvature, Threshold Dynamics, and Phase Field Theory on Finite Graphs.
MilanJournal of Mathematics
82, 1 (01 Jun 2014), 3–65. https://doi.org/10.1007/s00032-014-0216-8Amir Vaxman, Marcel Campen, Olga Diamanti, Daniele Panozzo, David Bommes, KlausHildebrandt, and Mirela Ben-Chen. 2016. Directional field synthesis, design, andprocessing. In
Computer Graphics Forum , Vol. 35. Wiley Online Library, 545–572.R. Viertel and B. Osting. 2019. An Approach to Quad Meshing Based on HarmonicCross-Valued Maps and the Ginzburg–Landau Theory.
SIAM Journal on Scien-tific Computing
41, 1 (2019), A452–A479. https://doi.org/10.1137/17M1142703arXiv:https://doi.org/10.1137/17M1142703W. Yu, K. Zhang, and X. Li. 2015. Recent algorithms on automatic hexahedral meshgeneration. In . 697–702.