Structured Regularization of Functional Map Computations
EEurographics Symposium on Geometry Processing 2019D. Bommes and H. Huang(Guest Editors)
Volume 38 ( ), Number 5
Structured Regularization of Functional Map Computations
Jing Ren , Mikhail Panine , Peter Wonka , and Maks Ovsjanikov KAUST, LIX, École Polytechnique, CNRS
Abstract
We consider the problem of non-rigid shape matching using the functional map framework. Specifically, we analyze a commonlyused approach for regularizing functional maps, which consists in penalizing the failure of the unknown map to commute withthe Laplace-Beltrami operators on the source and target shapes. We show that this approach has certain undesirable funda-mental theoretical limitations, and can be undefined even for trivial maps in the smooth setting. Instead we propose a novel,theoretically well-justified approach for regularizing functional maps, by using the notion of the resolvent of the Laplacian oper-ator. In addition, we provide a natural one-parameter family of regularizers, that can be easily tuned depending on the expectedapproximate isometry of the input shape pair. We show on a wide range of shape correspondence scenarios that our novelregularization leads to an improvement in the quality of the estimated functional, and ultimately pointwise correspondencesbefore and after commonly-used refinement techniques.
CCS Concepts • Computing methodologies → Shape analysis;
1. Introduction
Shape matching is a fundamental problem in geometry processingand computer graphics, with applications ranging from shape inter-polation [HRS ∗
16] to shape exploration [HWG14] and statisticalshape analysis [BRLB14].In this paper, we concentrate on the functional maps frameworkfor computing correspondences between shapes, which has provento be especially useful when dealing with non-rigid and especiallynear-isometric shape pairs. This approach, originally introducedin [OBCS ∗
12] and subsequently extended in multiple follow-upworks, e.g. [KBB ∗
13, RCB ∗
17, HWG14], is based on optimizingfor a linear mapping between function spaces defined on the shapes,which can be conveniently encoded as a matrix. One attractive fea-ture of this approach is that many desirable geometric objectivesfor the unknown pointwise correspondence, such as preservation ofgeodesic distances, can be conveniently encoded as constraints onthe matrices representing functional maps, and often lead to simpleconvex optimization problems.The majority of objectives when optimizing for a functional mapbetween a pair of shapes consist of a) preservation of some pre-computed descriptors and b) a functional map regularization term,based on promoting some desired global map properties. This lat-ter, regularization term is especially important for ensuring theoverall global consistency of the computed map.Common strategies for regularization include commutativitywith the Laplacian operator, exploited in the original article[OBCS ∗ Ours R e gu l a r i z i ng m a s k P o i n t w i s e m a p Figure 1:
Comparison of different masks. Standard: the Lapla-cian commutativity operator can be equivalently formulated as apenalty or regularizing mask. Slanted: the weight mask proposedin [RCB ∗
17] designed to promote a slanted structure. Ours: themask proposed in this paper based on the resolvent of the Laplacianoperator. The penalty masks of these three methods are visualizedin the first row for an example pair from the SHREC dataset, andthe corresponding optimized point-wise maps from the functionalmap pipeline are shown in the second row. partial correspondences [RCB ∗ area preservation by theunderlying pointwise map [OBCS ∗
12, KBB ∗
13, ROA ∗ c (cid:13) (cid:13)(cid:13)
13, ROA ∗ c (cid:13) (cid:13)(cid:13) a r X i v : . [ c s . G R ] S e p . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations objectives which would not be well-defined in the case of smoothsurfaces.In this paper, we introduce a novel way of regularizing functionalmaps, based on the concept of the resolvent of the Laplacian oper-ator. Unlike the original approach, based on commutativity withthe Laplacians, our method has a natural theoretically well-definedanalogue in the case of smooth surfaces, and is guaranteed to bebounded even in the full (infinite dimensional) basis case. More-over, our approach provides a simple one-parameter family of reg-ularizers, that can be tuned depending on the expected approximatedeviation from isometry in a given shape pair. Finally, we demon-strate on a wide range of challenging settings that our approachleads to a quantitative and qualitative improvement for the com-puted functional and eventually pointwise maps, without incurringany additional computational complexity (see Fig. 1 as an exam-ple).
2. Related Work
Shape matching in its full generality is an extremely well-studiedarea of geometry processing and computer graphics, and its fulloverview is beyond the scope of our paper. Therefore, below weonly review the most closely related methods based primarily onthe functional maps framework. We refer the interested readers torecent surveys including [VKZHCO11, TCL ∗
13, BCBB16] for anin-depth treatment of other shape matching approaches.
Functional Maps.
Our approach fits within the functional mapframework, which was originally introduced in [OBCS ∗
12] forsolving non-rigid shape matching problems, and extended signifi-cantly in follow-up works, including [KBB ∗
13, ADK16, KBBV15,RCB ∗ ∗ ∗ ∗
13, RCB ∗
17, BDK17] many natural properties on the un-derlying pointwise correspondences can be expressed as objec-tives on functional maps. Most notably, this includes: orthonor-mality of functional maps, which corresponds to the local area-preservation nature of pointwise correspondences [OBCS ∗ ∗
13, ROA ∗ ∗ pointwise products of func-tions, which corresponds to functional maps arising from point-to-point correspondences [NO17, NMR ∗ ∗
17, LRBB17]. Similarly, several otherregularizers have been proposed, including using robust norms andmatrix completion techniques [KBB ∗
13, KBBV15], exploiting therelation between functional maps in different directions [ERGB16],the map adjoint [HO17], and powerful cycle-consistency con-straints [HWG14] in the context of shape collections, among many others. More recently constraints on functional maps have beenintroduced to promote continuity of the recovered pointwise cor-respondence [PSO18], maps between curves defined on shapes[GBKS18], kernel-based techniques aimed at extracting more in-formation from given descriptor constraints [WGBS18], and an ap-proach for incorporating orientation information into the functionalmap infererence pipeline [RPWO18] among others.Among all of these, perhaps the most widely-used building blockfor regularizing functional maps, introduced in [OBCS ∗
12] and ex-tended in follow-up works, including [WHG13,RCB ∗ ∗ ∗ ∗ Optimal Transport.
We also note briefly that other commonly-used relaxations for matching problems include those based on op-timal transport, e.g. [SPKS16, MCSK ∗
17, VLR ∗ ∗ ∗
17, VLR ∗ kernels ratherthan preservation of e.g. geodesic distances. Our modification ofthe regularization term of functional maps can also be seen througha similar light, as we show that a more theoretically justified func-tional operator leads to an improvement of the overall structuralproperties of the map, which eventually result in more accuratefunctional and pointwise maps.
3. Background
Our work is based on the functional map representation and theestimation pipeline. Below we review the basic notions and themain steps for estimating a map between a pair of shapes using c (cid:13) (cid:13) . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations this framework, and refer the interested reader to a recent set ofcourse notes [OCB ∗
17] for a more in-depth discussion.
Basic Pipeline.
Given a pair of shapes, S , S represented as tri-angle meshes, and containing, respectively, n and n vertices, thegeneral pipeline for computing a map between them using the func-tional map representation, consists of the following main steps:1. Compute a small set of k << n and k << n basis func-tions on each shape. The most common choice consists in us-ing the first k eigenfunctions of the Laplace-Beltrami operatorof each shape, although other bases derived from the Hamil-tonian operator [CSBK18] and more localized basis functions[NVT ∗
14, MRCB18] have also been used.2. Compute a set of descriptor functions on each shape, that areexpected to be approximately preserved by the unknown map.Store their coefficients in the corresponding bases as columns ofmatrices A , A .3. Compute the optimal functional map C by solving the followingoptimization problem: C opt = arg min C E desc (cid:0) C (cid:1) + α E reg (cid:0) C (cid:1) (1)where the first term aims at the descriptor preservation: E desc (cid:0) C (cid:1) = (cid:13)(cid:13) C A − A (cid:13)(cid:13) , whereas the second term reg-ularizes the map by promoting the correctness of its overallstructural properties. As mentioned above, the most commonapproach consists of penalizing the failure of the unknown func-tional map to commute with the Laplace-Beltrami operators,which can be written as: E reg ( C ) = E comm (cid:0) C (cid:1) = (cid:13)(cid:13) C ∆ − ∆ C (cid:13)(cid:13) (2)where ∆ and ∆ are the Laplace-Beltrami operators of the twoshapes expressed in the respective bases. Here, and throughoutthe rest of our paper, unless stated otherwise (cid:107) · (cid:107) denotes thematrix Frobenius norm.4. Convert the functional map C to a point-to-point map, for ex-ample using an iterative approach, such as the Iterative ClosestPoint (ICP) in the spectral embedding, or using other more ad-vanced techniques [RMC15, EBC17].One of the attractive properties of this pipeline is that the func-tional map computation in step 3 leads to a simple (convex) leastsquares optimization with a relatively small number of unknowns,independent of the number of points on the shapes. This step hasbeen further extended e.g. using more powerful descriptor preser-vation constraints via commutativity [NO17], or using manifoldoptimization [KGB16] among many others (see also Chapter 3in [OCB ∗
4. Functional Map Regularization
Our main goals are to analyze the commonly-used functional mapregularization term, to bring attention to some of its theoreticallimitations, and to propose a novel regularizer with better theoret-ical properties, which lead to practical improvements. Therefore,in this in this section, we first consider the standard Laplacian-commutativity term, and then introduce our new regularizer basedon the resolvent operator of the Laplacian. Given descriptorGround-truth map G r ound - t r u t h O u r s S t a nd a r d Mask Functionalmap Pointwisemap
Figure 2:
Given a single pair of WKS descriptors, we optimize a100-by-100 functional map using the standard mask and our resol-vent mask and compare to ground truth. The converted point-wisemaps are shown on the right. We can see that the functional mapstemming from our resolvent mask has less noise than the func-tional map computed with the standard mask. Also, our resolventmask leads to a point-wise map with better quality.
As mentioned above, the Laplacian-commutativity term given inEq. (2) was first introduced to promote approximately isometriccorrespondences. If the functional map is expressed in the basisgiven by the eigenfunctions of the Laplace-Beltrami operator, andletting Λ and Λ represent vectors that store the eigenvalues ofthe Laplacians of shape S and S respectively, this term can beequivalently reformulated as: E comm (cid:0) C (cid:1) = (cid:13)(cid:13) C ∆ − ∆ C (cid:13)(cid:13) = (cid:13)(cid:13) C diag (cid:0) Λ (cid:1) − diag (cid:0) Λ (cid:1) C (cid:13)(cid:13) = (cid:13)(cid:13) C (cid:12) (cid:0) k Λ T (cid:1) − (cid:0) Λ Tk (cid:1) (cid:12) C (cid:13)(cid:13) = (cid:13)(cid:13)(cid:0) k Λ T − Λ Tk (cid:1) (cid:12) C (cid:13)(cid:13) (cid:44) E mask (cid:0) C (cid:1) = ∑ i j [ M LB ] i j [ C ] i j , (3)where (cid:12) is the entry-wise matrix multiplication, k is a k -dim all-ones vector and [ A ] i j denotes the ( i , j ) entry of a matrix A .In other words, the regularization term E comm can be written asa product between the matrix M LB , which we call the penalty mask matrix and the squares of the entries of the unknown functionalmap C . In the case of E comm , the entries of M LB are given as M LB ( i , j ) = ( Λ ( i ) − Λ ( j )) . Fig. 2 shows a heat map of this ma-trix (see the first row in the “Mask” column).Unfortunately, this simple formulation has certain fundamentalundesirable properties. In particular, in the case of smooth surfaces,Laplace-Beltrami operators are unbounded operators on square-integrable functions [MP49]. Consequently, in general, the energyterm (cid:107) C ∆ − ∆ C (cid:107) is undefined on smooth surfaces. As a sim-ple example, consider the situation where C = Id , the identityoperator, and the two surfaces are scaled versions of each other c (cid:13) (cid:13) . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations , k λ k TorusSphereWeyl Estimate
Figure 3:
The spectra of a sphere and a square torus of unit area,as well as the corresponding Weyl estimate. (i.e., ∆ = c ∆ for some constant c (cid:54) = (cid:107) C ∆ − ∆ C (cid:107) = | − c | (cid:107) ∆ (cid:107) is infinite.That being said, the ill-definiteness of (cid:107) C ∆ − ∆ C (cid:107) is nota mere question of scale. Recall that, by Weyl’s law [Dod81], largeLaplacian eigenvalues of surfaces can be estimated in terms of thesurface area S as follows: λ k ∼ π S k . (4)Thus, by rescaling the surfaces such that their areas match guar-antees that their eigenvalues have comparable asymptotic growth.This, however, is not sufficient to make (cid:107) C ∆ − ∆ C (cid:107) well-defined. This can be seen from another simple, if slightly arti-ficial, example. Consider a sphere ( M ) and a flat square torus( M ), both of unit area. While these surfaces are not diffeomor-phic, they are among the few whose Laplace-Beltrami eigenvaluescan be explicitly computed (see [Sau06a,Sau06b] among many oth-ers), which is why we consider them here. Said spectra, as well asthe corresponding Weyl estimate are illustrated in Fig. 3. For sim-plicity’s sake, we once again pick a functional map of the form C = Id . We illustrate the spectra, as well as the corresponding (cid:107) C ∆ − ∆ C (cid:107) on Fig. 4. Notice the divergence of the Lapla-cian commutator energy. On the same figure, we illustrate its pro-posed replacement, defined in the next section. Note its rapid con-vergence.In addition to being ill-defined in the continuous setting, theLaplacian commutativity mask has another significant problem.Namely, it penalizes the high frequency entries of functional mapsin the same way as the low frequency ones, despite the greater insta-bility of the former, and in this way fails to exhibit the funnel-likestructure observed in ground-truth functional maps (see e.g. Fig. 2above or Figure 4 in [OBCS ∗ . . . . k (cid:107) ∆ − ∆ (cid:107) . . . (cid:107) R ( ∆ ) − R ( ∆ ) (cid:107) Figure 4:
The Laplacian and Resolvent commutator energies forC = Id computed using the lowest k eigenvalues of the unitarea sphere and torus (see Fig. 3). Note the rapid convergence of (cid:107) C R ( ∆ ) − R ( ∆ ) C (cid:107) and the divergence of (cid:107) C ∆ − ∆ C (cid:107) .The spectra are rescaled with respect to the largest eigenvalue withk = , as described in Eq. (14) . to the Laplacian commutativity mask (labeled “Standard”), aver-aged over the same shape pairs. We perform the same computationfor the heuristic slanted mask introduced in [RCB ∗
17] and, finally,for our proposed resolvent-based mask defined in the next section.Notice that our mask better reproduces the funnel-like structure ofthe ground-truth maps.
In this section, we propose an alternative form for the Laplaciancommutativity regularizer that overcomes the problems identifiedin the previous section. The resulting regularizer is better from boththe theoretical and practical standpoints.The first issue of the original Laplacian-commutativity termarises from the fact that the Laplace-Beltrami operator is un-bounded. We thus propose to replace it with a meaningful boundedoperator, namely its resolvent. Below, we give a brief overview andrefer the reader to [RS80] for the detailed functional-analytic un-derpinnings of the following discussion.
Resolvent.
Let A : D → H be a densely-defined closed operatoron some Hilbert space H with domain D ⊂ H . Let ρ ( A ) ⊂ C bethe set of all complex numbers µ such that R µ ( A ) = ( A − µId ) − is defined and bounded. The set ρ ( A ) and the operator R µ ( A ) areknown as the resolvent set of A and the resolvent (operator) of A at µ , respectively. The resolvent set ρ ( A ) is the complement of thespectrum of A in the complex plane.The general idea of considering the resolvent of an unboundedoperator rather than the operator itself comes from the notion ofnorm-resolvent convergence, which is used to study the conver-gence of unbounded self-adjoint operators [RS80]. In that sense,our choice to use the resolvent of the Laplace-Beltrami operator isa natural one. Note that here, closedness is a technical conditionused in the definition of the resolvent. In particular, it is satisfied byself-adjoint operators, such as the ones we consider.Before proceeding further, we slightly generalize our problem.In what follows we will consider powers of the Laplacian ratherthan the Laplacian itself. Specifically, we will use ∆ γ for some γ > c (cid:13) (cid:13) . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations Ground-truth Standard mask Slanted mask Our mask
Figure 5:
The average squared ground-truth functional map andthree penalty masks over 100 FAUST non-isometric shape pairs. ∆ . As explained below, the parameter γ controls thefunnel-like structure of the resulting mask. Thus, this parametertakes on an important role in the numerical tests reported later inthis paper.Now since the Laplace-Beltrami operator ∆ is positive and self-adjoint, its spectrum is contained in the non-negative real line. Thesame also holds for ∆ γ . Thus, we are guaranteed that R µ ( ∆ ) is awell-defined bounded operator for any complex µ not in the non-negative real line.For our purposes, the resolvent of ∆ γ is expressed as R ( ∆ γ ) = (cid:0) ∆ γ − ( a + ib ) Id (cid:1) − , (5)where i is the imaginary unit and a , b ∈ R . Our proposal is to usethe resolvent of the Laplacian instead of the Laplacian in the com-mutator term. We thus define a new energy term, E resolvent (cid:0) C (cid:1) = (cid:13)(cid:13) C R (cid:0) ∆ γ (cid:1) − R (cid:0) ∆ γ (cid:1) C (cid:13)(cid:13) . (6)Before proceeding further, recall that we only need to considerbounded functional maps C , since functional maps that arise aspullbacks of diffeomorphisms are bounded. This last fact is shownin appendix A, for completeness. Our proposal is motivated by thefollowing result. Theorem 1 (Bounded Resolvent Commutativity)
Let C be abounded functional map. Then, in the operator norm, (cid:13)(cid:13) C R (cid:0) ∆ γ (cid:1) − R (cid:0) ∆ γ (cid:1) C (cid:13)(cid:13) op < ∞ . (7) Proof
The result directly follows from the fact that products andlinear combinations of bounded operators are bounded. As definedabove, R ( ∆ γ ) is bounded, which concludes the proof.As explained above, there is no analogous guarantee for the com-monly used Laplacian commutativity, since the Laplace-Beltramioperator itself is unbounded. In contrast, our resolvent-based com-mutativity is always well-defined, and leads to a commutator withbounded operator norm, even in the case of smooth surfaces.Notice that Theorem 1 holds in the operator norm, rather thanthe more convenient Frobenius norm, which we use to define theenergy in Equation (6). Our usage of the Frobenius norm can bejustified in two ways. First, a Frobenius norm version of Theorem1 holds for γ > /
2. This is shown in Lemma 2 of Appendix B.Thus, for γ > /
2, the usage of the Frobenius norm is perfectlyjustified in the smooth setting. For lower values of γ , we invoke thefact that, in practice, we work on a finite dimensional vector space.Since all norms are equivalent on finite dimensional vector spaces, we are free to replace the operator norm with the Frobenius one.In that case, we lose the guarantee that the energy E resolvent makessense in the full Laplace-Beltrami basis in the smooth setting. Yet,as shown in the experiments reported below, this does not seem tocause any issue.Similarly to the usual Laplacian commutativity term E comm , E resolvent can also be expressed as a mask matrix. This can be doneas follows. Notice that R ( ∆ γ ) is diagonal in the Laplacian eigenba-sis. Moreover, if ∆ has eigenvalues { λ n } ∞ n = , then, by matrix inver-sion, the eigenvalues { r n } ∞ n = of R ( ∆ γ ) will be given by: r n = λ γ n − a − ib = λ γ n − a ( λ γ n − a ) + b + b ( λ γ n − a ) + b i . (8)In practice, we choose a = b =
1, which yields r n = λ γ n ( λ γ n ) + + ( λ γ n ) + i . (9)Since the square of the Frobenius norm independently considers thereal and imaginary parts of a matrix, we can re-express Equation (6)in terms of two new matrices M Re and M Im , defined below. E resolvent (cid:0) C (cid:1) = (cid:13)(cid:13) M Re (cid:12) C (cid:13)(cid:13) + (cid:13)(cid:13) M Im (cid:12) C (cid:13)(cid:13) = ∑ i j [ M Re ] i j [ C ] i j + ∑ i j [ M Im ] i j [ C ] i j = ∑ i j (cid:16) [ M Re ] i j + [ M Im ] i j (cid:17) [ C ] i j (10)The matrices M Re and M Im correspond to the real and imaginaryparts of the eigenvalues of the resolvent, respectively. Explicitly,these matrices are given by M Re ( i , j ) = Λ ( i ) γ Λ ( i ) γ + − Λ ( j ) γ Λ ( j ) γ + M Im ( i , j ) = Λ ( i ) γ + − Λ ( j ) γ + M res ( i , j ) = M Re ( i , j ) + M Im ( i , j ) . (13)The split of M res into M Re and M Im will be revisited later, whenwe explore beyond the established theory and consider a maskconstructed from weighted combinations of these two matrices. As it stands now, we observe that in practice the mask definedin equation (13) decays too quickly as the Laplacian eigenvaluesgrow. In other words, the mask is not sufficiently sensitive to the c (cid:13) (cid:13)(cid:13)
1, which yields r n = λ γ n ( λ γ n ) + + ( λ γ n ) + i . (9)Since the square of the Frobenius norm independently considers thereal and imaginary parts of a matrix, we can re-express Equation (6)in terms of two new matrices M Re and M Im , defined below. E resolvent (cid:0) C (cid:1) = (cid:13)(cid:13) M Re (cid:12) C (cid:13)(cid:13) + (cid:13)(cid:13) M Im (cid:12) C (cid:13)(cid:13) = ∑ i j [ M Re ] i j [ C ] i j + ∑ i j [ M Im ] i j [ C ] i j = ∑ i j (cid:16) [ M Re ] i j + [ M Im ] i j (cid:17) [ C ] i j (10)The matrices M Re and M Im correspond to the real and imaginaryparts of the eigenvalues of the resolvent, respectively. Explicitly,these matrices are given by M Re ( i , j ) = Λ ( i ) γ Λ ( i ) γ + − Λ ( j ) γ Λ ( j ) γ + M Im ( i , j ) = Λ ( i ) γ + − Λ ( j ) γ + M res ( i , j ) = M Re ( i , j ) + M Im ( i , j ) . (13)The split of M res into M Re and M Im will be revisited later, whenwe explore beyond the established theory and consider a maskconstructed from weighted combinations of these two matrices. As it stands now, we observe that in practice the mask definedin equation (13) decays too quickly as the Laplacian eigenvaluesgrow. In other words, the mask is not sufficiently sensitive to the c (cid:13) (cid:13)(cid:13) . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations Figure 6:
The resolvent mask introduced in Eq. (13) with different γ . The darker the region, the larger the penalty. Here, ∆ and ∆ are near-isospectral, which explains the approximately symmetricstructure of the mask. higher frequencies. This is due to the scale introduced by the pa-rameter b found in the definition of the resolvent. There are a fewequivalent ways of addressing this issue.Our approach is as follows. We begin by computing the k lowesteigenvalues of the considered Laplacians. Then, both spectra arerescaled according to the rule Λ i (cid:55)−→ Λ i max (cid:0) max ( Λ ) , max ( Λ ) (cid:1) . (14)In other words, the spectra of the source and target shape arerescaled by the same factor, such that the largest considered eigen-value (over both spectra) becomes equal to 1. From Equation (8),one can see that this rescaling can be absorbed into a choice of b and a change in the weight of the resolvent energy in the over-all energy. Thus, conceptually, the spectral rescaling is equivalentto a choice of the parameter b used in the definition of the resol-vent. Alternatively, recall that rescaling the spectrum is equivalentto rescaling the surface. Consequently, this rescaling does not affectthe theoretical guarantees of the previous section.Returning to the previously discussed example of the sphere andtorus of unit area (see Fig. 3), we illustrate the resolvent energycomputed by the above procedure in Fig. 4. Note the rapid conver-gence of the resolvent commutator energy (in red) and the diver-gence of the Laplacian commutator energy (in blue). γ Having defined the way to compute the resolvent mask, we are nowready to explore the way in which tuning the parameter γ controlsits the funnel-like structure.Fig. 6 illustrates the resolvent mask for different values of γ . Thedarker the region, the more penalized the corresponding entry in thefunctional map C would be. Notice that the funnel-like structureof the mask changes with different values of γ .From Fig. 6, we see that the behaviour of the mask as a functionof γ can be separated into two regimes. The first corresponds to γ ∈ ( , ] . There, increasing γ results in a narrowing of the funnel-like shape of the mask. Thus, the larger γ is, the more the mask pe- nalizes functional maps that take eigenfunctions of ∆ to eigenfunc-tions of ∆ with distant eigenvalues. Correspondingly, small valuesof γ are more lax in that regard, allowing for maps between eigen-functions with quite different eigenvalues. The former choice of γ isappropriate when the shapes under consideration are approximatelyisometric, as one then expects the eigenfunctions and eigenvaluesof both surfaces to be roughly the same. The latter choice is soundfor shape pairs that are further away from isometry.The second regime corresponds to γ >
1. There, we observe aninversion of the funnel-like structure, with the reverse funnel shapebeing more and more pronounced as γ increases. This results in alow penalty for maps between low frequency eigenfunctions, whichdoes not respect the shape of the ground-truth maps (see Fig. 5).Later in this paper, we report empirically obtained optimal valuesfor γ on a benchmark dataset (see Fig. 14). In summary, we propose the following energy to compute a func-tional map for a pair of shapes, with a set of given descriptors: E (cid:0) C (cid:1) = α E desc + α E mult + α E orient + α E resolvent (15)where E desc is the descriptor-preserving term defined in Sec. 3, E mult is the descriptor-commutativity term defined as E mult = ∑ i (cid:13)(cid:13) C D mult1 i − D mult2 i C (cid:13)(cid:13) introduced in [NO17], E orient is theorientation-preserving term introduced in [RPWO18]. E resolvent isthe our new term introduced in Eq. (6) and discussed above.As mentioned above, when functional maps are expressed in theLaplace-Beltrami eigenbasis, this term can be written via a penaltyusing the mask matrix given in Eq. (13) .Before proceeding to a more extensive evaluation of the pro-posed energy, in Fig. 2 above, we provide an example of a func-tional map obtained using this energy and compare it to one ob-tained using the standard Laplacian commutator regularizer. Fig. 2also shows the resolvent mask with γ = .
5. Results
We tested the proposed approach using a MATLAB-based imple-mentation, which we adapt from the state-of-the-art functional mapapproach in [RPWO18]. Here, we first describe the benchmarkdatasets and the baseline methods.
Datasets.
We use the two datasets introduced in [RPWO18].These consist of shapes from the FAUST [BRLB14] andTOSCA [BBK08] datasets, which were remeshed so that the shapeshave different triangulations, and are no longer in one-to-one cor-respondence, making the matching more challenging and realistic.Specifically, we include 300 FAUST shape pairs and 284 TOSCApairs for evaluation.
Baselines.
To evaluate our new regularizer, we compare to the fol-lowing two masks: c (cid:13) (cid:13) . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations
50 100 150 200 2500 . . . . k : functional map size A v e r a g ee rr o r FAUST per-vertex measure standard slanted ours
50 100 150 200 2500 . . . k FAUST direct measure
Initialization with ICP 50 100 150 200 2500 . . . . . k TOSCA per-vertex measure
50 100 150 200 2500 . . . . k TOSCA direct measure
Figure 7:
Changing the functional map size. We randomly select 50 FAUST non-isometric pairs and 50 TOSCA isometric pairs. For eachof the pairs, we only use one pair of WKS [ASC11] descriptors. We then optimize for a functional map with different size ranging from 20to 250 using different Laplacian mask terms. The solid lines represent the initialization with different masks and the dashed lines are theICP-refined results. We report the per-vertex and the direct error measure. We can see that the proposed mask is much more stable than thestandard and the slanted mask as the functional map size increases. • "standard": the standard Laplacian-commutativity mask, whichis defined as M LB ( i , j ) = (cid:0) ∆ ( i ) − ∆ ( j ) (cid:1) as discussed before. • "slanted": the heuristic slanted diagonal penalty mask proposedin [RCB ∗ M ( i , j ) = exp (cid:18) − η (cid:113) i + j (cid:19) (cid:13)(cid:13)(cid:13)(cid:13) n (cid:107) n (cid:107) × (cid:16) ( i , j ) T − p (cid:17)(cid:13)(cid:13)(cid:13)(cid:13) where p = ( , ) T , and n = ( , r / k ) T is the line direction withslope r / k , where r is the estimated rank of the functional map,and k is the size of the square functional map. This weight matrixis originally applied to partial shape matching in [RCB ∗ η is set to the default value 0.03 as suggested in theoriginal paper. This mask is illustrated on Fig. 5. Note that itexhibits the desired funnel-like structure for the lower part of thespectrum, but not the upper part.We also compare three different settings: • Initialization . The wave kernel signatures (WKS) [ASC11] areused to construct E desc , E mult , and E orient in Eq. (15). Then weoptimize the functional map w.r.t. the energy defined in Eq. (15)with three different masks, namely, the standard, slanted, and ourresolvent mask. • ICP refinement . After initialization, we use ICP [OBCS ∗
12] torefine the computed functional maps. • BCICP refinement . After initialization, we use the recently pro-posed BCICP algorithm [RPWO18] to refine the computed func-tional and pointwise maps, using the open-source implementa-tion and default parameters provided by the authors.
Measurements.
In our experiments, we measured the quality of thefunctional maps and the recovered point-wise maps: • Point-wise maps . Since most shapes contain left-right symme-tries, which are indistinguishable for intrinsic methods, in eachdataset, we considered both the ground-truth direct and symmet-ric correspondences. To measure the accuracy of a computedmap, we used the following measures:– per-vertex error : for each vertex we accept the ground-truthdirect and symmetric correspondences and take the minimum Given descriptorGround-truth map o u r s s l a n t e d s t a nd a r d k = k = k = k = Figure 8:
Given one pair of WKS descriptors, as visualized on theleft top, we use different Laplacian mask terms to optimize for afunctional map with size k. The ground-truth map is visualized onthe bottom left. We can see that our mask is much more stable overdifferent size k. as the error of this vertex. This measure reflects the accuracyof the map regardless of the symmetry.– direct error : we compute the average per-vertex error to thedirect ground-truth correspondences only. This measure re-flects both the accuracy and the smoothness of the map. • Functional maps . We can evaluate a functional map by the qual-ity of its recovered point-wise map, or by measuring the penaltyfrom a given mask. Specifically, for a given functional map C and a mask matrix M , we can measure the total penalty as ∑ i j [ M ] i j [ C ] i j . For the functional maps pipeline, a set of corresponding descriptorsis given as input. We then optimize a k × k functional map by mini-mizing an objective function based on Eq. (15). Therefore, we haveto solve for k variables. If k is a smaller value, e.g., k <
50, the op- c (cid:13) (cid:13)(cid:13)
50, the op- c (cid:13) (cid:13)(cid:13) . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations Geodesic Error ( × − ) % c o rr e s pond e n ce FAUST standard (1) : 201.7slanted (2) : 188.2 ours (3) : 112.8 (1) + ICP : 116.5(2) + ICP : 141.3 (3) + ICP : 66.0 (1) + BCICP : 56.9(2) + BCICP : 116.5 (3) + BCICP : 47.6
Geodesic Error ( × − ) % c o rr e s pond e n ce TOSCA standard (1) : 213.0slanted (2) : 190.6 ours (3) : 125.7 (1) + ICP : 156.5(2) + ICP : 163.2 (3) + ICP : 90.2 (1) + BCICP : 81.8(2) + BCICP (3) + BCICP : 61.5
Figure 9:
For each dataset, we compare the quality of the functional maps between the standard mask (blue curves), the slanted mask (yellowcurves), and our resolvent mask (red curves) with the same set of parameters, where a single pair of descriptors is used to optimize a 100-by-100 functional map. The comparison is made in three different settings: comparing the initialization directly (dotted lines), comparingthe initialization with ICP refinement (dashed lines), and with the BCICP refinement (solid lines). Specifically, these curves are measured on300 FAUST shape pairs and 284 TOSCA pairs. The average direct errors are reported in the legends. timization problem is easier to solve since the number of variablesis small. However, if k is too small, the information that is encodedinto the optimized functional map is limited to the low-resolutionof the spectrum of the shapes, and hence it will be hard to transferdetailed information. On the other hand, if k is too large, solving theoptimization problem is potentially more time-consuming. Evenmore importantly, this optimization problem can become under-constrained when the number of variables exceeds the constraintsstemming from the input descriptors. In this case, we need effectiveregularizers to regularize the functional maps. In real experiments,the choice of the parameter k is a key hyper parameter.To quantify the stability and the effectiveness of the proposedresolvent mask, we conduct the following test: for each of the testpairs, we only use one pair of corresponding WKS descriptors. Wethen fix this descriptor pair and optimize for functional maps withdifferent sizes ranging from 20 to 250. We randomly select 50 pairsof FAUST and 50 pairs of TOSCA, and report the average errorover the tested shape pairs w.r.t. different functional map size inFig. 7. We can see that the standard mask fails to regularize thefunctional map with a large size: the average per-vertex error isthree or four times larger than ours. At the same time, the slantedmask has a better performance than the standard one in the per-vertex measure. However, the slanted mask has large direct errors,which suggests that the smoothness is not well preserved. We be-lieve this is due to the fact that the orientation-preserving regular-izer starts to fail to disambiguate the symmetry as the functionalmap size k increases. In this case, the slanted mask cannot helpthe orientation-preserving operator, while our mask can strengthenthe functionality of the orientation-preserving operator and leadsto maps with much lower per-vertex and direct error. In summary,even with limited constraints from a single pair of WKS descrip-tors, increasing the number of variables does not significantly affectthe performance of our mask. Fig. 8 shows an illustrative example.This test shows that our mask is much more stable and can betterregularize larger functional maps even in very challenging cases 1 2 5 10 200.050.1 A v e r a g e d i r ec t e rr o r standard (1)slanted (2)ours (3)(1) + ICP(2) + ICP(3) + ICP Figure 10:
Changing the number of input descriptors. We use dif-ferent numbers of descriptors to optimize a 100-by-100 functionalmap with different masks. The results are reported for 300 FAUSTshape pairs. As the number of descriptors increases, the results im-prove for all the masks. with little input information or constraints. Also, it suggests thatwith this new mask, we no longer need to tune the parameter k asmuch as needed by the standard mask to achieve a better result. In this experiment, we compare our resolvent mask to the stan-dard and the slanted mask on a larger FAUST and TOSCA datasetin three settings: the initialization, with ICP refinement, and withBCICP refinement. To make a fair comparison between differentmasks, the weights for different terms (the α i in Eq. (15)) are fixedacross different test pairs and different test masks.In this test, for each test pair, we use one pair of WKS descriptorsto optimize a 100-by-100 functional map. As reported in Fig. 9, ourmask leads to 40.4%, 43.1%, and 16.3% improvement w.r.t the set-tings of the initialization, with ICP refinement, and with BCICP re-finement respectively, over the best of the standard and the slanted c (cid:13) (cid:13) . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations Source Standard Slanted
Ours
Ground-truth I n iti a li za ti on + I C P Figure 11:
Example. Comparing the quality of the initial maps andafter ICP refinement via texture transfer. (One pair of descriptorsis used to optimize a 50-by-50 functional map)
Source Standard Slanted
Ours
Ground-truth
Figure 12:
Example. Comparing the quality of the ICP refinedmaps via color transfer. (One pair of descriptors is used to opti-mize a 100-by-100 functional map) mask on the FAUST dataset regarding the average direct error. Sim-ilarly, for the TOSCA dataset, ours achieves 34.2%, 42.3%, and24.8% improvement respectively over the best of the standard andthe slanted mask.Fig. 10 shows a test on the FAUST dataset, where we use dif-ferent number of carefully curated input descriptors based on WKS[ASC11] using parameters from [RPWO18] to compute 100-by-100 functional maps. As the number of input descriptors increases,more constraints are added to regularize the functional map. There-fore, for all the masks, and all the test settings, the results improve.Observe that when the number of variables is small (as shown inFig. 7 for small k ) or the number of descriptor constraints is large(as shown in Fig. 10 for large descriptor number), all masks per-form well since the problem is well-constrained. However, whenthis is not the case, our resolvent mask can still regularize thefunctional map better than the other two. For completeness, in Ap-pendix D we also include results with the BCICP refinement. Weremark that this refinement is very computationally and memoryintensive due, in part, to requiring all-pairs geodesic distances, butcan, as such, improve upon even very poor quality maps.Fig. 11 shows a qualitative example of using a single descriptorpair to optimize for a 100 ×
100 functional map. The first row showsthe quality of the initial maps with different masks, and the second Source Standard Slanted
OursFigure 13:
Example. Here we show a challenging pair of a lion anda cat and the results are refined by BCICP. We can see that, whenthe initial maps are in a low quality and the BCICP refinementfails to improve the initial maps, our new mask still gives a morereasonable map. (Twenty pairs of descriptors are used to optimizea 50-by-50 functional map) . . . . . γ A v e r a g ee rr o r Different γ (cid:107) M Re ( γ ) (cid:12) C (cid:107) F + (cid:107) M Im ( γ ) (cid:12) C (cid:107) F standardstandard + ICP oursours + ICP . . . . . . . w Different weight wM res = wM Im + ( − w ) M Re Figure 14:
Left: changing the parameter γ (as defined in Eq. (6) );Right: changing the relative weight w between the imaginary partand the real part (see the original definition in Eq. (13) ). The resultsare based on 50 TOSCA pairs. rows the quality of the corresponding maps refined by ICP. We cansee that our resolvent mask outperforms the other two masks. Alsothe quality of our initial map is close to the ICP refined map, whichshows that with our mask, we do not rely on the post-processingrefinement as much as the other two. Fig. 12 shows another ex-ample of the results refined by ICP. Fig. 13 shows a challengingnon-isometric example with results refined after BCICP. As shown in Fig. 5, the funnel-shape of our resolvent mask alignswell with the ground-truth functional map, which leads to a bet-ter performance over the standard and the slanted mask. To furtheranalyze the properties of our mask, we also conduct the followingexperiments: We test the range of the parameter γ , and the rela-tive weight between the complex and the real part to construct theresolvent mask. Moreover, we investigate the correlation betweenthe mask penalty added on a functional map and the correspondingrecovered point-wise map. In this section, we empirically tune the parameters in our resolventbased mask. First, we explore different values of γ . Recall that γ controls the funnel-like structure of the mask. Thus, it is expectedthat tuning γ can influence the functional map quality. c (cid:13) (cid:13)(cid:13)
Left: changing the parameter γ (as defined in Eq. (6) );Right: changing the relative weight w between the imaginary partand the real part (see the original definition in Eq. (13) ). The resultsare based on 50 TOSCA pairs. rows the quality of the corresponding maps refined by ICP. We cansee that our resolvent mask outperforms the other two masks. Alsothe quality of our initial map is close to the ICP refined map, whichshows that with our mask, we do not rely on the post-processingrefinement as much as the other two. Fig. 12 shows another ex-ample of the results refined by ICP. Fig. 13 shows a challengingnon-isometric example with results refined after BCICP. As shown in Fig. 5, the funnel-shape of our resolvent mask alignswell with the ground-truth functional map, which leads to a bet-ter performance over the standard and the slanted mask. To furtheranalyze the properties of our mask, we also conduct the followingexperiments: We test the range of the parameter γ , and the rela-tive weight between the complex and the real part to construct theresolvent mask. Moreover, we investigate the correlation betweenthe mask penalty added on a functional map and the correspondingrecovered point-wise map. In this section, we empirically tune the parameters in our resolventbased mask. First, we explore different values of γ . Recall that γ controls the funnel-like structure of the mask. Thus, it is expectedthat tuning γ can influence the functional map quality. c (cid:13) (cid:13)(cid:13) . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations .
02 0 .
05 0 . . . . . . M a s kp e n a lt y (cid:0) ∑ ij [ M ] ij [ C ] ij (cid:1) StandardSlantedOurs
Figure 15:
Correlation between the mask penalty on a functionalmap and the accuracy of the corresponding point-wise map. Fora TOSCA isometric pair, we sample 700 functional maps of size50-by-50 with different quality (normalized to the same scale). Wemeasure the penalty of different masks, as an indicator of the qual-ity of the functional map. We also measure the average geodesic er-ror (w.r.t. the direct measurement) of the recovered point-wise mapof the corresponding functional map. The correlation between thefunctional map penalty and the point-wise map quality is visualizedas a scatter plot with 700 samples. Compared to the standard andthe slanted mask, our resolvent mask applies a smaller penalty tofunctional maps with good underlying point-wise maps and moreheavily penalizes functional maps with bad underlying point-wisemaps.
On Fig. 14, we report results for 100 FAUST pairs, where theresult of the standard mask is colored blue, and ours is colored red.As is shown on the left, when γ lies between 0 and 1, our maskalways outperforms the standard mask over both the initialization(solid lines) and after the ICP refinement (the dashed lines).Note that the case γ = γ = γ > M Re and M Im to have different weights, ( − w ) and w for w ∈ [ , ] ,respectively. As shown in Fig. 14 (right), the convex shape of thered curve (where the weight w changes) demonstrates that both thereal part and the imaginary part contribute to the improvement overthe standard mask. Note that our mask with any convex combina-tion between the real part and the imaginary parts outperforms thestandard mask. In practice we always use the equal weight w = . To justify that our resolvent mask is a better regularizer than thestandard and the slanted ones, we also measure how the maskpenalty relates to the accuracy of the corresponding point-wisemap shown in Fig. 15. In this experiment, we generate 700 differ-ent point-wise maps with different levels of accuracy, then convert − − . . . per-vertex measure − − . . . . direct measureStandardSlantedOurs Figure 16:
Changing the relative weight α of the mask term. Wechoose 68 different values for α in the range of 0 and 10 for thistest. Table 1:
Average error on 284 TOSCA pairs of different masks,each mask using its own optimal weight from Fig. 16. Specifically,we set α to − , − , − for the standard, slanted, and ourmask respectively. Masktype Average error ( × − ) per-vertex measure direct measureIni + ICP + BCICP Ini + ICP + BCICPstandard 87.26 65.97 41.22 178.8 130.8 76.87slanted 64.66 50.76 37.45 166.9 131.5 100.0 ours 55.52 43.82 32.48 124.5 90.6 62.33 them to a functional map representation and measure the penaltyadded by different masks w.r.t. the average geodesic error computedfrom the pointwise maps. Each scatter point in Fig. 15 shows sucha test sample.We can observe that, compared to the standard and the slantedmask, the new mask induces a lower penalty on functional mapswith a good quality (i.e., smaller average geodesic error), and pe-nalizes a functional map with larger error more heavily. This furtherconfirms that using our resolvent mask is more likely to produce abetter functional, and ultimately better pointwise map. In our tests, we use γ = . w = . η = .
03 to construct the slantedmask. When comparing the three masks in the initialization setting,i.e., to optimize the energy defined in (15), the weights α i are setto the same values as reported in [RPWO18] but in a relative way(see Appendix D for more details). As for the comparison of theICP and BCICP refinement settings, the same default parametersand the same number of iterations are used for different masks.Fig. 16 shows the results of the changing the weight of the maskterm, i.e., α , while keeping the rest weights fixed on 50 TOSCAshape pairs. We can see that, for all different choices of the weight,our resolvent mask is always better than the standard mask. Whenthe weight is small enough, it seems the slanted mask can achievea faster decrease of the error than ours. However, the range of the c (cid:13) (cid:13) . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations Source Standard Slanted
Ours
Ground-truth I n iti a li za ti on + I C P Figure 17:
Example. Comparing the quality of the initial mapscomputed with different masks and after ICP refinement via texturetransfer. Each mask uses its own optimal weight from Fig. 16. effective weight is smaller than ours besides that the constructionis purely heuristic and lacks a theoretical justification.On the other hand, Fig. 16 also suggests that different maskshave their own preferred choice of the weight α . Therefore, in-stead of fixing α , we set it independently w.r.t. the optimal valueson a subset reflected in Fig. 16. Specifically, for the standard mask,we set α = − , for the slanted mask, we set α = − , andfor ours, we set α = − . The corresponding average errors onthe complete TOSCA dataset are reported in Table 1. We can seethat, even when the parameters are carefully tuned for the other twomasks, our resolvent mask still achieves the best accuracy. Fig. 17shows a qualitative example. Source S t a nd a r d S l a n t e d O u r s Figure 18:
Example. Comparing the quality of the initial mapscomputed with different masks of four SHREC shape pairs via tex-ture transfer.
To show the usefulness of our resolvent mask on non-isometricshape pairs, we test the 20 FourLeg shapes from the SHREC 2007dataset [GBP]. Specifically, we use 10 pairs of descriptors con-structed from 4 landmarks to optimize a 120-by-120 functionalmap using the standard, the slanted, and our resolvent mask. For afair comparison, we use the same parameters as the previous tests:we set γ = . η = .
03 for theslanted mask. Fig. 18 shows a qualitative example. We can see that . . . A v e r a g ee rr o r StandardSlantedOurs
Figure 19:
For 380 SHREC shape pairs, we plot the average errorfor each pair with different masks, where the result of the standard,slanted and our mask is colored in blue, yellow, and red respec-tively. Therefore, in the case where the red curve is below the blueor the yellow points, our mask leads to a lower error for these pairs.
Source Standard Slanted
OursFigure 20:
Failure case. Here we show a challenging shape pairfrom SHREC, where all the three masks fail to produce a good map. our resolvent mask gives the best initialization. We also measurethe accuracy of the 380 point-wise maps among the 20 shapes onthe given landmarks (21 ground-truth landmarks are given for eachshape). The average error for the standard, slanted, and our maskis 0.114, 0.122, and 0.107 respectively. Fig. 19 reports the averageerror for each shape pair. We observe that resolvent mask gives abetter result than the standard mask on 68.4% out of 380 pairs, andoutperforms the slanted mask on 69.5% pairs. The limited quanti-tative improvement of our mask over the other two masks is due to:(1) for the shape pairs where our mask significantly outperformsthe other two, e.g., as shown in Fig. 18, the map quality is mea-sured on only 21 landmarks, where the average error does not fullyreflect our improvement. (2) for some extremely challenging pairs,e.g., a failure case shown in Fig. 20, all the maps from differentmasks have a poor quality, and thus the “relative improvement” onthe average error is not informative. Moreover, we also computedthe average error in the case where the mask term is removed fromthe total energy, i.e., set α =
0. In this case, the average error overthe complete dataset is 0.112. We can see that, the standard andthe slanted mask can have a negative effective in this case, whileour resolvent mask still works to improve the map quality to someextent.
6. Conclusion, Limitations and Future Work
In this paper, we proposed a new regularizer, the resolvent Lapla-cian commutativity, for the functional map pipeline. We first ana-lyzed the limitations of the original Laplacian commutativity termand theoretically justified the effectiveness of our proposed newterm. This new regularizer can significantly improve the quality c (cid:13) (cid:13)(cid:13)
In this paper, we proposed a new regularizer, the resolvent Lapla-cian commutativity, for the functional map pipeline. We first ana-lyzed the limitations of the original Laplacian commutativity termand theoretically justified the effectiveness of our proposed newterm. This new regularizer can significantly improve the quality c (cid:13) (cid:13)(cid:13) . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations of the computed functional maps and the corresponding recoveredpoint-wise maps before or after refinement.However, our method also has several limitations that we wouldlike to overcome in future work. First, our proposed regularizer iswell justified on isometric and non-isometric shape pairs, but not onpartial shapes, where the ground-truth functional maps can have adifferent structure. Therefore, it would be interesting to extend theanalysis to partial shape pairs. Second, besides the funnel pattern,in our experiments, we also observed the slanted-diagonal structureof the ground-truth functional map of some non-isometric shapepairs as discussed in [RCB ∗ Acknowledgement
The authors would like to thank the anony-mous reviewers for their valuable comments and helpful sugges-tions. Parts of this work were supported by the KAUST OSRAwards No. CRG2017-3426 and CRG2018-3730, a gift from theNVIDIA Corporation, and the ERC Starting Grant No. 758800(EXPROTEA).
References [ADK16] A
FLALO
Y., D
UBROVINA
A., K
IMMEL
R.: Spectral general-ized multi-dimensional scaling.
International Journal of Computer Vi-sion 118 , 3 (2016), 380–392. 2[ASC11] A
UBRY
M., S
CHLICKEWEI
U., C
REMERS
D.: The wave kernelsignature: A quantum mechanical approach to shape analysis. In (2011), IEEE, pp. 1626–1633. 7, 9[BBK08] B
RONSTEIN
A. M., B
RONSTEIN
M. M., K
IMMEL
R.:
Nu-merical Geometry of Non-Rigid Shapes . Springer Science & BusinessMedia, 2008. 6[BCBB16] B
IASOTTI
S., C
ERRI
A., B
RONSTEIN
A., B
RONSTEIN
M.:Recent trends, applications, and perspectives in 3d shape similarity as-sessment. In
Computer Graphics Forum (2016), vol. 35, pp. 87–119.2[BDK17] B
URGHARD
O., D
IECKMANN
A., K
LEIN
R.: Embeddingshapes with Green’s functions for global shape matching.
Computers& Graphics 68 (2017), 1–10. 2[BRLB14] B
OGO
F., R
OMERO
J., L
OPER
M., B
LACK
M. J.: FAUST:Dataset and evaluation for 3D mesh registration. In
Proceedings IEEEConf. on Computer Vision and Pattern Recognition (CVPR) (Piscataway,NJ, USA, June 2014), IEEE. 1, 6[CSBK18] C
HOUKROUN
Y., S
HTERN
A., B
RONSTEIN
A. M., K
IMMEL
R.: Hamiltonian operator for spectral shape analysis.
IEEE transactionson visualization and computer graphics (2018). 3[Cut13] C
UTURI
M.: Sinkhorn distances: Lightspeed computation of op-timal transport. In
Advances in neural information processing systems (2013), pp. 2292–2300. 2[Dod81] D
ODZIUK
J.: Eigenvalues of the Laplacian and the heat equa-tion.
The American Mathematical Monthly 88 , 9 (1981), 686–695. 4[EBC17] E
ZUZ
D., B EN -C HEN
M.: Deblurring and denoising of mapsbetween shapes. In
Computer Graphics Forum (2017), vol. 36, WileyOnline Library, pp. 165–174. 2, 3 [ERGB16] E
YNARD
D., R
ODOLA
E., G
LASHOFF
K., B
RONSTEIN
M. M.: Coupled functional maps. In
3D Vision (3DV), 2016 FourthInternational Conference on (2016), IEEE, pp. 399–407. 2[GBKS18] G
EHRE
A., B
RONSTEIN
M., K
OBBELT
L., S
OLOMON
J.: In-teractive curve constrained functional maps. In
Computer Graphics Fo-rum (2018), vol. 37, Wiley Online Library, pp. 1–12. 2[GBP] G
IORGI
D., B
IASOTTI
S., P
ARABOSCHI
L.: Shape retrieval con-test 2007: Watertight models track. 11[HO17] H
UANG
R., O
VSJANIKOV
M.: Adjoint map representation forshape analysis and matching. In
Computer Graphics Forum (2017),vol. 36, Wiley Online Library, pp. 151–163. 2[HRS ∗
16] H
EEREN
B., R
UMPF
M., S
CHRÖDER
P., W
ARDETZKY
M.,W
IRTH
B.: Splines in the space of shells. In
Computer Graphics Forum (2016), vol. 35, Wiley Online Library, pp. 111–120. 1[HWG14] H
UANG
Q., W
ANG
F., G
UIBAS
L.: Functional map networksfor analyzing and exploring large shape collections.
ACM Transactionson Graphics (TOG) 33 , 4 (2014), 36. 1, 2[KBB ∗
13] K
OVNATSKY
A., B
RONSTEIN
M. M., B
RONSTEIN
A. M.,G
LASHOFF
K., K
IMMEL
R.: Coupled quasi-harmonic bases. In
Com-puter Graphics Forum (2013), vol. 32, pp. 439–448. 1, 2[KBBV15] K
OVNATSKY
A., B
RONSTEIN
M. M., B
RESSON
X., V AN - DERGHEYNST
P.: Functional correspondence by matrix completion.In
Proceedings of the IEEE conference on computer vision and patternrecognition (2015), pp. 905–914. 2[KGB16] K
OVNATSKY
A., G
LASHOFF
K., B
RONSTEIN
M. M.:MADMM: a generic algorithm for non-smooth optimization on mani-folds. In
Proc. ECCV (2016), Springer, pp. 680–696. 3[LRB ∗
16] L
ITANY
O., R
ODOLÀ
E., B
RONSTEIN
A. M., B
RONSTEIN
M. M., C
REMERS
D.: Non-rigid puzzles. In
Computer Graphics Forum (2016), vol. 35, Wiley Online Library, pp. 135–143. 2[LRBB17] L
ITANY
O., R
ODOLÀ
E., B
RONSTEIN
A. M., B
RONSTEIN
M. M.: Fully spectral partial shape matching. In
Computer GraphicsForum (2017), vol. 36, Wiley Online Library, pp. 247–258. 2[MCSK ∗
17] M
ANDAD
M., C
OHEN -S TEINER
D., K
OBBELT
L., A L - LIEZ
P., D
ESBRUN
M.: Variance-minimizing transport plans for inter-surface mapping.
ACM Transactions on Graphics (TOG) 36 (2017), 14.2[MP49] M
INAKSHISUNDARAM
S., P
LEIJEL
A.: Some properties of theeigenfunctions of the Laplace-operator on Riemannian manifolds.
Cana-dian J. Math 1 , 8 (1949). 3, 14[MRCB18] M
ELZI
S., R
ODOLÀ
E., C
ASTELLANI
U., B
RONSTEIN
M. M.: Localized manifold harmonics for spectral shape analysis.In
Computer Graphics Forum (2018), vol. 37, Wiley Online Library,pp. 20–34. 3[NMR ∗
18] N
OGNENG
D., M
ELZI
S., R
ODOLÀ
E., C
ASTELLANI
U.,B
RONSTEIN
M., O
VSJANIKOV
M.: Improved functional mappings viaproduct preservation. In
Computer Graphics Forum (2018), vol. 37, Wi-ley Online Library, pp. 179–190. 2[NO17] N
OGNENG
D., O
VSJANIKOV
M.: Informative descriptor preser-vation via commutativity for shape matching.
Computer Graphics Forum36 , 2 (2017), 259–267. 2, 3, 6[NVT ∗
14] N
EUMANN
T., V
ARANASI
K., T
HEOBALT
C., M
AGNOR
M., W
ACKER
M.: Compressed manifold modes for mesh processing.In
Computer Graphics Forum (2014), vol. 33, Wiley Online Library,pp. 35–44. 3[OBCS ∗
12] O
VSJANIKOV
M., B EN -C HEN
M., S
OLOMON
J.,B
UTSCHER
A., G
UIBAS
L.: Functional Maps: A Flexible Repre-sentation of Maps Between Shapes.
ACM Transactions on Graphics(TOG) 31 , 4 (2012), 30. 1, 2, 4, 7[OCB ∗
17] O
VSJANIKOV
M., C
ORMAN
E., B
RONSTEIN
M., R
ODOLÀ
E., B EN -C HEN
M., G
UIBAS
L., C
HAZAL
F., B
RONSTEIN
A.: Com-puting and processing correspondences with functional maps. In
ACMSIGGRAPH 2017 Courses (2017), SIGGRAPH ’17, pp. 5:1–5:62. 2, 3 c (cid:13) (cid:13) . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations [PSO18] P OULENARD
A., S
KRABA
P., O
VSJANIKOV
M.: Topologi-cal function optimization for continuous shape matching. In
ComputerGraphics Forum (2018), vol. 37, Wiley Online Library, pp. 13–25. 2[RCB ∗
17] R
ODOLÀ
E., C
OSMO
L., B
RONSTEIN
M. M., T
ORSELLO
A., C
REMERS
D.: Partial functional correspondence. In
ComputerGraphics Forum (2017), vol. 36, Wiley Online Library, pp. 222–236. 1,2, 4, 7, 12[RMC15] R
ODOLÀ
E., M
OELLER
M., C
REMERS
D.: Point-wise maprecovery and refinement from functional correspondence. In
Proc. Vi-sion, Modeling and Visualization (VMV) (2015). 3[ROA ∗
13] R
USTAMOV
R., O
VSJANIKOV
M., A
ZENCOT
O., B EN -C HEN
M., C
HAZAL
F., G
UIBAS
L.: Map-based exploration of intrinsicshape differences and variability.
ACM Transactions on Graphics (TOG)32 , 4 (July 2013), 72:1–72:12. 1, 2[RPWO18] R EN J., P
OULENARD
A., W
ONKA
P., O
VSJANIKOV
M.:Continuous and orientation-preserving correspondences via functionalmaps.
ACM Transactions on Graphics (TOG) 37 , 6 (2018). 2, 6, 7, 9,10, 14, 15[RS80] R
EED
M., S
IMON
B.:
Methods of Modern Mathematical PhysicsI: Functional Analysis , revised and enlarged ed. Academic Press, SanDiego, 1980. 4, 13[Sau06a] S
AUVIGNY
F.:
Partial Differential Equations 1: Foundationsand Integral Representations . Springer, Berlin, Heidelberg, 2006. 4[Sau06b] S
AUVIGNY
F.:
Partial Differential Equations 2: Functional An-alytic Methods . Springer, Berlin, Heidelberg, 2006. 4[SDGP ∗
15] S
OLOMON
J., D E G OES
F., P
EYRÉ
G., C
UTURI
M.,B
UTSCHER
A., N
GUYEN
A., D U T., G
UIBAS
L.: Convolutionalwasserstein distances: Efficient optimal transportation on geometric do-mains.
ACM Transactions on Graphics (TOG) 34 , 4 (2015), 66. 2[SPKS16] S
OLOMON
J., P
EYRÉ
G., K IM V. G., S RA S.: Entropic metricalignment for correspondence problems.
ACM Transactions on Graphics(TOG) 35 , 4 (2016), 72. 2[TCL ∗
13] T AM G. K., C
HENG
Z.-Q., L AI Y.-K., L
ANGBEIN
F. C., L IU Y., M
ARSHALL
D., M
ARTIN
R. R., S UN X.-F., R
OSIN
P. L.: Regis-tration of 3D point clouds and meshes: a survey from rigid to nonrigid.
IEEE TVCG 19 , 7 (2013), 1199–1217. 2[VKZHCO11] V AN K AICK
O., Z
HANG
H., H
AMARNEH
G., C
OHEN -O R D.: A survey on shape correspondence. In
Computer GraphicsForum (2011), vol. 30, pp. 1681–1707. 2[VLR ∗
17] V
ESTNER
M., L
ITMAN
R., R
ODOLÀ
E., B
RONSTEIN
A.,C
REMERS
D.: Product manifold filter: Non-rigid shape correspondencevia kernel density estimation in the product space. In
Proc. CVPR (2017),pp. 6681–6690. 2[WGBS18] W
ANG
L., G
EHRE
A., B
RONSTEIN
M. M., S
OLOMON
J.:Kernel functional maps. In
Computer Graphics Forum (2018), vol. 37,Wiley Online Library, pp. 27–36. 2[WHG13] W
ANG
F., H
UANG
Q., G
UIBAS
L. J.: Image co-segmentationvia consistent functional maps. In
Proceedings of the IEEE InternationalConference on Computer Vision (2013), pp. 849–856. 2[WLZT18] W
ANG
Y., L IU B., Z
HOU
K., T
ONG
Y.: Vector field maprepresentation for near conformal surface correspondence. In
ComputerGraphics Forum (2018), vol. 37, Wiley Online Library, pp. 72–83. 2
Appendix A:
Pullbacks are Bounded
Lemma 1
Let Φ : M → N be a diffeomorphism between con-nected compact oriented Riemannian manifolds. Then, the associ-ated pullback Φ ∗ is bounded as a map L ( N ) → L ( M ) . Proof
Let d M and d N denote the volume forms of M and N , respectively. Since C ∞ is dense in L , it is enough to show thatthere exists a constant B > f ∈ C ∞ ( N ) , (cid:90) M (cid:0) Φ ∗ f (cid:1) d M ≤ B (cid:90) N f d N . (16)By virtue of being a diffeomorphism, Φ is invertible. Using thepullback of Φ − , we can express the left-hand side of the desiredinequality as an integral over N .Begin by considering (cid:104)(cid:16) Φ − (cid:17) ∗ d M (cid:105) , the pullback of the volumeform of M . Since volume forms are top degree forms, there exists u ∈ C ∞ ( N ) such that (cid:104)(cid:16) Φ − (cid:17) ∗ d M (cid:105) = ud N (17) Φ is either orientation preserving or orientation reversing. In theformer case, u >
0. In the latter, u <
0. In either case, the left-handside of the desired inequality can be recast as (cid:90) M (cid:0) Φ ∗ f (cid:1) d M = (cid:90) N (cid:16)(cid:16) Φ − (cid:17) ∗ Φ ∗ f (cid:17) | u | d N = (cid:90) N f | u | d N . (18)It remains to bound this expression. (cid:90) N f | u | d N ≤ sup x ∈N | u ( x ) | (cid:90) N f d N (19)Since N is compact and u is continuous, the supremum is achievedand is finite. This concludes the proof. Appendix B:
Bounded Frobenius Norm for γ > / Definition 2 (Hilbert-Schmidt Norm)
Let A : H → H be a linearoperator between Hilbert spaces. Let A † denote the adjoint of theoperator A . Then, the Hilbert-Schmidt norm of A is given by: (cid:107) A (cid:107) HS = Tr (cid:16) A † A (cid:17) = ∞ ∑ i = (cid:104) e i , A † Ae i (cid:105) , (20)where { e i } ∞ i = is any orthonormal basis of H . This norm is alsoknown as the Schatten 2-norm.The following lemma is the main result of this appendix. Lemma 2
Let ∆ and ∆ be Laplacians on compact, connected,oriented surfaces M and M , respectively. Let C : L ( M ) → L ( M ) be a bounded operator. If γ > /
2, then: (cid:13)(cid:13) C R µ (cid:0) ∆ γ (cid:1) − R µ (cid:0) ∆ γ (cid:1) C (cid:13)(cid:13) HS < ∞ , (21)where µ is any complex number not on the non-negative real line. Proof
Operators with finite Hilbert-Schmidt norm are known asoperators of Hilbert-Schmidt class. It can be shown (see [RS80])that linear combinations of Hilbert-Schmidt class operators are of c (cid:13) (cid:13) . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations Hilbert-Schmidt class. Moreover, the product of a bounded opera-tor and a Hilbert-Schmidt class operator is also of Hilbert-Schmidtclass. Thus, it is sufficient to show that R µ ( ∆ γ ) has finite Hilbert-Schmidt norm for γ > / ∆ by { ψ k } ∞ k = and { λ k } ∞ k = , respectively. We assume that the eigenvalues are num-bered in the usual non-decreasing order. R µ ( ∆ γ ) is diagonal in the eigenbasis of ∆ and has eigenvalues1 / ( λ γ k − µ ) . Thus, the Hilbert-Schmidt norm of R µ ( ∆ γ ) is given by: (cid:107) R µ ( ∆ γ ) (cid:107) HS = ∞ ∑ k = (cid:104) ψ k , R µ ( ∆ γ ) † R µ ( ∆ γ ) ψ k (cid:105) = ∞ ∑ k = (cid:12)(cid:12) λ γ k − µ (cid:12)(cid:12) (22)We will establish the convergence of this series for γ > / ∑ ∞ k = k p . This series is well-known to converge if and only if p > λ k are non-negative and increasing with k towards ∞ , thefollowing inequality holds for all large enough k :14 λ γ k ≤ (cid:12)(cid:12) λ γ k − µ (cid:12)(cid:12) ≤ λ γ k . (23)By Weyl’s law, there exists a constant B > λ k ∼ Bk forlarge k [MP49]. Then, the inequality becomes14 B γ k γ ≤ (cid:12)(cid:12) λ γ k − µ (cid:12)(cid:12) ≤ B γ k γ . (24)Thus, by the comparison test, the series for (cid:107) R µ ( ∆ γ ) (cid:107) HS convergesif and only if γ > /
2. This concludes the proof.We conclude this appendix with two remarks on the above re-sult. First, note that we have shown a result slightly stronger thanrequired by the statement of the lemma. In fact, we only neededto show that the series for (cid:107) R µ ( ∆ γ ) (cid:107) HS converges if γ > /
2. The"only if" part was optional. We have done this to illustrate that theabove proof strategy is guaranteed to fail for γ ≤ /
2. Specifically,it is no longer sufficient to assume that C is merely a bounded op-erator. One also cannot simply assume C to be of Hilbert-Schmidtclass, as this rules out the important case of C = Id , which isbounded, but not Hilbert-Schmidt.As a final remark, note that an analogous proof strategy can beapplied to the Schatten p -norm, which can be seen as the l p gen-eralization of the Hilbert-Schmidt norm. There, the k th term of theseries would be 1 / | λ γ k − µ | p and convergence would be guaranteedfor γ > / p . Thus, one can find p large enough so that R µ ( ∆ γ ) hasa well-defined Schatten p -norm for any given γ > Appendix C:
Comparison to heat maskBesides our resolvent-based commutativity term, another natu-ral bounded option would be to use the commutativity with the heat operators, which are also bounded linear functional operators,which can act as regularizers on functional maps. This would lead
Figure 21:
Visualization of the heat mask with different time-scaleparameter T . − − − . . . . T A v e r a g ee rr o r heatstandardours Figure 22:
Changing the time-scale parameter T of the heat mask.Here we test the performance of the heat mask (yellow curve) withdifferent T and report the average direct error of 50 random FAUSTisometric shape pairs. The average error of the standard mask(blue dashed line) and the average error of our resolvent mask (reddashed line) are included for comparison. to the following mask matrix: M heat ( i , j ) = exp ( − T Λ ( i )) − exp ( − T Λ ( j )) , (25)where T is the scalar time parameter. Fig. 21 visualizes the heatmask with different values of T on an isometric shape pair. InFig. 22 we compared the performance of the heat mask with dif-ferent values of T ranging from 10 − to 500 on FAUST isometricshape pairs. We can see that, as a bounded operator, the heat maskhas a better performance than the standard mask. When T =
5, itgives the best performance among the tested values. The heat maskwith T = Appendix D:
Change the descriptor sizeFig. 23 shows the average direct error of 300 FAUST shapes withdifferent number of input descriptors. The descriptors are used tooptimize a 100-by-100 functional map for each test pair. The resultswith ICP/BCICP refinement are also included.We also compare to the results obtained directly using the codeand the dataset provided by the authors of [RPWO18] (see Table 2).The reproduced results that are reported in Table 2 are consistentwith the values reported in the Table 1-2 ("WKS + directOp +BCICP") in the paper, and Table 2-3 ("WKS + directOp") in the c (cid:13) (cid:13) . Ren, M. Panine, P. Wonka, & M. Ovsjanikov / Structured Regularization of Functional Map Computations A v e r a g e d i r ec t e rr o r standard (1)slanted (2)ours (3)(1) + ICP(2) + ICP(3) + ICP(1) + BCICP(2) + BCICP(3) + BCICP Figure 23:
The average geodesic error v.s. the number of descrip-tors used for functional map estimation on the FAUST dataset. Notethat our mask leads to better initialization, which results in bettermaps, even after the ICP and BCICP refinement.
Table 2:
We also compare our results to the exact settingof [RPWO18] on the same set of shape pairs of FAUST dataset. Theaverage direct geodesic error over 300 shape pairs are reportedand our mask leads to better results.
Avg. direct error ( × − ) Ini +ICP +BCICP[RPWO18] 175.5 121.2 58.2Ours 72.7 53.5 42.0supplementary materials of [RPWO18]. Note that [RPWO18] splitthe dataset into isometric and non-isometric categories, and herewe report the results altogether.Recall that the total energy to optimize is E (cid:0) C (cid:1) = α E desc + α E mult + α E orient + α E mask . We would like to highlight the factthat the parameters in [RPWO18] are not set in the same way asin our work. In particular, in [RPWO18] the standard Laplacianmask is used, and the weight α i are set to fixed values α ∗ i . While inour test, we used the proposed resolvent Laplacian mask, and theweight α i are set to α ∗ i / E i ( C ini ) , where E i ( C ini ) is the correspond-ing energy term acting on the initial functional map C ini .Our approach allows a better control over the relative contribu-tion of the different terms in the energy and we observed that ittypically works better in practice as well. Specifically, in our set-ting, the relative weight of each term α i E i ( C ) are fixed across dif-ferent shape pairs. In this case, if we change the mask construc-tion, we can conclude that the improvement indeed come from ourproposed resolvent mask. However, in the comparison to the ex-act setting of [RPWO18], since only α ∗ i is fixed over different testpairs, the relative weight of each term α ∗ i E i ( C ) can have differentscale since different test pairs may have different scale of eigenval-ues, descriptors and etc. Therefore, the improvement from Table 2is not obtained in a well controlled setting, since the improvementcan also come from the change of the relative weight of differentenergy terms as well as from our resolvent mask.Thus, in addition to the new Table 2, we also emphasize thatin Fig. 23 we provide a more fair and controlled comparisonto [RPWO18] in which 10 pairs of descriptors are used. Note that,in all the tests across the paper, we used the above discussed way to fix the relative weight of the mask term w.r.t. the rest terms to makesure the improvement indeed comes from our new resolvent mask. Appendix E:
Stability under remeshing and refinement n = 6890 n = 200 n=300 n=500 n=1000 n=3000 n=5000 n=6890SourceTarget S t a nd a r d O u r s Figure 24:
First row: The target shape is fixed, while the sourceshape is remeshed and downsampled to resolution ranging from200 to 5000 (the original target shape without remeshing is shownin the last column); Second row: we use the standard mask to opti-mize a 100-by-100 functional map with 3 descriptors. The recov-ered pointwise maps are visualized on the corresponding targetshape; Third row: similar to the second row but using our resol-vent mask.
Fig. 24 illustrate the stability of our resolvent mask underremeshing and refinement. Specifically, we fix the source shape,and remesh and downsample the target shape to different resolu-tions ranging from 200 to 5000. Note that the original source shape(first row, first column) and the original target shape (first row, lastcolumn) have the same triangulation. The downsampled meshes(using QSlim) are shown in the first row with the number of verticesreported in the above.We then compute a 100-by-100 functional map between thesource shape and the downsampled target shape using the stan-dard Laplacian mask (the corresponding point maps are shown inthe second row) and our resolvent mask (in the third row). We cansee that, our resolvent mask is more stable than the standard maskacross different mesh resolution and irregular/inconsistent triangu-lation. c (cid:13) (cid:13)(cid:13)