Je remercie chaleureusement toutes les personnes qui m'ont aidé durant mes travaux de recherche, notamment mes directeurs de thèse Maks Ovsjanikov et Antonin Chambolle, pour leurs intérêts et leurs soutiens durant ces trois années. J'aimerais remercier mon jury : Quentin Mérigot et Michael Bronstein qui ont rapporté cette thèse, Leonidas Guibas, Tamy Boubekeur et Boris Thibert qui ont accepté de faire partie du jury. J'aimerais aussi remercier la DGA pour avoir financé ces trois années de doctorat. Merci aux chercheurs rencontrés au cours de cette thèse et avec qui j'ai eu le plaisir collaborer : Justin Solomon, Omri Azencot et Mirela Ben-Chen.

Je tiens à remercier le personnel administratif du CMAP en particulier Nasséra Naar et Alexandra Noiret ainsi que Evelyne Rayssac au LIX sans qui cette soutenance n'aurait pas eu lieu.

Durant cette période sur plateau de Saclay il m'a été permis d'effectuer de nombreuses rencontres au CMAP, au LIX ainsi que parmi mes voisins de l'équipe Geometrica (ou DataShape ou tagada je ne sais plus). Je remercie donc chaleureusement toutes les personnes qui ont contribué à la bonne humeur et à l'atmosphère sympathique (mais studieuse) de ces laboratoires. Il serait malheureusement impossible de nommer toutes ces personnes ici sans faire d'impardonnables oublis. Fort de ce constat, j'espère que chacun pourra se reconnaitre dans les lignes qui suivent.

Merci à vous polytechniciens, polytechniciennes, centraliens, centraliennes, normaliens, normaliennes, ENStaiens, physiciens, physiciennes, statisticiens, statisticiennes, mathémagiciens, mathémagiciennes, adorateur de la mécanique kantique, propagandiste du transport optimal, annonciateur des méthodes de Stein, admirateur des distances de Wasserstein, renonciatrice des problèmes inverses, ordonnateur du permutahèdre, approximateur d'espace tangents, explicateur de la persistance homologique, promoteur en TDA, scruteur de graphs de Reeb, théoricien de l'action utilitariste, mapper-évangéliste, polygon mesh processeur, géomètre différentiel, flotteur par courbure moyenne, profiteur de quiver, zigzagueur en persistance, distance à la mesuriste, optimiste de forme, projeteur en simplex, peuple premier du bureau 2033, habitants de la deuxième ère du 2033, représentante de la minorité génoise, migrante portoricaine, festivaliers du film indépendant de Porquerolles, grand commandeur de la guerre des poireaux, décorateur en déchets animaliers, dompteur de moustiques, survivaliste primitiviste, lecteur de la Jaune et la Rouge, fan de Stupeflip, dessinateur de lapin compulsif, breton indépendantiste, touristes venus d'Ulm, grimpeur-bloqueur, gréviste CNTiste, militant de l'Humanuté.

of operators, is mapped to a low dimensional space with a simple Principal Component Analysis. Yet the embedding accurately depicts the degrees of freedom of the deformation. This remark along with interesting algebraic properties lead to the intuition that operators are more convenient to describe deformations than Euclidean embeddings.

In this thesis we attempt to use the shape difference operators for deformation synthesis and for exploring the space of deformable surfaces. This only prior work by Boscaini and al. in [15] reconstructs an embedded surface from an operator representing the intrinsic structure of the shape. However, this approach is limited as it does not take into accounts extrinsic information. We explore different ways of encoding extrinsic curvature information as operators with theoretical guarantees of completeness.

Chapter 1

Introduction

Operator Representation in Geometry Processing

A challenging issue in geometry processing is the wide range of data structures. Discrete geometry is encoded in various ways mostly depending on the acquisition procedure and the application domain. Raw data from a scanner are stored as noisy point-clouds; physically-based simulation is typically done on simplicial complexes; animators and architectural geometry practitioners traditionally work with quad meshes. However, standard problems in geometry processing, such as shape matching, deformation transfer and deformation analysis, should be solved regardless of the representation. Most of the time each community provides their own solutions with respect to their native data which often do not allow cross-data analysis. This leads to the necessity of defining a framework agnostic to the underlying structure.

The basic idea of the functional representation framework, first introduced in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF], is to reach for the most common denominator of geometry analysis, namely real-valued functions. Functions have been historically used and well-studied in many contexts providing a large literature to rely on. For example differential operators have been defined on triangle and tetrahedral meshes (Finite Element Method [START_REF] Brenner | The mathematical theory of finite element methods[END_REF], Discrete Exterior Calculus [START_REF] Hirani | Discrete exterior calculus[END_REF]), polygonal surfaces [START_REF] Alexa | Discrete laplacians on general polygonal meshes[END_REF] and point-clouds (meshfree methods based on Radial Basis Functions or Moving Least Squares [START_REF] Gross | Point-based graphics[END_REF]) to name just a few. Moreover many of these methods come with convergence and stability analysis. Consequently functional operators can be implemented in many situations while studying a single continuous theory. As an example in [START_REF] Crane | Geodesics in Heat: A New Approach to Computing Distance Based on Heat Flow[END_REF] the authors propose an algorithm to compute geodesic distances based on the Laplace-Beltrami operator allowing them to exhibit results on radically different data structures. This flexibility is also one key of the success of spectral methods and global descriptors such as Global Point Signature [START_REF] Rustamov | Laplace-beltrami eigenfunctions for deformation invariant shape representation[END_REF], Heat Kernel Signature [START_REF] Sun | A concise and provably informative multi-scale signature based on heat diffusion[END_REF] and Wave Kernel Signature [START_REF] Aubry | The wave kernel signature: A quantum mechanical approach to shape analysis[END_REF].

The baseline idea developed in this thesis is to represent standard tools for deformable surfaces (diffeomorphisms, intrinsic distortions, tangent and extrinsic vector fields) as operators acting on functions, following the intuition that their discretizations have been well-studied and are easily comparable across data structures. The starting point of the thesis is the functional map representation for diffeomorphisms, and shape differences, representing intrinsic distortions. We propose in Chapter 5 a novel method, robust to noise and large deformations, to solve shape matching problems with functional maps. Chapters 6, 7 and 8 are providing different answers to an open problem in this framework: the conversion of operator-based representation to point-based representation. Namely, given a functional map we propose an algorithm to convert it to a continuous vertex-to-vertex assignment or given shape difference operators we suggest solutions to find the underlying deformed surfaces.

Functional Maps

The space of diffeomorphisms is challenging to handle in both continuous and discrete settings. In the simple case of two triangle meshes, the simplest representation consists of assigning a vertex on the first shape to a vertex on the second shape. However when dealing with different number of vertices this makes, for example, the computation of its inverse difficult. A more complete description is therefore needed. The most straightforward answer is to consider a piecewise linear map with the possibility of mapping vertices inside triangles. This way it is possible to compute the map and its inverse consistently. Such a description of discrete diffeomorphism has been used in practice [START_REF] Aigerman | Seamless surface mappings[END_REF][START_REF] Aigerman | Orbifold tutte embeddings[END_REF] but suffers some drawbacks as it cannot be stored compactly and often leads to challenging optimization problems.

In this context representing the mapping ϕ as the linear composition operators f → f • ϕ provides a flexible alternative globally agnostic to change of connectivity and easily defined on various types of data. Functional maps are a simplified and compact representation of correspondences as they are often computed between reduced function bases. They can be stored as matrices and analyzed using standard linear algebra tools. Follow-up works (e.g. [START_REF] Rodolà | Dense non-rigid shape correspondence using random forests[END_REF][START_REF] Rodolà | Partial functional correspondence[END_REF][START_REF] Kovnatsky | Functional correspondence by matrix completion[END_REF][START_REF] Pokrass | Sparse modeling of intrinsic correspondences[END_REF]) focus on computations of functional maps between shapes with unknown correspondences by minimizing various kinds of distortion with guiding functional constraints. However, functional constraints are often unreliable in case of large or non-isometric deformation and do not fully describe the mapping. We propose a supervised learning method to identify the most trustworthy set of constraints and the part of the surface on which the functional map is accurate. Another challenge of the operator representation is the problem of obtaining a point-to-point map from a functional map. The current methods [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF][START_REF] Rodola | Point-wise map recovery and refinement from functional correspondence[END_REF] do not ensure that resulting map is continuous and may lead to flaws in the map. To overcome this issue we restrict the conversion to continuous maps by using an adequate representation of vector fields as operators.

Shape Differences

Representing metric distortion can be challenging when dealing with heterogeneous data. The shape difference operators [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] take advantage of the functional map representation to provide a coordinate free encoding of the metric distortion up to isometric deformation. Those operators act linearly on functions and separate local area change and angle modification, even if those changes are noisy and not easy to define directly across different data structures, function correspondences induce a notion of smoothness up to a certain level-of-details.

Shape difference operators were originally presented as analysis tools since they can be used for comparison and composition of deformations across shapes. Interestingly the space of shape differences seems "flat" and well behaved. For instance in Figures 7.7 and 7.9 (see also Figure 3 in [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF]) a collection of surfaces, represented as a collection remain: the existence of internal symmetries and the absence of global constraints linking nearby points lead to discontinuities in the resulting map. In Chapter 6 we provide a method to restrict the conversion process to continuous maps. While this problem is challenging with standard analysis tools, it becomes much more tractable with the operator representation. The outcome is a provably smooth conversion of a functional map using vector field flow.

Overall, the main contribution of this first part is a complete pipeline for shape matching problems in the functional map framework: from the computation of reliable maps to a final point-to-point correspondence.

While the first part is aimed at finding similarities across near-isometric shapes, in the second part, our goal is to quantify and manipulate differences in the form of non-isometric distortion. In particular we are interested in a complete representation of surfaces as operators.

Deformation to Shape

Shape differences enjoy many interesting algebraic properties. In particular they allow the transfer, composition and addition (in special cases) of deformations. Thus, exploring and creating new intrinsic deformation is relatively easy in this representation as it only requires matrix manipulations. The only missing step to produce embedded surfaces is a conversion from the functional characterization to a surface embedding. Namely, given a base shape and shape differences, find the deformed surface. This problem is closely related to recovering the shape from its Laplacian operator which is known to be possible theoretically but for which very few computational methods are available. This thesis attempts to fill this hole with two different approaches.

The first development in Chapter 7 starts with the remark that the discrete metric on triangle meshes (edge length) can be recovered from shapes differences by solving a set of linear equations. It follows that in the case of noisy and limited information the metric can be recovered through two convex optimization problems. However, intrinsic information alone is not sufficient to recover an embedded surface as multiple solutions exist. We investigate the possibility of using offset surfaces to encode curvatures. The outcome is a complete characterization of triangle meshes through functional operators allowing us to deform surfaces according to shape difference operators and a corresponding continuous theory.

The second development in Chapter 8 introduces a functional characterization of extrinsic vector fields. To do so we introduce a unique shape difference accountable for all intrinsic changes. Then we consider the operator associated to an infinitesimal displacement. This way deformation fields are characterized by the distortion they induce on the metric. Interestingly they enjoy similar properties as shape differences (composition, transfer) enabling significantly simpler deformation reconstruction compared to the previous method. We provide theoretical proofs of informativeness in both continuous and discrete cases.

This thesis constitues a first step toward using functional representation for characterization of deformable surfaces and deformation synthesis.

Chapter 2

Introduction en français

L'une des difficultés rencontrées dans le traitement de forme 3D est la large gamme de structures de données. Les formes 3D et autres géométries discrètes sont codées de diverses manières dépendant principalement de la méthode d'acquisition et du domaine d'application. Les données brutes provenant d'un scanner sont stockées sous forme de nuages de points souvent bruités ; les simulations physiques s'effectuent habituellement sur des complexes simpliciaux ; les modeleurs 3D ainsi que la communauté de la géométrie architecturale travaillent traditionnellement avec des maillages quadrangulaires. Cependant, de nombreux problèmes standards liés au traitement de la géométrie, tels que la mise en correspondance de formes, le transfert de déformations et l'analyse de déformations, doivent être résolus indépendamment de la représentation utilisée. La plupart du temps chaque communauté fournit ses propres solutions vis-à-vis de ses données de base qui, souvent, ne sont pas adaptables aux autres types de représentation. Cela nous conduit à la nécessité de définir un cadre d'analyse indépendant des structures sous-jacentes.

L'idée à l'origine de la représentation fonctionnelle, d'abord introduite dans [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF], est de revenir au dénominateur commun dans l'analyse de la géométrie, à savoir les fonctions à valeurs réelles. Les fonctions définies sur la géométrie ont été historiquement utilisées et étudiées de manière approfondie dans de nombreux contextes fournissant une abondante littérature sur laquelle se reposer. Par exemple, les opérateurs différentiels ont été définis sur des maillages triangulaires et tétraédriques (Méthode des Eléments Finis [START_REF] Brenner | The mathematical theory of finite element methods[END_REF], Discrete Exterior Calculus [START_REF] Hirani | Discrete exterior calculus[END_REF]), sur des maillages polygonaux [START_REF] Alexa | Discrete laplacians on general polygonal meshes[END_REF] et des nuages de points (méthodes "sans maillage" basées sur des fonctions de base radiale ou Moving Least Squares [START_REF] Gross | Point-based graphics[END_REF]) pour n'en citer que quelques-uns. De plus, beaucoup de ces méthodes vont de paire avec une analyse de leur convergence et de leur stabilité. Par conséquent les opérateurs fonctionnels peuvent être mis en oeuvre dans de nombreuses situations tout en se référant à une théorie unique. A titre d'exemple, dans [START_REF] Crane | Geodesics in Heat: A New Approach to Computing Distance Based on Heat Flow[END_REF] les auteurs proposent un algorithme pour calculer les distances géodésiques grâce à l'opérateur de Laplace-Beltrami. Ils obtiennent ainsi des résultats sur des structures de données radicalement différentes. Cette flexibilité est également l'une des clés du succès des méthodes spectrales et des descripteurs globaux tels que Global Point Signature [START_REF] Rustamov | Laplace-beltrami eigenfunctions for deformation invariant shape representation[END_REF], Heat Kernel Signature [START_REF] Sun | A concise and provably informative multi-scale signature based on heat diffusion[END_REF] et Wave Kernel Signature [START_REF] Aubry | The wave kernel signature: A quantum mechanical approach to shape analysis[END_REF].

L'idée de base développée dans cette thèse est de représenter des outils standards pour les surfaces déformables (difféomorphismes, distorsions intrinsèques, champs de vecteurs tangents et extrinsèques) par des opérateurs agissant sur des fonctions, suivant l'intuition que leurs discrétisations ont été préalablement étudiées et sont facilement comparables quelles que soient les structures de données considérées. Le point de départ de cette thèse est la représentation fonctionnelle des difféomorphismes appelée Application Fonctionnelle (ou "Functional Map" en anglais), et des distorsions intrinsèques nommées Opérateurs de Déformation (ou "Shape Differences" en anglais). Nous présentons dans le Chapitre 5 une nouvelle méthode, robuste aux bruits et aux grandes déformations, pour trouver des correspondances entre surfaces en utilisant une approche fonctionnelle. Les Chapitres 6, 7 et 8 fournissent des réponses différentes à un problème ouvert dans le cadre fonctionnel : la conversion de la représentation à base d'opérateurs vers une représentation ponctuelle. Plus précisément, étant donné une application fonctionnelle, nous proposons un algorithme pour la convertir en une affectation continue de sommet à sommet ou à partir d'un opérateur de différence de forme, nous proposons des solutions pour retrouver la surface déformée sous-jacente.

Applications fonctionnelles

L'espace des difféomorphismes est difficile à appréhender que ce soit dans le cadre continu ou discret. Considérons deux maillages triangulaires, la représentation la plus simple consisterait à attribuer à chaque sommet de la première forme un sommet de la seconde. Toutefois, lorsque l'on considère deux maillages avec un nombre différent de sommets le calcul de l'application inverse devient difficile. Une description plus complète est donc nécessaire. Une possibilité serait de considérer une application linéaire par morceaux permettant d'assigner des sommets à l'intérieur des triangles. De cette façon, il est possible de calculer le difféomorphisme et son inverse de façon cohérente. Une telle description discrète a été mise en oeuvre notamment dans [START_REF] Aigerman | Seamless surface mappings[END_REF][START_REF] Aigerman | Orbifold tutte embeddings[END_REF] mais possède quelques inconvénients car elle ne peut-être stockée de manière compacte et conduit souvent à des problèmes d'optimisations non-convexes difficiles à appréhender.

Dans ce contexte, représenter l'application ϕ qui lie deux formes géométriques par l'opérateur de composition f → f • ϕ fournit une alternative globalement indifférente au changement de connectivité et facilement définie sur différents types de données. Les applications fonctionnelles sont des représentations simplifiées et compactes des correspondances car elles sont souvent calculées entre des bases de fonctions de tailles réduites. Elles peuvent être stockées sous forme de matrices et analysées à l'aide d'outils d'algèbre linéaire standards. Les travaux qui ont suivis (par exemple [START_REF] Rodolà | Dense non-rigid shape correspondence using random forests[END_REF][START_REF] Rodolà | Partial functional correspondence[END_REF][START_REF] Kovnatsky | Functional correspondence by matrix completion[END_REF][START_REF] Pokrass | Sparse modeling of intrinsic correspondences[END_REF]) se concentrent sur des calculs d'applications fonctionnelles entre surfaces dont la mise en correspondance est inconnue. Pour ce faire ils minimisent différents types de distorsions intrinsèques auxquelles s'ajoutent des contraintes fonctionnelles. Cependant, ces contraintes fonctionnelles sont souvent peu fiables en cas de déformations importantes ou non-isométriques et ne décrivent les correspondances entre les surfaces que de façon incomplète. Nous proposons donc une méthode d'apprentissage supervisée pour identifier l'ensemble des contraintes le plus digne de confiance et la partie de la surface sur laquelle l'application fonctionnelle est la plus précise. Un autre défi de cette représentation fonctionnelle est le problème de l'obtention d'une application point-à-point à partir d'un opérateur. Les méthodes actuelles [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF][START_REF] Rodola | Point-wise map recovery and refinement from functional correspondence[END_REF] ne garantissent pas que l'application résultante soit continue et peuvent conduire à des défauts de continuité. Pour surmonter ces problèmes, nous limitons la conversion à des applications continues en utilisant une représentation adéquate des champs de vecteurs par des opérateurs.

Opérateurs de déformation

Représenter les distorsions de la métrique peut être difficile lorsqu'il s'agit de traiter des données de provenances variées. Les opérateurs de déformation [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] tirent parti de la représentation fonctionnelle pour fournir un encodage des distorsions intrinsèques faisant abstraction du plongement de la surface dans l'espace ambiant. Ces opérateurs agissent linéairement sur les fonctions et font la distinction entre les modifications de volumes et l'évolution des angles. Ces changements sont locaux, facilement perturbables et difficiles à définir correctement selon les différentes structures de données. Toutefois les opérateurs de déformation parviennent à les rendre stables et facilement analysables.

Les opérateurs de déformation ont été initialement présentés comme des outils d'analyse, car ils peuvent être utilisés pour comparer et composer des déformations. Il est intéressant de remarquer que l'espace généré par ces opérateurs semble "plat" et stable. Par exemple les Figures 7.7 et 7.9 (voir aussi Figure 3 dans [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF]) représentent une collection de surfaces, interprétée comme une collection d'opérateurs, est envoyée dans un espace euclidien de petite dimension grâce à une simple analyse en composantes principales. Ce plongement décrit précisément les degrés de liberté des déformations. Cette remarque ainsi que des propriétés algébriques intéressantes conduisent à l'intuition que les opérateurs sont plus commodes pour décrire des déformations intrinsèques que l'analyse directe d'un plongement dans l'espace euclidien.

Dans cette thèse, nous tentons d'utiliser les opérateurs de déformation pour synthétiser et explorer l'espace des surfaces déformables. L'unique travail préalable par Boscaini et al. dans [START_REF] Boscaini | Shape-fromoperator: Recovering shapes from intrinsic operators[END_REF] retrouve le plongement d'une surface à partir d'un opérateur représentant la structure intrinsèque de la forme. Cependant, cette approche est limitée car elle ne prend pas en compte les informations extrinsèques. Nous explorons différentes façons de stocker les informations données par la courbure en utilisant des opérateurs fonctionnels.

Résumé de la thèse

Cette thèse porte sur le traitement de surfaces déformables dans le cadre d'une représentation fonctionnelle. L'étude est divisée en deux parties qui coïncident avec l'analyse de deux types d'opérateurs : les applications fonctionnelles et les opérateurs de déformation.

Des surfaces aux déformations

Le but de cette partie est de trouver des correspondances entre deux formes données. Nous nous intéressons en particulier aux informations délivrées par la structure intrinsèque de ces formes. Pour rendre le problème plus abordable, nous nous plaçons dans un modèle restrictif en ne considérant que des déformations quasi-isométriques.

Trouver des correspondances entre formes 3D pour des déformations quasi-isométriques fut présentée comme la première application de la représentation fonctionnelle par [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] et est toujours la direction de recherche privilégiée dans ce domaine. L'une des raisons de ce succès étant que ce problème, difficile de prime abord, devient un système des moindres carrés dans cette nouvelle représentation. La formulation la plus basique prend en compte la distance à l'isométrie et des correspondances fonctionnelles, souvent obtenues à partir de descripteurs globaux (HKS, WKS) restant stables même après de légères déformations. Cependant, l'application résultante peut avoir des artefacts dus à de mauvaises correspondances fonctionnelles, à des perturbations ou à une interpolation spectrale trompeuse. De nombreuses variantes ont été proposées pour améliorer la qualité des résultats [START_REF] Rodolà | Dense non-rigid shape correspondence using random forests[END_REF][START_REF] Kovnatsky | Functional correspondence by matrix completion[END_REF][START_REF] Pokrass | Sparse modeling of intrinsic correspondences[END_REF]) afin de permettre des correspondances partielles [START_REF] Rodolà | Partial functional correspondence[END_REF] ou de tirer parti de cohérences cycliques dans des réseaux de difféomorphismes [START_REF] Huang | Functional map networks for analyzing and exploring large shape collections[END_REF].

Dans le Chapitre 5 nous prenons une approche légèrement différente : au lieu de modifier la formulation initiale, nous tentons d'obtenir les résultats optimaux à partir des données initiales. Les contraintes fonctionnelles, potentiellement peu fiables, nous obligent à trier et à pondérer ces correspondances. Par ailleurs, l'application fonctionnelle résultante n'est fiable que sur un sous-espace de fonctions qu'il nous appartient d'identifier. Notre principale contribution est un algorithme d'apprentissage supervisé capable de calculer un ensemble de pondérations évaluant conjointement l'utilité de chaque descripteur et capable de déterminer les fonctions les plus fidèlement transférées par l'application.

Une fois calculées les applications fonctionnelles offrent un autre défi qui consiste à les convertir en une représentation ponctuelle, utilisée dans la plupart des applications. L'article original propose une méthode qui peut être comprise comme l'algorithme Iterative Closest Point appliqué au domaine spectral. Bien qu'il ait été amélioré par des articles postérieurs [START_REF] Rodola | Point-wise map recovery and refinement from functional correspondence[END_REF], un problème majeur demeure : l'existence de symétries internes et l'absence de contraintes globales reliant les points proches les uns des autres conduisent à des discontinuités dans l'application résultante. Dans le Chapitre 6 nous fournissons une méthode pour restreindre ce processus de conversion aux seules applications continues. Ce problème, difficile en utilisant des outils d'analyse standards, devient beaucoup plus abordable en utilisant une représentation fonctionnelle. Ce problème, difficile en utilisant des outils d'analyse standards, devient beaucoup plus abordable en utilisant une représentation fonctionnelle. Le résultat est la conversion d'une application fonctionnelle vers une application continue en utilisant le flot d'un champ de vecteurs.

Ainsi, la contribution principale de cette première partie est une procédure complète pour résoudre les problèmes de mise en correspondance de formes dans le cadre d'une représentation fonctionnelle : du calcul d'une application fonctionnelle fiable à une correspondance ponctuelle entre surfaces.

Alors que la première partie vise à trouver des similarités entre des formes quasiisométriques, dans la deuxième partie notre objectif est de quantifier et de manipuler des distorsions non-isométriques de la métrique. En particulier, nous nous intéressons à une représentation complète des surfaces par des opérateurs fonctionnels.

Synthèse de formes et représentation fonctionnelle

Les opérateurs de déformation jouissent de nombreuses propriétés algébriques intéressantes. En particulier, ils permettent le transfert, la composition et l'addition (dans des cas particuliers) de déformations. Ainsi, l'exploration et la création de nouvelles déformations intrinsèques sont aisées dans cette représentation car elles ne nécessitent que des manipulations de matrices. L'unique étape manquante est la conversion d'une caractérisation fonctionnelle des déformations en un plongement dans l'espace euclidien. Plus précisément, étant donné une forme de base et un opérateur de déformations, il nous faut retrouver la surface déformée. Ce problème est étroitement lié à la reconstruction d'une forme à partir de son laplacien, connu pour être théoriquement possible bien que peu d'algorithmes soient actuellement disponibles pour le résoudre. Cette thèse tente de combler ce vide par deux approches différentes.

Le premier développement, présent dans le Chapitre 7, commence avec la remarque que la métrique discrète d'un maillage triangulaire (longueur des arêtes) peut être récupérée à partir des opérateurs de déformations, et ce, en résolvant un système d'équations linéaires. Il en résulte que dans le cas d'informations bruitées et limitées la métrique peut être récupérée grâce à deux problèmes d'optimisations convexes. Cependant, l'information intrinsèque n'est pas suffisante pour obtenir le plongement d'une surface car il existe plusieurs solutions à ce problème. Dans ce Chapitre nous étudions donc la possibilité d'utiliser des surfaces parallèles pour encoder la courbure. Il en ressort une caractérisation complète des maillages triangulaires par des opérateurs fonctionnels ainsi que des garanties théoriques de reconstruction.

Le deuxième développement inclut dans le Chapitre 8 introduit une caractérisation fonctionnelle des champs de vecteurs extrinsèques. Pour cela, nous introduisons un opérateur de déformation unique, responsable de tous les changements de la métrique. Ensuite, nous considérons l'opérateur associé à un déplacement infinitésimal. De cette façon, les champs de déformation sont caractérisés par la distorsion qu'ils induisent sur la métrique. Ils jouissent de propriétés similaires aux opérateurs de déformations (composition, transfert), tout en permettant une reconstruction des déformations nettement plus simple que la méthode précédente.

Cette thèse constitue une première étape vers l'utilisation de la représentation fonctionnelle pour la caractérisation des surfaces déformables et la synthèse de déformation.

Publications liées à cette thèse

Cette thèse est en partie basée sur les publications suivantes :

• E. Corman, M. Ovsjanikov and A. Chambolle, Supervised descriptor learning for non-rigid shape matching, in ECCV 2014 Workshops, Part IV, Springer International Publishing, 2014. [START_REF] Corman | Supervised descriptor learning for non-rigid shape matching[END_REF] • E. Corman

Context in Differential Geometry and its Discrete Equivalents

This chapter introduces the main mathematical tools that will be used throughout the thesis. A brief overview of differential geometry focused on intrinsic geometry is provided. We put emphasis on the role of the metric tensor and the pullback metric in the definitions of differential operators. This chapter also provides a brief remainder of the Finite Element Method applied to the discretization of common differential operators. This chapter can be skipped by readers familiar with the continuous theory.

Differential Geometry

The study of surface deformations is central in this thesis. Therefore some notions and terminology from differential geometry are recalled here. The aim of this section is to introduce the notions of metric tensor and more importantly of pullback metric by a diffeomorphism. The pullback metric accounts for the intrinsic distortion underlying a mapping between shapes. This will be a key tool for our analysis of surface deformations as it allows a classification of mappings according to their impact on the local structure.

To represent diffeomorphisms and distortions as operators on functions, we will need some basic differential operators such as gradients and Laplacians. We will focus on the change of variables which is equivalent to a change of metric in differential operators. This remark allows us to construct operators provably informative of the intrinsic structure and explain the leading role of the Laplace-Beltrami operator in geometry processing.

This basic introduction to differential geometry is greatly inspired by [START_REF] Morita | Geometry of differential forms[END_REF][START_REF] Rosenberg | The Laplacian on a Riemannian Manifold[END_REF]. The Einstein summation convention is sometimes used to lighten the notations. Repeated indices imply a summation, i.e. g ij f j means j g ij f j . Moreover when dealing with tensors upper indices denote an element of the inverse matrix A ij = (A -1 ) ij .

Manifold

In practice we will most of the time deal with surfaces in R 3 . However, operations on intrinsic informations (e.g. transport, composition of distortion) may lead to metric tensors which do not correspond to embedded surfaces and are only meaningful for general manifolds. Thus, it is important to start with abstract definitions.

A topological manifold M is a space such that each point has a neighborhood diffeomorphic to R n . To introduce differential calculus, it is necessary to always refer to the Euclidean space where derivatives are well-defined. To do so we introduce a local parametrization of M called an atlas. An atlas on M is a countable collection of pairs (U α , ψ α ) named charts where U α are sets covering M and ψ α : U α → R n are homeomorphisms acting as local parametrizations.

Since the sets U α cover the entire manifold some points maybe found in two different charts. If U α and U β overlap then one can change the coordinate system by means of the transition map:

ψ βα = ψ β • ψ -1 α : ψ α (U α ∩ U β ) → ψ β (U α ∩ U β ).
Note that ψ βα maps an open set of R n to an open set of R n , so at a point p in U α ∩ U β one can decompose this mapping in terms of coordinates (x 1 (p), . . . , x n (p)) with respect to U α and (y 1 (p), . . . , y n (p)) with respect to U β . Two systems of coordinates are linked by:

y i (p) = ψ βα (x 1 (p), . . . , x n (p)) i .
A transition map and two coordinate systems are represented in Figure 3.1. This formulation is convenient as it allows us to use standard differential calculus in R n to study the properties of a curved surface. In particular the Jacobian matrix of ψ βα will appear in coordinate change for vector fields (Section 3.1.3).

The charts introduce a local structure on the manifold. However, to have a differential manifold we need to ensure that the local properties are consistent among themselves. Definition 3.1.1 A topological manifold M equipped with an atlas S = (U α , ψ α ) is called differential manifold if all its coordinate changes ψ βα are C ∞ maps.

U α U β M α β βα x 1 x 2 y 2 y 1

Tangent Vectors and Local Coordinates

Tangent vectors are a first step to introduce local intrinsic structure on a manifold. Let f : M → R be a real valued function on a differential manifold. If, for all coordinate systems (U, ψ) in an atlas,

f • ψ -1 : ψ(U ) ⊂ R n → R is a C ∞ function on ψ(U ) then f is by definition a C ∞ function on M .
Now let (U, ψ) be a coordinate system around a point p ∈ M and (x 1 , . . . , x n ) the coordinate functions with respect to U . The derivative of f in the direction x i denoted ∂ i f is naturally defined as:

∂ i f = ∂(f • ψ -1 ) ∂x i • ψ ∈ C ∞ (M ).
In case of a surface in R 3 , the gradient of a function is tangent to this surface. Moreover, the set of all possible directions for function derivatives defines the set of all tangent vectors. By analogy, the application V p : C ∞ (M ) → R acting on differentiable functions at point p ∈ M :

V p (f ) = n i=1 V i (p)∂ i f | p , (3.1)
is understood as a vector tangent to the manifold. The local coordinates of the vector V i (p) ∈ R depend on the chart (U, ψ) and correspond to the derivative in the canonical directions ∂ ∂x i p . We denote T p M , the set of all tangent vectors at point p, defined as the set of all linear combinations of the basis

∂ ∂x i p i
. The canonical basis satisfies the property ∂x j ∂x i = δ ij . A vector field on M is a smooth assignment of a tangent vectors V p ∈ T p M at each point p. Namely V is a C ∞ vector field if:

V (f ) ∈ C ∞ (M ), ∀f ∈ C ∞ (M ).
We denote the set of smooth vector fields V(M ).

Mappings

Along this thesis we will study relations between manifolds through mappings. Let M and N be two differential manifolds, ϕ : N → M is a map between them. Since the map is not real-valued we have to take into account the parametrizations from M and N to introduce a differential structure.

Definition 3.1.2 The map ϕ : N → M is a C ∞ map if for all charts (U, ψ) on N and (V, ϑ) on M , the functions ϑ • ϕ • ψ -1 : R n → R n are of class C ∞ (R n ).
Local coordinate systems (x 1 (p), . . . , x n (p)) around p ∈ N and (y 1 (ϕ(p)), . . . , y n (ϕ(p))) around q = ϕ(p) ∈ M are linked by the mapping ϕ:

y i (ϕ(p)) = ϑ • ϕ • ψ -1 (x 1 (p), . . . , x n (p)) i . (3.2)

Differential map

A C ∞ map not only describes pointwise correspondences but also map tangent vectors from one surface to another. The distortion of a tangent space is determined by the Jacobian of the mapping dϕ p .

Using the chain rules for derivatives, Equation (3.2) yields a change of coordinates formula or pushforward of tangent vector fields:

dϕ p ∂ ∂x i = n j=1 ∂y j ∂x i ∂ ∂y j , (3.3)
where dϕ denote the differential map or Jacobian of the transformation. Therefore, a map not only relates points but also tangent spaces. The differential map defines a linear mapping between the tangent space on N at point p and the tangent space on M at point ϕ(p):

dϕ p : T p N → T ϕ(p) M.

Categories of mappings

Mappings are categorized according to the properties of the differential map:

• Immersion: At each point p ∈ N the Jacobian is injective. It implies that the image of the manifold ϕ(N ) can intersect itself so it may not be a submanifold of M .

• Embedding: ϕ is immersion and also a homeomorphism from N onto ϕ(N ). In this case no self-intersection can occur and the image ϕ(N ) is a submanifold of M .

• Diffeomorphism: ϕ : N → M is a bijection and its inverse ϕ -1 : M → N is a C ∞ map. These requirements are stronger than for embeddings, in particular it implies that the Jacobian is bijective and that N, M have the same intrinsic dimension.

Riemannian Manifold

Metric tensor

The metric tensor extends the notion of scalar product in R n to tangent vectors. It is important to recall that general metric tensors can be arbitrarily chosen and are not necessarily induced by the embedding of a surface. The choice of the metric determines the specifying of the tangent plane, thus, specifying distances on the manifold and many differential operators.

Definition 3.1.3 A Riemannian metric is a family of bilinear forms defined at each point p of the differential manifold M by:

g p : T p M × T p M → R.
As Euclidean scalar products the metric should be symmetric, positive definite at each point. Moreover, a Riemannian metric varies smoothly on the manifold: for all differentiable vector fields V, U on M , p → g p (V p , U p ) is a smooth function from M to R.

Given a system of local coordinates (x 1 , . . . , x n ) around p, the set ∂ ∂x 1 , . . . , ∂ ∂xn is a basis of the tangent space T p M . So, the metric can be written as a matrix of R n×n :

g ij (p) := g p ∂ ∂x i p , ∂ ∂x j p .
Thus, the metric tensor can be interpreted as a Gramian matrix of tangent vectors. Geometrically the square root of the determinant is the volume of the parallelotope defined by the basis ∂ ∂x i . Therefore, the quantity det(g) can be understood as a local volume or area on the manifold.

Given two tangent vectors V, U ∈ T p M their scalar product is computed as a usual bilinear form within the coordinate system:

g p (V, U ) = ij V i g ij U j .

Surface in R 3

An important example that will be frequently used in this thesis is the case where M is an embedded surface of R 3 . In this configuration, the ambient Euclidean scalar product ., . naturally induces a metric on the surface. Given a chart (U, ψ) in the neighborhood of p ∈ M , the mapping

F = ψ -1 : ψ(U ) ⊂ R 2 → U ⊂ M ⊂ R 3 maps an open subset of R 2
to an open subset of R 3 . Therefore ∂F ∂x i (p) is a vector of R 3 tangent to the surface at p and can be identified with a basis of T p M . The induced metric in this basis is:

g ij (p) = ∂F ∂x i (p), ∂F ∂x j (p) . (3.4)
Many local coordinates expressions can be replaced by vectors in the ambient space. For instance the scalar product g(V, U ) can be replaced by V , Ū where V is the tangent vector V expressed in the global coordinate system. To avoid heavy notations we will consider the change of coordinates implicit when using the ambient scalar product.

Geodesic distance

The metric is a local quantity determining lengths of vectors. It induces a global structure on the manifold enabling us to measure distances between points. The geodesic distance is defined as the shortest path between two points on a manifold.

The application c : [0, 1] → M is a parametrized curve on M . Its velocity c ′ is a tangent vector whose length can be computed with the metric g. The arc length L(c) of the curve c is defined by similarity with the Euclidean case:

L(c) := 1 0 g c(t) (c ′ (t), c ′ (t))dt.

The pullback metric

The pullback metric is a key tool for analyzing metric distortion. Namely, given an immersion ϕ : N → M between the differential manifolds M and N , we are able to express g M , the metric tensor on M , as a distorted metric on N . This new tensor encodes the distortion of the tangent space by the underlying deformation ϕ.

Definition 3.1. [START_REF] Anguelov | Scape: shape completion and animation of people[END_REF] The pullback of a metric on M , denoted ϕ ⋆ g M , is defined as the action of the tensor on the pushforward tangent vectors:

(ϕ ⋆ g M ) p (V, U ) = g M ϕ(p) (dϕ(V ), dϕ(U )), ∀V, U ∈ T p N.
Given a coordinate system, Equation (3.3) leads to a local expression of the pullback metric:

(ϕ ⋆ g M ) ij (p) = ∂y k ∂x i g M kl (ϕ(p))
∂y l ∂x j .

The Jacobian of the immersion distorts the bilinear form in a symmetric way. If ϕ is an immersion, the tensor ϕ ⋆ g M is a metric on N .

Three types of mappings are defined according to the intrinsic properties they preserve (local areas, angles, geodesic distances). Figure 3.3 shows examples of deformation of the plane.

• Area preserving map. A diffeomorphism is said to be area preserving if local areas are preserved:

det(ϕ ⋆ g M ) p = det(g N ) p , ∀p ∈ N.
• Conformal map. A diffeomorphism ϕ : N → M is said to be conformal if there exists a positive function h ∈ C ∞ (M ), sometimes called the conformal factor, such that:

h(q)(ϕ ⋆ g M ) p = g N p ∀p ∈ N.
Such deformations preserve angles between curves but affect local areas. A classical example of such deformation are the Möbius transformations which span all the bijective conformal maps from the sphere to itself.

• Isometry. A conformal and area preserving diffeomorphism is an isometry. In this case, the pullback metric agrees with the metric on N :

(ϕ ⋆ g M ) p = g N p ∀p ∈ N.
A map is an isometry if and only if geodesic distances are preserved, namely d M (ϕ(p), ϕ(q)) = d N (p, q), as proven by the Myers-Steenrod theorem [START_REF] Myers | The group of isometries of a riemannian manifold[END_REF].

Integral on Manifold

Area preserving maps are naturally related to integral on manifold. In fact, a change of variables in integral is equivalent to using the measure induced by the pullback metric [START_REF] Morita | Geometry of differential forms[END_REF].

Volume form

First, we introduce the volume form dµ g defined locally by:

dµ g = det(g) dx 1 ∧• • •∧dx n .
The volume form induces a Borel measure µ g :

µ g (U ) = U dµ g , ∀U ⊂ M.
Note that µ will abusively stand for the local area det(g) instead of the measure.

Given a chart (U α , ψ α ) and a coordinate system (x 1 , . . . , x n ) the integral is then defined through its mapping to the Euclidean space:

Uα f (p) dµ g (p) = ψα(Uα) (f det(g)) • ψ -1 α (x 1 , . . . , x n ) dx 1 . . . dx n .

Pushforward measure

Let ϕ : N → M be a diffeomorphism between two manifolds and g N a metric on N . The pushforward measure ϕ ⋆ µ g N is a measure on M standing for the change of variables:

M f d(ϕ ⋆ µ g N ) := N f • ϕ dµ g N .
Applying the change of coordinate formula for volume form [START_REF] Morita | Geometry of differential forms[END_REF] leads to Proposition 3.1.5. Proposition 3.1.5 The volume form induced by the pushforward measure is the volume form corresponding to the pullback metric:

d(ϕ ⋆ µ g N ) = dµ (ϕ -1 ) ⋆ g N .
When ϕ is area preserving, integrals are left untouched by a change of variables:

M f dµ g M = N f • ϕ dµ g N , ∀f integrable.

Inner product

Let f, g : M → R be smooth functions on M . The space of square integrable functions with respect to the volume form dµ g is defined as:

L 2 (M, g) = f : M → R integrable : M f 2 dµ g < ∞ .
This space is equipped with the inner product f, g

(M,g) L 2 := M f g dµ g inducing the norm f (M,g) L 2 := f, f (M,g) L 2 1 2
. Mention of the metric and the manifold is sometimes omitted to avoid heavy notation

Gradient, Divergence and Laplacian

Gradient, divergence and Laplacian are three fundamental differential operators to study intrinsic geometry. Their definitions entirely rely on the underlying metric tensor. Moreover, as we will see later, they also fully specify the metric structure of the manifold. For this reason, the Laplacian is a key tool in geometry processing and computer graphics.

In particular it has been used for: fairing and compression (e.g. [START_REF] Karni | Spectral compression of mesh geometry[END_REF]), geodesic computation (e.g. [START_REF] Crane | Geodesics in Heat: A New Approach to Computing Distance Based on Heat Flow[END_REF]) and global descriptors (e.g. [START_REF] Reuter | Laplace-beltrami spectra as 'shape-dna'of surfaces and solids[END_REF][START_REF] Sun | A concise and provably informative multi-scale signature based on heat diffusion[END_REF][START_REF] Aubry | The wave kernel signature: A quantum mechanical approach to shape analysis[END_REF]) to name just few applications. Two properties explain this success: its invariance under isometric deformation and its eigenfunctions comparable to a Fourier basis on manifolds.

For a more comprehensive study of the Laplace-Beltrami operator in differential geometry we refer the interested reader to [START_REF] Rosenberg | The Laplacian on a Riemannian Manifold[END_REF].

Laplace-Beltrami operator

The directional derivative V (f ) is the action of a tangent vector on a function. It can also be thought of as the scalar product between a fixed vector and the function gradient. Thus, the gradient can be defined as the tangent vector satisfying:

V (f ) = g(V, ∇f ), ∀V ∈ V(M ).
Given a local coordinate system we can write

V i ∂ i f = V i g ij g jk ∂ k f
, so the components of the gradient are given by: (∇f

) i = g ij ∂ j f .
The gradient depends on the inverse metric tensor. As our main topic is intrinsic geometry, it is interesting to relate the pushforward gradient with the pullback metric. Lemma 3.1.6 Let ϕ : N → M be a diffeomorphism between two differential manifolds and g N a metric on N . The pushforward of the gradient is equal to the gradient defined by the pullback metric:

dϕ p ∇ g N (f • ϕ) = ∇ (ϕ -1 ) ⋆ g N f (ϕ(p)).
Proof Given local coordinates (x 1 , . . . , x n ) around p ∈ N and (y 1 , . . . , y n ) around q = ϕ(p) ∈ M , a direct computation yields:

dϕ p ∇ g N (f • ϕ) = dϕ p g ij N ∂ ∂x j (f • ϕ) ∂ ∂x i = g ij N ∂ ∂x j (f • ϕ) ∂y k ∂x i ∂ ∂y k = g ij N (ϕ -1 (q)) ∂y l ∂x j ∂f ∂y l (q) ∂y k ∂x i ∂ ∂y k = ((ϕ -1 ) ⋆ g N ) ij (q) ∂f ∂y j (q) ∂ ∂y i = ∇ (ϕ -1 ) ⋆ g N f (q).
This simple lemma expresses the idea that instead of using the differential map to transfer gradients of functions from N to M , one can transfer intrinsic information and then compute gradients locally with respect to this metric. Computing a the differential map may be challenging when considering meshes with different triangulation but transferring metrics, as we will see later, is more tractable since it is can be done by pullback of the Laplacian.

For manifolds without boundaries the Laplacian can be defined through the integration by parts identity:

M g(∇g, ∇f ) dµ g = - M g∆f dµ g . (3.5)
The Laplacian is therefore a symmetric and negative definite operator on smooth functions. The expression in local coordinates relies heavily on the metric:

∆f = 1 det(g) ∂ i det(g)g ij ∂ j f .

Eigen-decomposition

The Laplace-Beltrami operator has a countable set of eigenfunctions:

∆ g φ i = λ i φ i ,
where the eigenvalues are negative and decreasing

0 = λ 0 ≥ λ 1 ≥ • • • ≥ λ i .
Moreover, if M is a smooth manifold the eigenfunctions are infinitely differentiable functions. Moreover they provide an orthonormal basis of L 2 (M, g) ordered from low to high frequency, as the total integral of the gradient equals the absolute value of the eigenvalue. Meaning that if M f 2 dµ g < ∞ then f can be written has a linear combination of the functions {φ i } i≥0 :

f = lim k→+∞ k i=0 f, φ i (M,g) L 2 φ i .
In the computer graphics literature, the basis {φ i } i≥0 is often understood as a generalization of Fourier basis on Manifolds. An example is shown in Figure 3.4. They are sometimes referred to as manifolds harmonics [START_REF] Taubin | A signal processing approach to fair surface design[END_REF].

Commutativity with isometries

The Laplacian heavily depends on the metric. Like gradients, change of variables in Laplace-Betrami operators is equivalent to changing the metric tensor. Lemma 3.1.7 Let ϕ : N → M be a diffeomorphism between manifolds and g a metric on N , then: Proof To lighten the notations we denote g ⋆ the pullback metric (ϕ -1 ) ⋆ g. Given the local coordinates (x 1 , . . . , x n ) around p ∈ N and (y 1 , . . . , y n ) around q = ϕ(p) ∈ M , the Laplace-Beltrami operator can be transported from N to M using Lemma 3.1.6:

∆ g (f • ϕ) = (∆ (ϕ -1 ) ⋆ g f ) • ϕ, ∀f ∈ C ∞ (M ).
N g • ϕ∆ g (f • ϕ) dµ g = - N g p (∇ g (g • ϕ), ∇ g (f • ϕ)) dµ g (p) = - N g p dϕ -1 ϕ(p) (∇ g ⋆ g) , dϕ -1 ϕ(p) (∇ g ⋆ f ) dµ g (p) = - M g ϕ -1 (q) dϕ -1 q (∇ g ⋆ g) , dϕ -1 q (∇ g ⋆ f ) dµ g ⋆ (q) = - M g ⋆ (∇ g ⋆ g, ∇ g ⋆ f ) dµ g ⋆ = M g∆ g ⋆ f dµ g ⋆ = N g • ϕ (∆ g ⋆ f ) • ϕ dµ g .
The fundamental lemma of calculus of variations leads to the equality between the functions.

As a consequence of Lemma 3.1.7, the Laplace-Beltrami operator commutes with isometries. Moreover, the converse also holds: the Laplacian determines the metric up to isometries. This property is well-known in the computer graphics community and makes this operator very attractive for shape analysis. For instance, Proposition 3.1.8 ensures that the Heat Kernel Signature [START_REF] Sun | A concise and provably informative multi-scale signature based on heat diffusion[END_REF] and Wave Kernel Signature [START_REF] Aubry | The wave kernel signature: A quantum mechanical approach to shape analysis[END_REF] fully encode intrinsic information.

Proposition 3.1.8 Let ϕ : N → M be a bijection between manifolds then it is an isometry if and only if:

∆ g N (f • ϕ) = (∆ g M f ) • ϕ, ∀f ∈ C ∞ (N ).
Proof We prove the necessary and sufficient conditions separately:

• Assuming that ϕ is an isometry and using Lemma 3.1.7 lead to the equalities

∆ g M f = ∆ (ϕ -1 ) ⋆ g N f = ∆ g N (f • ϕ) • ϕ -1 .
• Now suppose that the Laplacians commute with the map ϕ. The operators on both manifolds are therefore isospectral and the eigenfunctions of one operator are eigenfunctions of the other after composition with ϕ. It follows that the heat kernel [START_REF] Rosenberg | The Laplacian on a Riemannian Manifold[END_REF], defined as h t (p, q) = i≥0 exp(tλ i )φ i (p)φ i (q), is preserved by the map:

h N t (p, q) = h M t (ϕ(p), ϕ(q)), ∀p, q ∈ N, ∀t > 0.
This equality holds true if and only if ϕ is an isometry, as detailed in [START_REF] Sun | A concise and provably informative multi-scale signature based on heat diffusion[END_REF].

Discretization

To be defined correctly, differential operators require smooth surfaces. Obviously, in a discrete setting smoothness has to be redefined since most of the data (scanner, MRI) are point clouds or meshes. We will limit our theoretical study to triangular meshes.

Although in practice the functional map framework is stable across data representations, all of the results presented in this thesis assume that the data comes in the form of triangle meshes. Therefore, in this section we show how differential operators are approximated using the standard Finite Element Method.

Discrete Manifold

Triangle meshes are composed of three sets (X , E, F) creating a discrete manifold. The set of vertices X whose elements x i ∈ R n are points in space. The set of edges E linking two vertices together forming an undirected graph. If a pair (i, j) belongs to E then the edge vector is denoted e ij = x j -x i where the order of the indices define the orientation.

We require an additional structural element that is the graph should be composed only by triangles such that an edge is shared by only two triangles. A triangle is defined by the indices of three vertices (i, j, k) belonging to the set F. Each side consists of an edge, so {(i, j), (j, k), (k, i)} belongs to E. The notations are summarized in Figure 3.5a. Meshes are thought of a piecewise linear approximations of a smooth manifold. When refined by adding sample points, the approximation error should decrease and eventually go to zero. On such discrete structure it seems natural to consider only piecewise linear functions. So, a discrete function f is defined as an assignment of a scalar value per vertex f : X → R and extended linearly inside triangles. For a point p inside the triangle (x i , x j , x k ) the interpolation is expressed using the hat basis composed of the functions B i defined in Figure 3.5b:

f (p) = f i B i (p) + f j B j (p) + f k B k (p).
The measure µ(T ) denotes the area of triangle T and ℓ ij = x i -x j the edge length. These two quantities will be used to describe the intrinsic structure of a discrete manifold.

x j x i x k x l α ij β ij T ilj w ij `ik T ijk `il `jl `jk (

Finite Element Method and Cotangent Weight Formula

For a short but efficient introduction to finite elements for computer graphics the reader may enjoy [START_REF] Botsch | Polygon Mesh Processing[END_REF] Chapter 3.

Inner product

Once the interpolation is written with the hat basis we can define the discrete inner product ., . L 2 by computing the integration of piecewise linear functions. Functions on vertices are elements of R X so the inner product is a symmetric definite matrix A, called mass matrix, of size |X | × |X |. The element A ij is then the inner product between B i and B j :

A ij = T ∈F T B i (p)B j (p)dp.
For computational efficiency it is sometimes preferable to use the lumped mass matrix. This matrix is diagonal whose elements are j A ij .

Gradient

Piecewise-linear functions are only differentiable inside triangles. Therefore, the discrete gradient takes as input a function defined at vertices and outputs a vector per triangle. By linearity, it is enough to compute the gradient of the hat functions to find the gradient of any function [START_REF] Botsch | Polygon Mesh Processing[END_REF]:

∇f (x) = 1 2µ(T ijk ) R 90 • E T ijk f k -f j f i -f j , (3.6)
where R 90 • is the counterclockwise rotation by 90 • around the normal and E T ijk = x j -x i , x k -x j contains a local basis of edges.

Laplacian operator

The discretization of the Laplace-Beltrami operator L is obtained by discretizing the inner product M ∇f, ∇g dµ represented by the stiffness matrix W ∈ R |X |×|X | and is deduced from W by integration by parts (Equation (3.5)) so L = A -1 W . The discrete inner product follows the discrete gradient in Equation 3.6:

g ⊤ W f = T ∈F 1 4µ(T ijk ) g k -g j g i -g j ⊤ E ⊤ T E T f k -f j f i -f j . (3.7)
After rearranging the terms the matrix becomes a graph Laplacian:

(W f ) i = j∼i w ij (f i -f j ), w ij =          1 8µ(T ijk ) (ℓ 2 ij -ℓ 2 jk -ℓ 2 ki ) (i, j) ∈ E + 1 8µ(T ijl ) (ℓ 2 ij -ℓ 2 jl -ℓ 2 li ), -k∼i w ki , i = j 0,
else.

(3.8)

The cotangent weight formula, a more popular expression of the Laplacian in the computer graphics community, is obtained by remarking that cot

α ij = 1 4µ(T ijk ) (-ℓ 2 ij + ℓ 2 jk + ℓ 2 ki ).
So the weight at an edge (i, j) ∈ E reads:

w ij = - 1 2 (cot α ij + cot β ij ). (3.9)
where α ij and β ij are the angles opposite of the edge.

Eigenfunctions

Given the discrete Laplacian, finding the eigen-decomposition is done by solving the generalized eigenvalue problem:

W Φ = AΦΛ, subject to Φ ⊤ AΦ = I,
where Λ is the matrix of eigenvalues. Since W is symmetric and A is positive definite, there exists a solution with non-negative eigenvalues. The matrices A and W are sparse so this problem is solved efficiently by standard sparse solvers.

Discrete Local Coordinates

Although the Finite Element Method has a nice interpretation as an approximation of a function space by piecewise linear functions, it does not relate immediately to the continuous concept of metric and local charts. In this section we make one link between those two worlds.

The metric tensor expressed as a matrix of edge lengths was first introduced in 1961 in [START_REF] Regge | General relativity without coordinates[END_REF] by Regge. It was originally intended for the study of discrete general relativity. This field is often referred as Regge calculus. For an accessible introduction see [START_REF] Khavari | Regge Calculus as a Numerical Approach to General Relativity[END_REF].

Discrete gradient

Triangle faces provide a straightforward analog of tangent spaces as each triangle can be easily mapped to a 2D plane. Of course, unlike the continuous case, the transitions between coordinates are not smooth. In this section, we are going to construct the gradient operator from directional derivatives.

The directional derivative of a piecewise linear function f in the direction of an edge is, by definition, given by the classical finite difference formula:

∇f, e ij e ij = f j -f i e ij .
We denote E T = (e ij , e jk ) a particular coordinate system at a given triangle. Inside a triangle the derivative of f can only be computed along the edges, therefore we store them in the matrix ∂f = f j -f i f k -f j . Given a vector V in the discrete tangent space and its decomposition in local coordinate V = E T α, the directional derivative is expressed through a matrix multiplication:

V (f ) = α ⊤ ∂f = ∇f, V .
This formulation is the discrete counterpart of Equation (3.1). In the direction of the local basis we have the equality: E ⊤ T ∇f = ∂f . To obtain a gradient we need to invert E T knowing that gradients are tangent vectors. To do so we are going to use the metric induced by the ambient space defined by:

g T := E ⊤ T E T ,
as suggested by Equation (3.4). This leads to the following defining Proposition.

Proposition 3.2.1 Given a local coordinate system E T = (e ij , e jk ) in triangle T , the gradient in local coordinate reads g -1 T ∂f and in the global coordinate system:

∇f = E T g -1 T ∂f.
In this discussion we merely computed the derivative of a piecewise linear function therefore it should correspond to the Finite Element Method.

Discrete metric

The metric tensor possesses several properties that will prove useful in Chapters 7 and 8 when converting functional representation of deformation into embedded surfaces. Proposition 3.2.2 Given a local frame E T = (e ij , e jk ) the metric tensor g T = E ⊤ T E T is symmetric positive definite and its determinant is such that:

det(g T ) = 4µ(T ) 2 .
Moreover the metric tensor can be rewritten exclusively in term of edge lengths:

g T = 1 2 2ℓ 2 ij ℓ 2 ki -ℓ 2 ij -ℓ 2 jk ℓ 2 ki -ℓ 2 ij -ℓ 2 jk 2ℓ 2 jk .
Proof The discrete metric is by definition symmetric positive definite for non-degenerate triangles. The computation of the determinant follows Heron's formula.

The second statement comes from the fact that the edges sum to zero e ij +e jk +e ki = 0 so:

ℓ 2 ki = e ij + e jk , e ij + e jk = ℓ 2 ij + ℓ 2 jk + 2 e ij , e jk .
Equivalence with Finite Element Method

The discretization using the Finite Element Method, Equation (3.6), and using the local coordinates, Proposition 3.2.1, are the same. First let's recall this simple property satisfied by every invertible symmetric matrix in R 2×2 :

g -1 T = 1 4µ(T ) 2 D ⊤ g T D,
where

D = 0 1 -1 0 . (3.10)
Starting from Proposition 3.2.1 and using Equation (3.10) lead directly to the conclusion:

∇f = E T g -1 T ∂f = 1 2µ(T ) ( 1 2µ(T ) E T D ⊤ E ⊤ T )E T D∂f = - 1 2µ(T ) R 90 • E T D∂f.
The inner product matrix W also enjoys a formulation in terms of local coordinate equivalent to Equation (3.7):

g ⊤ W f = T ∈F g j -g i g k -g j ⊤ g -1 T f j -f i f k -f j µ(T ).

Extension to simplicial complexes

One advantage of Regge calculus is that metric tensors and gradients can be easily defined for higher dimensional simplicial complexes. Moreover, they are still computable (but not necessarily meaningful) for non-manifold triangulations while being equivalent to the Finite Element Method in nicer scenarios.

x 1 x 2

x 3

x 4 E S 3 For example in Figure 3.6 a k-simplex S k is defined by k + 1 vertices (x 1 , . . . , x k+1 ). The first k edges form the local basis:

E S k = (x 2 -x 1 , . . . , x k+1 -x 1 ). The discrete metric tensor g S k = E ⊤
S k E S k is expressed locally using the edge lengths:

g ij = ℓ 2 1i , i = j 1 2 (ℓ 2 ij -ℓ 2 1i -ℓ 2 1j ), i = j. , 1 ≤ i, j ≤ k.
The volume of the simplex is accessible through the determinant of the metric by µ(S k ) = det(g S k )/k!. The gradient inside the simplex is now computed with the formula:

∇f = E S k g -1 S k   f 2 -f 1 . . . f k -f 1   .
Given a pure simplicial k-complex S, the inner product matrix W can be assembled using the following rule:

g ⊤ W f = S k ∈S   g 2 -g 1 . . . g k -g 1   ⊤ g -1 S k   f 2 -f 1 . . . f k -f 1   µ(S k ).

Chapter 4

Operator Representation

This chapter reviews two tools of the operator based representation framework: functional maps [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF], representing diffeomorphisms, and shape differences [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF], encoding the intrinsic deformation. The next chapters of this thesis explore the possibilities offered by these representations therefore this introduction attempts to give a complementary background to the reader and highlight open problems in this framework. We focus the presentation on the theoretical properties of these operators and how they are preserved after discretization. Although most of the results discussed here were already provided in the original papers, some are novel and not mentioned in follow-up works. In particular, Theorem 4.1.1 gives a characterization of composition operators never used in computer graphics and Propositions 4.2.3 and 4.2.4 describe the precise conditions in which compositions and summations are allowed for shape differences.

Functional Map

Real-valued functions are easily generated on many types of shapes (coordinate functions, descriptors) plus they enjoy interesting algebraic structure as they can be added, multiplied and composed with maps. To take advantage of this structure one can represent a diffeomorphism ϕ : N → M as a composition operator C ϕ acting on real-valued functions:

C ϕ : L 2 (M ) -→ L 2 (N ) f -→ f • ϕ.
Instead of studying directly the deformation induced by ϕ, composition operators assess the changes in the space of square-integrable functions L 2 defined on the manifolds M and N . These operators have been studied and used in many fields including dynamical systems (Koopman operator theory [START_REF] Henrion | Applied koopmanism[END_REF]) and functional analysis [START_REF] Singh | Composition Operators on Function Spaces[END_REF]. The immediate property of this point of view is that a composition acts linearly with respect to functions. In the discrete settings, composition are therefore represented by matrices, meaning that the space of diffemorphisms can be analyzed with standard linear algebra tools. Of course, all linear operators do not represent a composition (for example f → 0 is not a composition) but it suggests a convexification of the space of diffeomorphisms as linear mappings. In computer graphics composition operators were first introduced in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] under the name of functional maps. They were originally used in intrinsic shape matching problems and led to many follow-up works. Two directions are now prevalent: either to make better use of the structure of the space (e.g. [START_REF] Rodolà | Partial functional correspondence[END_REF][START_REF] Rodolà | Dense non-rigid shape correspondence using random forests[END_REF][START_REF] Ovsjanikov | Shape matching via quotient spaces[END_REF]) or to extend the functional representation to other quantities (e.g. metric distortion [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF], tangent vector fields [START_REF] Azencot | An operator approach to tangent vector field processing[END_REF], fluid mechanics [START_REF] Azencot | Functional fluids on surfaces[END_REF]).

In this section the term composition operator refers to the composition with an underlying map and should not be mistaken for functional map designating a linear operator acting on real-valued functions.

Mathematical Properties

All results presented here were first stated in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] except for Theorem 4.1.1. Although this fundamental property is well-known in operator theory [START_REF] Singh | Composition Operators on Function Spaces[END_REF], it has not been used in a computer graphics context. A more comprehensive study and careful computations can be found in [START_REF] Singh | Composition Operators on Function Spaces[END_REF].

Characterization

A functional map is a linear operator between function spaces. A first question is what characteristic properties a linear operator should satisfy to be a composition operator.

The following theorem provides an answer through a single property of the composition: the composition of a multiplication is the multiplication of the compositions.

C(f g) = C(f )C(g), ∀f, g ∈ L 2 (M ).
A sketch of the proof is formally reproduced here as it helps building an intuition (for a complete version see [START_REF] Singh | Composition Operators on Function Spaces[END_REF] Theorem 2.1.13). The necessary condition is obviously satisfied by a composition operator. The sufficient condition is proven in three steps. The first part is to show that the image of an indicator function on M is an indicator function on N . Second, C must define a homomorphism φ between algebras of measurable sets and third φ induces a diffeomorphisms ϕ.

Let χ U be the indicator function of a measurable set U ⊂ M . The necessary condition of Theorem 4.1.1 yields:

C(χ U ) = C(χ 2 U ) = C(χ U )C(χ U ) = (C(χ U )) 2 .
Therefore C(χ U ) is also an indicator function. Any measurable set U ⊂ M is mapped to a set V ⊂ N through C, defining an map between sets by: φ(U ) = V .

Let U 1 and U 2 be two disjoint subsets of M , then:

C(χ U 1 ∪U 2 ) = C(χ U 1 + χ U 2 ) = C(χ U 1 ) + C(χ U 2 ) = χ V 1 + χ V 2 .
Thus the image of a union of sets by φ is the union of the images. Similarly φ preserves intersections and differences between sets. So it can be shown that φ is a homomorphism between σ-algebras.

Finally, it follows, using a technical property of Borel sets (Theorem 2.1.12 in [START_REF] Singh | Composition Operators on Function Spaces[END_REF]), that there exists a measurable transformation ϕ such that φ(U ) = ϕ -1 (U ).

Algebraic properties

Many algebraic properties of the space of diffeomorphisms are naturally transferred to the space of composition operators. In a slightly simpler form, however, since we are dealing with linear operators.

The most basic property we can hope for is composition of mappings. Obviously, a composition of two mappings is also a composition operator since

C ϕ (C ψ (f )) = f •ψ •ϕ = C ψ•ϕ (f ).
Thus, a composition of mappings becomes a composition of linear operators and we have proved the following: Proposition 4.1.2 Let ϕ : N → M and ψ : M → O be mappings between manifolds then their composition is the composition of the composition operators:

C ϕ • C ψ = C ψ•ϕ .
As a consequence of Proposition 4.1.2, an invertible mapping ϕ becomes an invertible composition operator by noting that

C ϕ • C ϕ -1 (f ) = f . Corollary 4.1.3 Let C ϕ : L 2 (M ) → L 2 (
N ) be a composition operator then C ϕ is invertible if and only if ϕ is invertible. Moreover, we have:

C -1 ϕ (f ) = C ϕ -1 (f ), ∀f ∈ L 2 (N ).

Intrinsic information

Composition operators contain valuable information about the intrinsic structure of the manifold. For example Proposition 3.1.5 states that the pushforward measure describes the area distortion of the underlying deformation. Introducing the composition operator in the change of variables formula leads to the equality:

C ϕ (f ), C ϕ (g) M L 2 = M f g det((ϕ -1 ) ⋆ g N ) det(g M ) dµ g M . (4.1)
Thus the linear operator

C ⋆ ϕ • C ϕ (f ) : L 2 (M ) → L 2 (M )
, where the adjoint operator

C ⋆ ϕ is defined by C ⋆ ϕ (f ), g L 2 (M ) = f, C ϕ (g) L 2 (N )
, accounts for the multiplication of a function by the ratio of local areas. Thus, C ⋆ ϕ • C ϕ is an operator representation of the area distortion which will reappear in Section 4.2.

Proposition 3.1.8 completes the intrinsic information given by composition operators: composition operators representing isometries should commute with the Laplacian. • ϕ is area preserving if and only if:

C ⋆ ϕ (C ϕ (f )) = f, ∀f ∈ L 2 (M ).
• ϕ is an isometry if and only if:

C ϕ (∆ N f ) = ∆ M C ϕ (f ), ∀f ∈ L 2 (M ), ∆ N f ∈ L 2 (N ).
Theorem 4.1.4 shows that composition operators fully characterize intrinsic geometry. However, as such, we cannot use them to compare metrics of different manifolds since by definition compositions start and end at different function spaces. For intrinsic distortion analysis we will prefer the shape difference operators introduced in Section 4.2.

Discrete Functional Maps

In the discrete setting the space L 2 is replaced by an approximation with piecewise linear functions -a standard choice in FEM. A function is then a vector in R |X | assigning a real value per vertex. Therefore, functional maps will be understood as matrices. Let (φ M 1 , . . . , φ M |X | ) be an orthonormal basis of such a space (i.e. φ ⊤ i Aφ M j = δ ij where A is the mass matrix), then any function f is uniquely defined by the coefficients

a M i = f ⊤ A M φ M i in the decomposition f = |X | i=1 a M i φ M i .
The functional map is represented by a matrix transporting coefficients on M to coefficients on N :

C ϕ (f ) = |X | j=1 a M j C ϕ (φ M j ) := |X | i,j=1 C ij a M j φ N i .
The elements of C ϕ are by definition the coefficients of the functions C ϕ (φ M j ) when decomposed in the orthonormal basis φ N . This leads to the scalar products:

C ij := (φ N i ) ⊤ A N P M N φ M j , 1 ≤ i, j ≤ |X |.
where A N is the mass matrix an N and P M N assigns vertices on M to vertices on N . Let Φ M be the matrix containing all the basis functions φ M i , the discrete functional map reads:

C M N = Φ ⊤ N A N P M N Φ M . (4.2)
Representing a given vertex-to-vertex assignment by a functional map is usually done in two steps: This representation is flexible as the choice of basis can be adapted to each application. However, this is still an ongoing research topic. Some authors for example suggest a compactly supported basis [START_REF] Neumann | Compressed manifold modes for mesh processing[END_REF]. In this thesis, we will mostly focus on the popular choice of the basis of eigenfunctions of the Laplace-Beltrami operator.

Discrete properties

Every matrix does not represent a valid functional map since it has to satisfy Proposition 4.1.1, namely the composition of a product should be a product of the compositions. In fact, the product of two piecewise linear functions is not piecewise linear. So, this property cannot hold exactly when considering classic Finite Elements. This problem will not be discussed further in this thesis even thought it could improve the operator representation of mappings.

Nevertheless many of the other properties are still satisfied exactly at the condition that the basis span the entire space of discrete functions, namely Φ ⊤ AΦ = I and ΦΦ ⊤ A = I. Using Equation (4.2) the following operations hold:

• A composition of mappings becomes a multiplication of matrices as suggested by Proposition 4.1.2: The main theorem in [START_REF] Zeng | Discrete heat kernel determines discrete Riemannian metric[END_REF] proves that the stiffness matrix uniquely defines the edge length up to global scaling. This strong result is mostly due to the rigidity of triangle meshes. In particular,

C P N C M P = (Φ ⊤ N A N P P N Φ P )(Φ ⊤ P A P P M P Φ M ) = Φ ⊤ N A N P P N P M P Φ M = C M N . • Proposition 4.1.3 remains true since C N M = Φ ⊤ M A M P -1 M N Φ N is the inverse of the functional map C M N = Φ ⊤ N A N P M N Φ M . • Theorem 4.1.
C M N Φ ⊤ M W M Φ M = Φ ⊤ N W N Φ N C M N if
and only if M and N have identical triangle inner angles.

In practice a basis spanning the full space is not desirable as it requires |X | functions when shapes can contain several dozen thousand vertices. Besides, if the discrete manifolds have different connectivity, an accuracy up to a triangle is not relevant in many applications. It is often a better choice to use a small (compared to the number of vertices) family of smooth functions so the functional map is resilient to discretization noise.

Reduced basis

The eigenfunctions of the Laplace-Beltrami operator form a basis of the space of squareintegrable functions. Moreover, they are naturally ordered according to their smoothness since they are minimizers of the Dirichlet's energy [START_REF] Vallet | Spectral geometry processing with manifold harmonics[END_REF].

Therefore, we set our reduced basis to be the first k eigenfunctions of the Laplace-Beltrami operator. We are not considering the full function space but only a reduced space of smooth functions. This way we are compressing the information into a k N × k M matrix by keeping only the low frequency structure of the diffeomorphism. Nevertheless, this basis has several interesting properties regarding the representation of nearly isometric deformations. Since the Laplacian commutes with isometries, nearly isometric maps are represented by nearly diagonal matrices. Furthermore, a recent paper by Rodolà and al. [START_REF] Rodolà | Partial functional correspondence[END_REF] remarks that even when large parts of a shape are missing the functional map has a slanted-diagonal structure, making this representation useful for challenging problems such as partial matching. Regarding algebraic operations on functional maps, all the discrete properties mentioned above still hold but only in a reduced space of functions. In the reduced basis of the LB eigenfunctions the functional map acts as a low-pass filter on M and on N . Two failure cases can be identified depending on k N and k M the size of the bases. First scenario, the function do not lie in the span of Φ M so the projection onto the basis distorts the function. The projection, however, is well transferred to N . Second scenario, f is represented by the family Φ M but the transfer induces high frequency distortions lying outside of the span of Φ N . Increasing k N solves this issue. Figure 4.2 illustrates the effect of varying k M and k N on function transfer. In conclusion, a function f is accurately transferred if f and P M N f are well represented in their respective basis. We can characterize this subspace by:

f ∈ R |X | : f ∈ span(Φ M ), P M N f ∈ span(Φ N ) .
The problem of identifying the subspace of stable functions is tackled in Chapter 5.

k M = 5 k M = 50 k M = 100 k N = 5 k N = 50 k N = 100 k M k N f f • ϕ C M N f wrt k N , k M L 2 relative error

Other problems in the functional map framework

In this section we have only considered the problem of representing a given mapping between surfaces. The study of approximating an unknown functional map with minimal distortion using functional correspondences is covered in Chapter 5.

One of the challenge in this framework is to convert a functional map into a point-topoint map used in many applications. In [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] the authors propose a way to solve this problem by transferring highly localized functions. However, those functions are not well represented in a reduced basis introducing noise and destroying the smoothness of the recovered map. We propose a way to address this issue in Chapter 6.

Shape Differences

Characterization and operations on deformations are a fundamental tools in geometry processing with many applications, including deformation design, shape search and the organization of shape collections. The composition operator introduced in the previous section provides a lot of information on the intrinsic structure of a diffeomorphism as shown in Theorem 4.1.4. However, composition is not a convenient tool to compare mappings or metrics. For example, if M and N 1 , N 2 are shapes, by definition the functional maps ending points are respectively L 2 (N 1 ) and L 2 (N 2 ). As such it seems difficult to compare the deformation between M and N 1 with the deformation from M to N 2 without a map linking N 1 to N 2 . Moreover, an isometric deformation does not change the metric but is represented by a functional map different from identity. Thus, we need a representation of intrinsic deformation invariant under composition with isometries.

The shape differences address some of these issues. They are operators acting on functions and representing the intrinsic distortion induced by a given functional map. They are also symmetric allowing comparison between deformations and invariant under isometric composition.

Definition

Introduced by [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF], the shape difference operators describe a shape deformation by considering the change of two inner products between functions f, g M L 2 := M f g dµ and f, g M H 1 0 := M g(∇f, ∇g) dµ. Namely, given a pair of shapes M, N and a functional map C ϕ : L 2 (M ) → L 2 (N ) the authors introduce the area-based and conformal shape difference operators D A : L 2 (M ) → L 2 (M ) and D C : H 1 0 (M ) → H 1 0 (M ) respectively, as linear operators acting on (and producing) real-valued functions on M implicitly via the following equations:

f, D A (g) M L 2 := C ϕ (f ), C ϕ (g) N L 2 ∀f, g ∈ L 2 (M ) f, D C (g) M H 1 0 := C ϕ (f ), C ϕ (g) N H 1 0 ∀f, g ∈ H 1 0 (M ). (4.3)
By H 1 0 (M ) we denote the function space of square integrable functions with L 2 (M ) gradients and zero integrals:

H 1 0 (M ) = f ∈ L 2 (M ) : M ∇f 2 dµ < +∞, M f dµ = 0 .
This space seems natural when studying the conformal shape difference since D C maps any constant function to zero. When equipped of the scalar product ., .

L 2 + ., . H 1 0 , H 1 0 is a Hilbert space. The bilinear form (f, g) → f, g M H 1 0
is continuous and coercive thanks to the Wirtinger's inequality [START_REF] Brezis | Functional analysis, Sobolev spaces and partial differential equations[END_REF]. Moreover, for a given g in

H 1 0 (M ) the linear form f → C ϕ (f ), C ϕ (g) N H 1 0
is continuous assuming that C ϕ represents a composition with a diffeomorphism. All conditions of the Lax-Milgram theorem [START_REF] Brezis | Functional analysis, Sobolev spaces and partial differential equations[END_REF] are satisfied therefore there exists a unique D C (g) satisfying Definition 4.3. The existence and uniqueness of the area-based shape difference are guaranteed by the same argument.

It is important to note that the definition of shape differences does not require a composition operator but merely a functional map (i.e. a linear map between function spaces). The properties derived in the following, however, assumes there exists an underlying diffeomorphism ϕ : N → M .

Fundamental Properties

As shown by [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] those operators fully characterize intrinsic deformations. More precisely, the following equivalences link properties of the mapping to properties of the operators. Theorem 4.2.1 Given a pair of surfaces M, N and a diffeomorphism ϕ : N → M with the functional representation C ϕ , the followings hold for all functions f in H 1 0 (M ):

• ϕ is area-preserving if and only if

D A (f ) = f , • ϕ is conformal if and only if D C (f ) = f , • ϕ is an isometry if and only if D A (f ) = f and D C (f ) = f .
The last item is a direct consequence of the first two. The first item has already been proven by Theorem 4.1.4. Besides Equation (4.1) provides a direct representation of the operator D A :

D A (f ) = det((ϕ -1 ) ⋆ g N ) det(g M ) f, ∀f ∈ L 2 (M ). (4.4)
The proof of the second item requires a technical lemma proven in [START_REF] Schumacher | Conformal maps and p-dirichlet energies[END_REF] and reproduced here has it will be useful later on (in particular in Chapter 8). Lemma 4.2.2 generalizes the fundamental lemma of calculus of variations to symmetric tensors. Lemma 4.2.2 Let (M, g) be a Riemannian manifold (possibly with Lipschitz boundary) and let A be a differentiable symmetric 2-tensor field on M . Then A is the null tensor if and only if:

M g(∇f, A∇g) dµ g = 0, ∀f, g ∈ C ∞ (M ).
The proof of the last point follows after a few computations.

Proof of Theorem 4.2.1 Let g N be the metric on N and let g ⋆ denote (ϕ -1 ) ⋆ g N the pullback metric on M . Using Lemma 3.1.6 and the definition of the pullback metric, we have:

C ϕ (f ), C ϕ (g) N H 1 0 = N g N p dϕ -1 ϕ(p) (∇ g ⋆ f ) , dϕ -1 ϕ(p) (∇ g ⋆ g) dµ g (p) = M g N ϕ -1 (q) dϕ -1 q (∇ g ⋆ f ) , dϕ -1 q (∇ g ⋆ g) dµ g ⋆ (q) = M g ⋆ (∇ g ⋆ f, ∇ g ⋆ g) dµ g ⋆ .
The scalar product g ⋆ (∇ g ⋆ f, ∇ g ⋆ g) is rearranged to g M ∇ g M f, A∇ g M g where

A ij = (g ⋆ ) ik g kj and g M is the metric on M . The definition of the pushforward measure leads to the equation:

C ϕ (f ), C ϕ (g) N H 1 0 = M g M ∇ g M f, λA∇ g M g dµ g M , λ = det(g ⋆ ) det(g M ) . (4.5)
We are now ready to prove the equivalence.

• Suppose that ϕ is conformal then g M = λg ⋆ . In Equation (4.5) it implies that λA ij = δ ij so D C is the identity operator.

• Suppose that D C is identity then using Definition (4.3) and Equation (4.5) lead to:

M g M ∇ g M f, (Id -λA)∇ g M g dµ g N = 0, ∀f, g.
Thanks to Lemma 4.2.2 we conclude that g M = λg ⋆ .

Shape differences are interpreted as ratio of intrinsic quantities. The area-based operator is the ratio between the local area of the pullback metric with the local area of the reference metric. A similar explanation applies to the conformal operator, Equation

B ij = det((ϕ -1 ) ⋆ g N ) det(g M ) ((ϕ -1 ) ⋆ g N ) ik g M kj . (4.6)
Thus, D C describes the ratio of the pullback metric and the reference metric weighted by the local area. The tensor B is immune to multiplication of the metric by a positive function, so this measurement of distortion cannot be used to assess area changes. In that sense D A and D C are orthogonal and complementary when describing intrinsic deformations.

The shape difference operators are defined by the mean of a diffeomorphism but Equations (4.5) and (4.6) rely only weakly on it. So these operators can be extended to (and should be thought as) comparison of an arbitrary metric, not necessarily induced by an embedded surface, with the reference metric.

Algebraic Properties

The results of this section, namely Propositions 4.2.3 and 4.2.4, are novel and do not appear in the original paper [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] or any follow-up papers.

Linear combination

Unlike composition operator there is no known defining property equivalent to Theorem 4.1.1. However, Equations (4.4) and (4.6) provide a sense of the structure of the shape difference space.

Equation (4.4) tells us that the area-based shape difference is a multiplication operator associated to a positive function. Thus multiplying D A by a positive scalar represents a global scaling of the target metric. Moreover adding two area-based shape differences results in a linear operator representing a shape difference. However, the underlying metric may not correspond to an embedded surface.

For the conformal operator the situation is even simpler. From Equation (4.6) we take out that this operator is characteristic of a 2-tensor with determinant one. In general linear combinations do not produce a valid conformal shape difference.

Besides those algebraic considerations, the space of shape differences is vast as it represents any metric distortion. At the moment, constraining shape differences to represent the metric of an embedded surface is an open question.

Deformation support

The addition between operators is only possible in case of non-overlapping deformations. As a consequence of Theorem 4.2.1, the set Ω ⊂ M supports a deformation represented by the operator D (area or conformal-based) if and only if:

D(f ) = f, ∀f ∈ H(M ) : Suppf ∩ Ω = ∅,
where Suppf = {p ∈ M : f (p) = 0} denotes the support of a real-valued function on M . The following proposition specifies necessary conditions for summation of deformations.

Proposition 4.2.3 Let D 1 , D 2 : H(M ) → H(M ) be two shape difference operators with respective supports Ω 1 and Ω 2 such that Ω 1 ∩ Ω 2 = ∅, then the linear operator:

D + (f ) = D 1 (f ) + D 2 (f ) -f, ∀f ∈ H(M ),
is the shape difference operator corresponding to the deformation D i on Ω i and to an isometry elsewhere.

Proof The proof is restricted to the case of conformal shape difference as it is very similar for the area-based one. Using Equation (4.6) the operator D + is separated into a sum of three integrals:

f, D + (g) M H 1 0 = Ω 1 ∇f, B 1 ∇g dµ + Ω 2 ∇f, B 2 ∇g dµ + Ω c 1 ∪Ω c 2 ∇f, ∇g dµ.
Since the support of these integrals are disjoint they can be written as a single integral:

f, D + (g) M H 1 0 = M ∇f, B∇g dµ,
where B is the smooth tensor on M defined by:

B(p) =    B 1 (p), p ∈ Ω 1 B 2 (p), p ∈ Ω 2 Id(p), p ∈ Ω c 1 ∪ Ω c 2 .

Composition of diffeomorphisms

Another possible operation on shape difference operators is the composition of two (possibly overlapping) deformations. As the property introduced here applied to both types of operators we will denote D the operator corresponding to either D A or D C and H the function space L 2 or H 1 0 . The following proposition was proved in a restricted discrete setting for one operator in the original article [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF].

N M P ϕ • φ φ ϕ L 2 (N ) L 2 (M ) L 2 (P ) C ϕ•φ C φ C ϕ
D ϕ•φ = D ϕ • C -1 ϕ • D φ • C ϕ .
Proof The proof relies only on Definition (4.3). Let f, g be a pair functions in H(M ):

f, D ϕ•φ (g) M H = C ϕ•φ (f ), C ϕ•φ (g) N H = C ϕ (f ), D φ (C ϕ (g)) P H = f, D ϕ C -1 ϕ D φ (C ϕ (g)) M H .
Proposition 4.2.4 confirms the intuition that the shape differences are invariant under composition with isometries. Suppose that the diffeomorphism φ is an isometry then Theorem 4.2.1 tells us that D φ is the identity operator on H(P ). Using the composition formula we conclude that the distortion of the composition is equal to the distortion induced by ϕ alone:

D ϕ•φ = D ϕ • C -1 ϕ • D φ • C ϕ = D ϕ .
The shape difference operators associated to the inverse mapping ϕ -1 can be computed from D ϕ using the composition formula:

D ϕ -1 •ϕ = D ϕ -1 • C ϕ • D ϕ • C -1
ϕ . Formally, we have:

D ϕ -1 = C ϕ • (D ϕ ) -1 • C -1 ϕ .

Distortion pullback

An attractive property of this framework is the possibility to transport a distortion from a shape to another. Given a shape difference on shape P and a functional map between M and P , we would like to describe an equivalent distortion of the metric on M . Intuitively, the operator

D ⋆ = C -1 ϕ • D φ • C ϕ
is a good candidate for the pullback of the shape difference D φ encoding the metric deformation from g P to (φ -1 ) ⋆ g N . It follows from Proposition 4.2.4 that D ⋆ is characteristic of the same distortion but after a pullback on M , namely it compares (ϕ -1 ) ⋆ g P with ((ϕ

• φ) -1 ) ⋆ g N .
This operator pullback suggested various experiments in the original shape difference article [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF]. For example the authors find similar deformations across different shapes collections. Given a collection of cats in different positions the pullback shape differences are used to identify similar poses in a collection of lions.

Discrete Shape Differences

In the discrete setting, a functional map is defined by Equation (4.2). Thus, given orthonormal bases Φ M and Φ N , Definition (4.3) yields the matrices

D A , D C ∈ R k M : D A = C ⊤ M N C M N , D C = Φ ⊤ M W M Φ M -1 C ⊤ M N Φ ⊤ N W N Φ N C M N .
We recover the interpretation of shape differences as ratios. The area-based shape difference is the ratio between mass matrices and the conformal operator compares stiffness matrices.

When the basis is restricted to the first eigenfunctions of the Laplacian, the conformal shape difference has a slightly simpler expression ([106] Section 5 Option 2):

D C = Λ -1 M C ⊤ M N Λ N C M N ,
where Λ is the diagonal matrix storing the eigenvalues of the Laplacian. When represented by matrices the shape differences can be analyzed by standard linear algebra tools. For example, the eigenvalues offer a quantification of the distortion. According to Theorem 4.2.1, an eigenvalue close to one is associated to an eigenfunction highlighting isometric areas. Figure 4.4 shows the most distorted function in the simple deformation scenario.

Discrete properties

Theorem 4.2.1 is satisfied by the discretization. More precisely, in Chapter 7, we show that one can recover the discrete metric (edge lengths) from both shape differences by solving two linear systems of equations. In fact, this chapter proposes a solution to the broader problem of finding embedded surfaces from operators.

The algebraic operations of addition and composition described by Propositions 4.2.3 and 4.2.4 hold true for a smooth subspace of functions. This is a direct consequence of Section 4.1.2 showing that the algebraic properties of functional maps are preserved after discretization.

As for functional maps, basis reduction can generate artifacts depending on the basis sizes k M on M and k N on N . A low k N leads to a bad representation of hight frequency distortions since some functions in Φ M are not represented in Φ N after transfer. So the shape differences unjustly map to zero some functions. It follows that eigenvalues (and eigenfunctions) lower the one are less reliable than eigenvalues greater than one. A low k M implies that Φ M contains only low frequency functions so the shape differences tend to smooth hight frequency deformation. Overall Figure 4.5 shows that a tradeoff has to be found between k M and k N to represent correctly the distortion at a certain level of detail.

Organization of the Thesis

The rest of the thesis is organized in two parts. The first part is dedicated to the computation and analysis of functional maps, in particular applied to shape matching Source Target problems. In Chapter 4 we have assumed that the underlying diffeomorphism was given as a point-to-point correspondence so the functional map is an approximation of a composition operator. Chapter 5 tackles the converse problem: given two shapes we try to compute the functional map that fits a particular deformation model, namely the one that minimizes the intrinsic distortion. For this purpose we provide a supervised learning algorithm able to select features and identify functions that jointly produce the best functional map. To complete our tour on shape matching problems we move on to converting a functional map to a point-to-point representation. Although this problem was already considered in the original article, it gave no guarantee of recovering a continuous map. Chapter 6 suggests a solution based on repairing a given diffeomorphism by a vector field flow.

k M = 5 k M = 50 k M = 100 k N = 5 k N = 50 k N = 100 k M k N f D C (f ) D C (f ) wrt k M , k N L 2 relative error
The second part explore intrinsic deformation through the prism of shape differences. As proven in the previous section shape differences completely characterize the metric in the continuous setting. Chapter 7 expents the analysis to triangular meshes by proposing an algorithm to recover edge lengths from operators. Furthermore, shape differences are blind to some curvature deformations hence giving an incomplete description of embedded surfaces. A possible solution explored in Chapter 7 is to encode the metric of an offset surface leading to a coordinate free description of triangular meshes by four operators.

Shape difference operators were introduced for analysis purposes. However, their interesting algebraic structure makes them suitable for exploration and synthesis of new deformations. In Chapter 7 a first step is made in the direction of deforming shapes using an operator-based representation with one drawback: to obtain embedded surfaces we need first to extract the metric information and then reconstruct the mesh by solving an Multidimensional scaling problem which is computationally expensive and error prone. Chapter 8 looks at this problem from another angle and studies a characterization of deformation fields up to isometric deformation. For this, we first introduce a unified shape difference that fully characterizes all intrinsic changes. Then we consider the operator associated to an infinitesimal displacement, and characterize deformation fields by the distortion they induce on the metric. Interestingly these operators enjoy similar properties as shape differences (composition, transfer), which enable significantly simpler deformation reconstruction compared to the previous method, and enable novel deformation synthesis applications. We provide theoretical proofs of informativeness in both continuous and discrete settings.

Part I

Shape to Deformation

This first part is dedicated to the computation and analysis of functional maps. Namely, we use feature selection to improve solutions of shape matching problems in Chapter 5 and then propose an algorithm to convert functional maps to continuous point-to-point maps in Chapter 6.

Chapter 5

Supervised Descriptor Learning for Non-Rigid Shape Matching

In this chapter, we present a novel method for computing correspondences between pairs of non-rigid shapes. Unlike the majority of existing techniques that assume a deformation model, such as intrinsic isometries, a priori and use a pre-defined set of point or part descriptors, we consider the problem of learning a correspondence model given a collection of reference pairs with known mappings between them. Our formulation is purely intrinsic and does not rely on a consistent parametrization or spatial positions of vertices on the shapes. Instead, we consider the problem of finding the optimal set of descriptors that can be jointly used to reproduce the given reference maps. We show how this problem can be formalized and solved for efficiently by using the recently proposed functional maps framework. Moreover, we demonstrate how to extract the functional subspaces that can be mapped reliably across shapes. This gives us a way to not only obtain better functional correspondences, but also to associate a confidence value to the different parts of the mappings. We demonstrate the efficiency and usefulness of the proposed approach on a variety of challenging shape matching tasks.

Introduction

Finding high quality correspondences is a key component in many tasks including statistical shape analysis [START_REF] Hasler | A statistical model of human pose and body shape[END_REF], deformation transfer [START_REF] Sumner | Mesh-based inverse kinematics[END_REF] and interpolation (morphing) [START_REF] Kilian | Geometric modeling in shape space[END_REF] among others. While a number of efficient techniques have been proposed to address the problem of rigid alignment [START_REF] Tam | Registration of 3d point clouds and meshes: A survey from rigid to nonrigid[END_REF], the problem of general non-rigid shape matching remains difficult. Most existing methods for finding correspondences between non-rigid shapes rely on an a priori deformation model, which specifies the space of "reasonable" maps between shapes. Perhaps the most popular and widely used such model is that of approximate intrinsic isometries [START_REF] Bronstein | Generalized multidimensional scaling: a framework for isometry-invariant partial surface matching[END_REF][START_REF] Mémoli | On the use of Gromov-Hausdorff Distances for Shape Comparison[END_REF], where the mapping is assumed to preserve geodesic distances between all pairs of points on the shapes. A more general possibility is to consider conformal deformations, which are only assumed to preserve angles [START_REF] Lipman | Mobius voting for surface correspondence[END_REF][START_REF] Kim | Blended intrinsic maps[END_REF] or to parameterize the space of possible maps using a fixed deformation model [START_REF] Zhang | Deformation-driven shape correspondence[END_REF]. Although these techniques can produce good results when the deformation satisfies the a priori model, they can fail badly as soon as even moderate deviations from the model are introduced. This is especially critical since many natural deformations, such as articulated motion of humans or animals are known to induce potentially significant geodesic distortion [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF]. Incorporating the possibility for such distortion into a deformation model is challenging especially using a purely axiomatic (theoretical) approach.

Rather than trying to devise a theoretical deformation model capable of adapting to known deformations, several communities have tackled this challenge by using a datadriven approach, where the space of "reasonable" maps or deformations is learned from a set of examples, e.g. [START_REF] Cootes | Training models of shape from sets of examples[END_REF]. Since obtaining example deformations is often significantly easier than devising a unified theoretical deformation model, such an approach allows the resulting techniques to remain flexible yet efficient in the particular settings where they are applied. Most data-driven approaches for devising a deformation model, however, rely heavily on a consistent parametrization of the deformation domain (e.g. on a fixed grid in Euclidean space), and perform statistical analysis on the positions of vertices of the shapes [START_REF]Active shape models-their training and application[END_REF][START_REF] Bookstein | Landmark methods for forms without landmarks: morphometrics of group differences in outline shape[END_REF][START_REF] Dryden | Statistical shape analysis[END_REF][START_REF] Gao | A review of active appearance models, Systems, Man, and Cybernetics, Part C: Applications and Reviews[END_REF]. When computing correspondences between pairs of surfaces in 3D, such parametrization is often unavailable and moreover, shapes can undergo severe deformations which are difficult to capture using purely extrinsic approaches. In this chapter, we propose a purely intrinsic method for exploiting prior correspondence information between pairs of shapes to find better correspondences between a reference shape and a new previously unseen instance. Rather than doing the learning over, e.g., the positions of the vertices on the shapes, we propose to find the optimal set of point descriptors that can be jointly used to reproduce the given reference maps. While such an optimization is, in general, very complicated, since even to evaluate how well the descriptors can reproduce a given map would require a full solution of the shape matching problem, we show how this problem can be formalized and solved for efficiently by using the recently proposed functional maps framework [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF]. Moreover, we demonstrate how to extract the functional subspaces that can be mapped reliably across shapes. This gives us a way to not only obtain better functional correspondences, but also to associate a confidence value to the different parts of the mappings. Our approach is also quite general since it can be used as a preprocessing step of other methods using functional maps [START_REF] Ovsjanikov | Shape matching via quotient spaces[END_REF][START_REF] Huang | Functional map networks for analyzing and exploring large shape collections[END_REF][START_REF] Azencot | An operator approach to tangent vector field processing[END_REF] in order to improve the quality of the results and help to handle difficult deformation. Note that in this chapter we focus on the shape matching problem which is the most developed application of the functional maps.

Related Work

Non-rigid shape matching is a very-well developed area and its complete overview is beyond the scope of this chapter (see, e.g., [START_REF] Bronstein | Numerical geometry of non-rigid shapes[END_REF][START_REF] Van Kaick | A survey on shape correspondence[END_REF] for recent surveys of this field). We therefore concentrate on the work directly related to ours, namely near-isometric shape matching with special emphasis on approaches that utilise prior knowledge for establishing correspondences between pairs of shapes. The vast majority of techniques for non-rigid shape matching implicitly make use of a deformation model for finding correspondences between geometric shapes. Perhaps the most common model in the context of intrinsic (i.e., not relying on vertex positions and not assuming approximate alignment) approaches is approximate isometries, introduced by Bronstein et al. [START_REF] Bronstein | Generalized multidimensional scaling: a framework for isometry-invariant partial surface matching[END_REF] and Mémoli [START_REF] Mémoli | On the use of Gromov-Hausdorff Distances for Shape Comparison[END_REF]. This model has been used by a large number of methods, (e.g., [START_REF] Huang | Non-rigid registration under isometric deformations[END_REF][START_REF] Tevs | Isometric registration of ambiguous and partial data[END_REF][START_REF] Ovsjanikov | One point isometric matching with the heat kernel[END_REF][START_REF] Sahillioğlu | Coarse-to-fine combinatorial matching for dense isometric shape correspondence[END_REF][START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] among many others) that all assume that the sought correspondences must approximately preserve pairwise geodesic distances. Another set of approaches is based on a more relaxed model, conformal mappings, used by, e.g., [START_REF] Lipman | Mobius voting for surface correspondence[END_REF][START_REF] Kim | Blended intrinsic maps[END_REF] where only angles are assumed to be preserved. Other techniques, such as the one used by Zhang et al. [START_REF] Zhang | Deformation-driven shape correspondence[END_REF] explicitly deform a shape using a fixed deformation model to find correspondences between non-rigid shapes. All of these approaches use a model given a priori to find correspondences, which can be problematic if the real deformations do not agree with the given model. Interestingly, it has recently been observed [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] that even articulated motion of humans can induce noticeable isometric distortion, which could explain some of the difficulties encountered by previously proposed techniques. In contrast, other works have proposed to learn an appropriate deformation model from a set of examples, and then use this model for shape matching. Perhaps the best-known example of this approach are Active Shape and Active Appearance Models [START_REF] Cootes | Training models of shape from sets of examples[END_REF][START_REF]Active shape models-their training and application[END_REF] and their variants (see, e.g., [START_REF] Gao | A review of active appearance models, Systems, Man, and Cybernetics, Part C: Applications and Reviews[END_REF]) used widely in Computer Vision. In a similar vein, techniques in Statistical Shape Analysis [START_REF] Dryden | Statistical shape analysis[END_REF] use the distribution of positions of pre-specified landmark points in 3D to learn a statistical deformation model over which inference can be made. Related techniques are commonly used in medical imaging and Morphometrics [START_REF] Bookstein | Landmark methods for forms without landmarks: morphometrics of group differences in outline shape[END_REF] and in Geometry Processing communities, e.g. [START_REF] Anguelov | Scape: shape completion and animation of people[END_REF][START_REF] Hasler | A statistical model of human pose and body shape[END_REF] among many others. However, all of these methods assume the existence of a common domain over which learning can be made, and which most often is done using vertex coordinates of either landmark points or all points on a fixed reference shape. In the context of intrinsic shape matching, where shapes lack labeled landmark points and can undergo severe deformations, vertex coordinates are often not relevant, limiting the applicability of such techniques. Rather than relying on vertex positions, recent methods have considered using derived properties such as point or triangle descriptors for learning. Thus, Kalogerakis et al. [START_REF] Kalogerakis | Learning 3d mesh segmentation and labeling[END_REF] and Van Kaick et al. [START_REF] Van Kaick | Prior knowledge for part correspondence[END_REF] have proposed using a set of example shapes to train classifiers for part segmentation and labeling, which can then be used to establish partlevel correspondences. Similarly, Chen et al. [START_REF] Chen | Schelling points on 3d surface meshes[END_REF] explore the predictive power of various descriptors for detecting distinctive landmark (schelling) points identified by users. These methods, while similar to ours in learning on the level of descriptors do not, however, specifically address the shape matching problem. Perhaps most closely related to ours are recent works by Litman et al. [START_REF] Litman | Learning spectral descriptors for deformable shape correspondence[END_REF] and Rodolà et al. [START_REF] Rodolà | Dense non-rigid shape correspondence using random forests[END_REF], where the authors use a set of examples to learn the most informative descriptors that are used directly in the context of shape matching. Our approach is fundamentally different, however, since rather than trying to identify descriptors that can distinguish different points, we propose to find the optimal descriptor set that can be used to jointly produce the entire map across shapes. We thus avoid the problem of obtaining consistent correspondences present in these approaches (and obtained during post-processing), since consistency is incorporated directly in the learning stage. Crucially, we use the recently proposed functional map representation [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] that allows us to formulate the learning problem purely intrinsically, while permitting to directly control and optimize for the influence of descriptors on the quality of the final map.

Goals Given a collection of (training) shapes with known correspondences our goal is to identify the most informative descriptor set that can be used to solve the non-rigid shape matching problem on new (test) instances. Besides we want to learn where are the most stable correspondences.

Consistent Functional Maps

Our method is based on the functional map representation introduced in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF]. In this section, we give a brief overview of the representation and the method used in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] to construct a functional map for a given pair of shapes. While our method is general, throughout the chapter we assume that all shapes are represented as triangle meshes, and all functions are expressed as vectors in the basis of the eigenfunctions of the Laplace-Beltrami operator. This basis needs to be computed beforehand on each shape. The objective is to output a uniquely defined functional map.

Functional Map Representation

The functional map representation is based on the observation that given two surfaces M 0 and M i , a point-to-point map ϕ i : M i → M 0 induces a map between function spaces

C i : L 2 (M 0 ) → L 2 (M i ), where L 2 (M )
is the set of square integrable functions defined on the surface M . The functional map C i is defined by composition with ϕ i as

C i f = f • ϕ i .
The operator C i is a linear transformation and given a basis it can be represented as a matrix in the discrete setting. This matrix can be easily computed if the map ϕ is known. The basic method described in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] approximates the functional map C i using a set of linear constraints. The first type of constraints is given by a set of pairs of functions, which we refer to below as "probe functions", that are expected to be preserved by the deformation. The second is a regularization term coming from the deformation model. This leads to the least squares problem:

X i = arg min C CG 0 -G i 2 F + α C ⊙ W 2 F , (5.1)
where . F denotes the Frobenius norm. The use and meaning of each term will be detailed in the following paragraphs.

Probe Functions

The probe functions can be represented by two matrices G 0 and G i , where each pair of corresponding columns represents a pair of functions g 0 , g i such that C i g 0 ≈ g i is expected to hold for the unknown C i . In practice we normalize the corresponding column so that each column has the same L 2 norm and g 0 , g i contain the coefficients of the probe functions in a given basis (e.g. LB eigenfunctions) . Thus, the functional map should verify CG 0 ≈ G i . In the context of isometric matching the probe functions are given by classical descriptors, such as the HKS [START_REF] Sun | A concise and provably informative multi-scale signature based on heat diffusion[END_REF], or WKS [START_REF] Aubry | The wave kernel signature: A quantum mechanical approach to shape analysis[END_REF].

Regularization

In addition to the probe function constraints, the authors of [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] have proposed a regularization using the assumption that the deformation is nearly isometric. This assumption is equivalent to the commutativity of C i with the Laplace-Beltrami operator, namely

C i ∆ 0 = ∆ i C i .
In the discrete setting the eigenfunctions of the Laplace-Beltrami operator are used as function basis. Thus, this equation can be written as C i ⊙ W = 0 where "⊙" denotes the component-wise multiplication and the matrix W is defined by W kl = λ i k -λ 0 l with λ i k the k th eigenvalue of the Laplace-Beltrami operator on the surface M i . Other assumptions on the deformation model can be used. For example in [START_REF] Kovnatsky | Coupled quasi-harmonic bases[END_REF] the authors regularize the shape mathcing problem by imposing that the map should be area-preserving implying that the matrices C i are a-orthonormal as explained in Chapter 4. However, this assumption is weaker and leads to a non-convex problem more challenging to incorporate in our supervised learning algorithm.

Uniqueness of the solution. In practice the eigenvalues of the Laplace-Beltrami operator of two different shapes are always numerically different except for the zero eigenvalues. Thus, the only zero coefficient of W is W 11 which weights the coefficient C 11 of the functional map. Since the corresponding eigenfunctions are constant, C 11 maps the constant functions of L 2 (M 0 ) into the constant functions of L 2 (M i ). The coefficient C 11 should always be one. Therefore, since W is non zero everywhere, the solution of (5.1) is unique without any assumptions on G 0 . However, the probe functions G 0 and G i are in practice composed of symmetric functions, so we cannot hope to have information about the antisymmetric functions. Nevertheless, the problem (5.1) still has a unique solution that might simply map the antisymmetric functions poorly.

Main Challenge

In the original article [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] the probe functions are assumed to be given, so how to choose them was not discussed. As mentioned in introduction, this choice can already be challenging. For example in Figure 5.1a the smoothed Gaussian curvature computed on two different meshes provides a decent functional correspondence. At the same time, in Figure 5.1b the logarithm of the Gaussian curvature, while intrinsic in theory, does not result in a useful correspondence. One option to identify the best descriptors would be to simply find the most stable probe functions in the example (training) set, by learning spectral descriptors [START_REF] Litman | Learning spectral descriptors for deformable shape correspondence[END_REF] for example. However, some descriptors (e.g., the constant function) can be stable without at all being informative. More importantly, however, as can be seen from Equation (5.1), the descriptors influence the resulting functional maps X i jointly. As an example, if a correspondence is described by several probe functions the resulting functional map will tend to respect this constraint while other meaningful correspondences will be arbitrarily put aside due to their low redundancy. So picking the best descriptors independently will not necessarily result in high quality maps. Thus, the key idea developed in this chapter is to introduce weights for probe functions, over which learning can be done. As explained below, the probe function constraint will be replaced by:

CG 0 D -G i D 2
F where D is an unknown diagonal matrix of weights.. The weights D will be optimized so that the weighted descriptors are jointly as informative as possible. This will allow us to improve the quality of the functional maps and to extract the most stable functional subspaces.

Algorithm Outline

We propose a two-step method described in following two sections and summarized in Figure 5.2. Given collection of shapes, we learn the most informative set of weights D by solving an optimization problem. We then extract a functional basis whose components are ordered by quality of correspondence. When given a previously unseen shape, we use this information to compute a high-quality functional map using the optimal weights and to discard the badly mapped functions by reducing the functional space.

Selection of the Best Functional Correspondences

The idea developed here is to assign a weight to each pair of probe functions. These weights can then be tuned according to their consistency in the matching. Since a priori there is no reason to choose one probe function over another, we propose to learn the optimal weights given a training set of shapes. As input we need a set of n triangulated meshes with known correspondences representing the same object undergoing a set of deformations. Our main assumption is that the optimal weights on the probe functions should be stable across the shapes in the collection. Thus, if we are given a new deformation of the same shape, the learned weight should also select the consistent probe functions. The output of our algorithm will be a set of weights for the probe functions, which, as we will show below, can then be used to find correspondences between new, unseen shape instances.

G i Y D ⋆

Probe Functions Weight Optimization Basis Learning

Training Shapes

D ⋆ Learning G N +1 X p N +1 X N +1 (D ⋆ ) Probe Functions Weighted Functional Map Space Reduction New Shape Testing C i Functional Maps Learnt Weight X i (D ⋆ ) C i Functional Maps Y Learnt Basis
Section 4 Section 3

Weighting the probe functions

As described above, our main idea is to introduce a set of weights on the functional correspondences to measure their usefulness in finding a relevant map by using a diagonal matrix D. For a given weight matrix D, the linear constraints given by the probe functions become CG 0 D = G i D. We can then define the function X i (D), which maps a given sets of weights to the corresponding functional map, via the solution of the optimization problem:

X i (D) = arg min C CG 0 D -G i D 2 F + α C ⊙ W 2 F (5.2)
We choose here to fix α and tune D. We could also try to tune all the parameters (α and D) but the coefficients would be defined up to a multiplicative constant and C(D) may no longer be well-defined when α is equal to zero. Since all the functional maps start from the reference shape M 0 , this shape obviously plays special role in our method. Ideally we would like to take as reference the most "average" shape of the collection. Following this idea, a simple procedure is presented in [START_REF] Shapira | Cross-collection map inference by intrinsic alignment of shape spaces[END_REF] to find the shape of the collection which minimizes the average isometric distortion. However, in our experiments we chose the symmetric standing pose as reference. Note that an interesting future work would be to use a more complex graph of correspondences between the shapes of a collection and impose cycle consistency as proposed in [START_REF] Huang | Functional map networks for analyzing and exploring large shape collections[END_REF].

As discussed in the previous section X i (D) is well-defined and differentiable. Note that the weight matrix D has a global effect on X i (D), and the problem (5.2) cannot be separated in terms of the individual components of D.

Finding the best weights

Our goal is to find a set of weights which will provide a relevant set of functional maps. The main hypothesis is that the weights should be stable across deformations that belong to the same category. That is, the most relevant probe functions for the training shape pairs will also be well-suited to find a functional map to new unseen shapes in that category.

Learning from a given collection. We assume that we are given a collection of n nearly isometric deformations of the same object with known functional maps C i . The optimal weights D ⋆ are the ones that produce an approximation X i (D) that is closest to the ground truth C i . Thus, we want to solve the following optimization problem:

D ⋆ ∈ arg min D n i=1 X i (D) -C i , (5.3)
where the sum is over the set of given training maps C i . Interestingly, the optimization problem 5.3 is closely related to optimal control problems as the functional maps X i (D) are themselves solutions of an optimization problem. Similar setting appears in [START_REF] Litman | Learning spectral descriptors for deformable shape correspondence[END_REF] where a feedback loop is used to achieve a supervised dictionary learning. Note that the choice of the norm is important. We would like the functional map X i (D) to match C i over as-large-as possible functional subspace. This is equivalent to minimizing the rank of the difference X i (D) -C i . Thus, the important quantities are the singular values of the differences

X i (D) -C i .
The naive choice of the squared Frobenius norm is not well suited for our problem since

A 2 F = σ(A) 2 2
where σ(A) is the vector containing the singular values of the matrix A. Therefore this norm would give a large weight to the biggest singular values, which correspond to the worst-matched functional subspaces. Among these subspaces is the space of antisymmetric functions that we have no hope of mapping since the probe functions give us very little information about this subspace. At the same time, the small singular values have little influence on the minimization whereas they are the ones we would like to optimize.

The choice of the norm. To tackle the rank minimization problem we choose the following norm which is a regularization of the l 0 -norm:

A ǫ = n i=1 σ(A) 2 i σ(A) 2 i + ǫ . (5.4)
Note that the problem (5.3) is differentiable as long as . ǫ is differentiable. The gradient can be computed efficiently using the Jacobian matrix of the singular values as expressed in [START_REF] Papadopoulo | Estimating the jacobian of the singular value decomposition: Theory and applications[END_REF]. In practice, we solve this optimization problem using a standard L-BFGS algorithm. The choice of ǫ can have a big impact on the results. In fact since we are using a gradient descent method the big singular values are in the flat part of . ǫ therefore their gradient will be granted a small weight. On the contrary the singular values in the slope will have a big influence on the minimization. So with an ǫ too small only a few singular values will be minimized but with an ǫ too large many singular values will not be well minimized. We chose the parameter ǫ such that at the initialization 80 percent of the singular values satisfy

σ 2 i σ 2
i +ǫ ≤ 0.9.

Basis function extraction

Since the probe functions can give redundant information in some shape parts and incomplete information in others, the resulting functional map will map some subspaces of L 2 (M 0 ) with more confidence than others. Using a collection of shapes we would like to extract the most stable subspaces.

For this purpose we propose to use the learned optimal weights D and the resulting estimated functional maps X i (D) and to identify stably mapped functional subspaces by comparing X i (D) to the reference maps C i . The input here is the same as in the previous section. We need n shapes with known ground truth functional maps C i and a set of consistent probe functions G i . The output will be Y an orthonormal basis of L 2 (M 0 ) ordered with decreasing confidence. As we demonstrate in Section 5.5, in most cases this order remains stable even for maps that are estimated to previously unseen shapes.

Identifying stable subspaces

The best mapped function

y 0 ∈ L 2 (M 0 ) is such that X i y 0 is the closest to C i y 0 for all i.
Such function is solution of the problem:

y 0 ∈ arg min y∈L 2 (M 0 ), y =1 n i=1 (X i -C i )y 2 F
We can then iteratively define an orthonormal basis of L 2 (M 0 ) ordered by decreasing accuracy in the mapping, by solving the following problem:

y n+1 ∈ arg min y∈L 2 (M 0 ), y,y j =0 ∀j≤n n i=1 (X i -C i )y 2 F
Such a basis can be efficiently computed by considering the singular value decomposition of the matrix:

B =   X 1 -C 1 . . . X n -C n   = U ΣV t .
It is well-known that y j must be equal to singular vectors corresponding to the j th smallest singular value of B. We can, therefore, form a new orthonormal basis Y of L 2 (M 0 ) composed of the singular vectors of B by increasing singular values. This allows us to quantify the quality of the mapping of a functional subspace by looking at the singular values of B: the smaller the singular values are, the better the mapping.

Functional map to a test shape using a reduced basis

Now if we are given a new unseen shape M n+1 that does not belong to the training set, we first compute its probe functions and store them in a matrix G n+1 . We then compute the functional map X n+1 by using the previously solved for weight matrix D. Finally, since we know that X n+1 contains some badly mapped subspaces (for example the antisymmetric functions), by using Y p the p first column of Y , we compute the reduced map X p n+1

X p n+1 = X n+1 Y p : L 2 (M 0 ) ∩ L 2 (Im(Y p )) → L 2 (M n+1 ).

Experimental Results

Functional correspondences

The probe functions used to solve the problem in Eq. ( 5.2) are given by various descriptors computed on each shape:

• Heat Kernel Signature [START_REF] Sun | A concise and provably informative multi-scale signature based on heat diffusion[END_REF] for multiple values of t

• Wave Kernel Signature [START_REF] Aubry | The wave kernel signature: A quantum mechanical approach to shape analysis[END_REF] at three different energies

• Gaussian and Mean Curvature

• Logarithm of the absolute value of Gaussian and Mean Curvature

• Mesh Saliency [START_REF] Lee | Mesh saliency[END_REF] The HKS, WKS and Mesh Saliency are computed at various scales to ensure a wide variability. We process the curvature functions to obtain a family of descriptors. Since the curvatures can have very high peaks we take the logarithm of their absolute value to put more weight on the small curvatures areas. The family of functions is then created by considering the solution of the Heat Diffusion Equation at various times when each function is used as initial heat distribution over the surface.

Isometric Shape Matching

TOSCA Dataset. Our method is a tool which can be used in addition to other methods using the functional framework in order to improve the approximation of the functional maps. As several methods using functional maps [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF][START_REF] Ovsjanikov | Shape matching via quotient spaces[END_REF] have been shown to be more efficient than the state-of-the-art methods, we compare our trained maps to the baseline "original" method described in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF].

We have evaluated our method on the shape matching benchmark TOSCA [START_REF] Bronstein | Numerical geometry of non-rigid shapes[END_REF]. For each shape class we use all the available shapes for training, except one for testing and we choose the standard undeformed pose as shape M 0 . We compare three ways of weighting the probe functions: a single weight for all the functions, a weight per category of descriptors and one weight per probe function.

For all the experiments we express all functions in the basis given by the first 50 eigenfunctions of the Laplace-Beltrami operator. We compute 50 probe functions divided in 9 categories (WKS is divided in three categories with tree different energies) of descriptors. We take 5 functions per category except for the Mesh Saliency where 10 functions are computed. Since all the shapes in TOSCA have an internal symmetry, we cannot hope to recover the entire functional map, and thus Eq. (5.4) is a reasonable choice of norm.

The experiments follow the pipeline shown in Figure 5.2. First we learn the optimal weights and extract the ordered basis using the training set of shapes. Second we are given an unknown shape. We use the optimal weights to compute the functional map and the extracted basis to suppress the badly mapped function subspaces. The L-BFGS algorithm used to solve the optimization problem in Eq. ( 5.3) is initialized with the naive functional maps solution of (5.1) with α = 10 -3 . We compare all of the functional maps and subspaces computed with our method to the baseline "naive" map, obtained using the identity matrix D, which correspond to the original method described in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF].

Performance The proposed approach was implemented in MATLAB. Note that the number of vertices of each shape has no effect on the performance since all the functions are expressed in a reduced functional basis. The most time consuming task in our pipeline is the training part which requires to solve a difficult non-linear optimization problem (5.3). The processing cost is dominated by the computation of the gradient of the energy, which is done by solving two linear systems for each shape of the training set at each iteration. However, the contributions of each shape to the gradient are independent so this can be done in parallel. The learning process with a training set of 10 shapes took about 45 min on an Intel i7 processor without parallelization.

Optimal Weights Figure 5.3 (left) shows the weights obtained after solving the problem in Eq. ( 5.3) with a training set composed of 9 cats. To demonstrate the importance of weighting the probe functions on the quality of the functional map, we study the distribution of the singular values of the difference X n+1 (D) -C n+1 for the different learned weights. In Figure 5.3 (right) each curve depicts the percentage of singular values below the threshold given on the x-axis. For the perfect map, all singular values would be zero. As can be seen, the functional maps with the optimal weights have a bigger concentration of small singular values than the naive functional map. Therefore there exists a bigger functional subspace on which these functional maps provide a good approximation of the ground truth. Note that the naive map has no small singular values and is indeed a very bad approximation.

Stable subspaces From the naive maps and functional maps with optimal weights, we extract four function bases ordered by decreasing stability. The most stable functions for each case are shown in Figure 5.4. Even the most stable functions from the naive maps are not mapped very accurately since they are very bad approximation of the groundtruth. For the other bases the functions seem consistent with the information we would expect from descriptors as HKS and WKS: a distinction between flat area (body) and salient area (legs, tail, head). Note that even with only one weight we are able to retrieve meaningful stable areas. We also evaluate the extracted functional basis by computing the difference between the ground truth map and our approximation on the unseen shape:

ǫ i = (X n+1 (D) -C n+1 )y i 2 , (5.5)
where y i is the i th function of the extracted basis. We compare this error for the three weighting strategies with the naive map in Figure 5.5 (left). The extracted basis was ordered by decreasing quality on the training set. Note that this order is still preserved on the unknown shape for the 25 first functions of the basis. Most of all we are able to identify the worst mapped subspaces, which can be safely removed. Despite only estimating the functional maps on a subspace of the full functional space, we converted them to point-to-point correspondences using the method described in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF]. Figure 5.5 right compares the quality of the point-to-point correspondences before and after reducing the space dimension from 50 to 25. We obtain better results with our learned weight than with the naive map. For the weighted maps, the reductions perform better or similarly than the full maps. Thus our basis extraction manages to identify correctly the most stable subspaces. For the naive map our reduction space strategy fails as there is no well-mapped subspace.

Non-Isometric Shape Matching

Until now we have assumed the deformation to be nearly isometric. Our algorithm to find the optimal set of probe functions and the extraction of the most stable subspaces do not contain any explicit knowledge of the type of deformation. In fact, this assumption is only used to construct the least-squares problem (5.1). Which means that our framework can be adapted to any kind of deformation model as long as we have a consistent way of computing a functional map form probe functions and a reduced basis that is approximately stable.

To test this case, we consider a man or a gorilla and 12 women in different poses from the TOSCA dataset. The ground truth functional maps are computed using a thousand user-picked correspondences. The meshes have different number of vertices and different connectivity. We train our method on 10 poses and use the last for testing.

The resulting problem is very noisy for two reasons. First, the ground truth functional maps are computed from sparse correspondences, and therefore can be inaccurate on some functional subspace. However, the use of a collection of shape for training allows us to remove this noise. Second, since all the probe functions used are designed for isometric deformation, few are going to contain useful information.

In order to introduce a wide variability of functions we pre-compute on each shape 50 basis functions, 310 probe functions and we put a weight on each probe function. The optimization algorithm is initialized as in the previous experiment. Figure 5.7 (left) shows the most stable functions learned from the training maps. Each of these functions is also mapped to a previously unseen pose using the "ground truth" map converted to a functional map. Note that the functions are badly transferred, due to the incompatibility of the LB basis and the noise in the input maps. Compare this with Figure 5.7 (right) where the probe functions have been weighted using our method. The stable functions indicate the head, the hands and the feet to be the most stable area. Besides, these functions are correctly mapped on a new shape using a computed functional map with the learned weights. Clearly, the fact that we use a collection helps removing the noise of the input data. Thus, our method is able to correctly identify the most stable functional subspace under mild assumptions on the underlying deformation.

Conclusion

In this work, we presented a method to learn the most informative descriptors for non-rigid shape matching, from a given set of shape correspondences. Our method is purely intrinsic and allows us to obtain high quality consistent correspondences to new, unseen shapes, and to identify the most reliably mapped functional subspaces. The approach is flexible and can potentially be applied to scenarios that lack a good theoretical deformation model, as demonstrated by meaningful even non isometric deformation. One of its weaknesses is a relatively high cost for the training due to the non-convex nature of the energy. In the future, we plan to explore more efficient optimization strategies.

In this chapter we proposed an efficient method to compute functional maps. However in many practical scenarios, one would like to obtain pointwise correspondences between the given shapes. Thus, in the following chapter we present a method for converting functional maps to point-to-point correspondences, that enforces continuity in the resulting map.

Stable functions before training. Stable functions after training. Chapter 6

Continuous Matching via Vector Field Flow

In this chapter we present a new method for non-rigid shape matching designed to enforce continuity of the resulting correspondence. Our method is based on the functional map representation, which allows efficient manipulation and inference but often fails to provide a continuous point-to-point mapping. We address this problem by exploiting the connection between the operator representation of mappings and flows of vector fields. In particular, starting from an arbitrary continuous map between two surfaces we find an optimal flow that makes the final correspondence operator as close as possible to the initial functional map. Our method also helps to address the symmetric ambiguity problem inherent in many intrinsic correspondence methods when matching symmetric shapes. We provide practical and theoretical results showing that our method can be used to obtain an orientation preserving or reversing map starting from a functional map that represents the mixture of the two. We also show how this method can be used to improve the quality of maps produced by existing shape matching methods, and compare the resulting map's continuity with results obtained by other operator-based techniques.

Introduction

Computing correspondences or mappings between 3D shapes is one of the key building blocks in many areas of digital geometry processing, including deformation transfer [START_REF] Sumner | Mesh-based inverse kinematics[END_REF], shape interpolation (morphing) [START_REF] Kilian | Geometric modeling in shape space[END_REF] and statistical shape analysis [START_REF] Hasler | A statistical model of human pose and body shape[END_REF] among many others. This problem is particularly challenging in the case of shapes undergoing non-rigid deformations, where the notion of the optimal map may be difficult to define and optimize for. As mentioned in the previous chapter, most of the successful global methods proposed to find correspondences between pairs of non-rigid shapes in the recent years have relied on a variant of the conformal [START_REF] Wang | Conformal geometry and its applications on 3d shape matching, recognition, and stitching, Pattern Analysis and Machine Intelligence[END_REF][START_REF] Lipman | Mobius voting for surface correspondence[END_REF][START_REF] Kim | Blended intrinsic maps[END_REF] or fully isometric [START_REF] Bronstein | Generalized multidimensional scaling: a framework for isometry-invariant partial surface matching[END_REF][START_REF] Tevs | Isometric registration of ambiguous and partial data[END_REF][START_REF] Sahillioğlu | Coarse-to-fine combinatorial matching for dense isometric shape correspondence[END_REF][START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] deformation models, which assume that either the angles or the geodesic distances between pairs of points are approximately preserved by the mapping. Although such models have very appealing theoretical properties, using them directly can often lead to difficult non-linear, non-convex optimization problems [START_REF] Bronstein | Generalized multidimensional scaling: a framework for isometry-invariant partial surface matching[END_REF]. Therefore, most recent works in this direction have concentrated on finding a low-dimensional parameterization of the space of mappings, that allows for efficient optimization techniques (e.g. [START_REF] Lipman | Mobius voting for surface correspondence[END_REF][START_REF] Ovsjanikov | One point isometric matching with the heat kernel[END_REF][START_REF] Brunton | A lowdimensional representation for robust partial isometric correspondences computation[END_REF]). Among such low-dimensional representations of the space of correspondences, one particularly appealing approach is based on the framework of functional maps [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF], which consider mappings as linear operators between the corresponding function spaces. One of the weaknesses of the functional map representation, however, is that by representing mappings as correspondences between functions, it requires an additional post-processing step to obtain a point-to-point map after computing the optimal functional map. The basic approach for this conversion step, proposed in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] and used in most follow-up works, assigns points by considering the mapping between the corresponding Dirac deltafunctions. Since each delta-function is mapped independently, however, this approach can (and most often does) introduce significant artifacts and discontinuities into the final point-to-point mapping (see the first two columns of Figure 6.6). This makes the resulting correspondences unusable in settings that require continuity of the mapping, such as texture transfer. Additional pair-wise terms can potentially be introduced in the conversion procedure, but this would require creating variables for points with potentially very expensive consistency constraints, which very quickly loses the appeal of the functional map framework, and reduces to direct optimization. In this context, we propose a novel method for converting a functional map to a point-topoint map, which guarantees continuity and does not rely on any pairwise consistency constraints, making it computationally efficient. Our main idea is to represent the target point-to-point map as a composition of an arbitrary continuous map between the two surfaces and a flow associated with an unknown vector field on one of them. By relying on the recently proposed operator representation of vector fields [START_REF] Azencot | An operator approach to tangent vector field processing[END_REF], we show that the optimal vector field can be computed efficiently entirely within the functional map framework, and the computation of the final map requires a single discretization of vector field advection. We also employ the supervised learning technique presented in Chapter 5 that not only helps to obtain better functional maps but also helps to identify functional subspaces where the map is reliable, which significantly helps to improve the final point-to-point map.

Our method also helps to address the symmetric ambiguity problem inherent in many intrinsic correspondence methods when matching symmetric shapes. We provide practical and theoretical results showing that our method can be used to obtain an orientation preserving or reversing map starting from a functional map that represents the mixture of the two. Finally, we test our method on a shape collection and show that we can produce maps that are both continuous and have smaller geodesic distortion compared to the results obtained by existing techniques.

Related Work

Below we concentrate on the recent works that are directly related to ours, consisting of methods for global near-isometric shape matching with special emphasis on approaches that guarantee the continuity of the resulting maps. As mentioned in the introduction, most of the existing techniques for non-rigid shape matching use a deformation model for finding correspondences between 3D shapes. The two most common models in this setting include approximate intrinsic isometries and conformal mappings. The former model, which was originally introduced by Bronstein et al. [START_REF] Bronstein | Generalized multidimensional scaling: a framework for isometry-invariant partial surface matching[END_REF] and Mémoli [START_REF] Mémoli | On the use of Gromov-Hausdorff Distances for Shape Comparison[END_REF] assumes that pairwise geodesic distances are approximately preserved by the deformation. The first works that use this assumption lead to continuous maps by design, but result in very challenging optimization problems that are difficult to solve with more than a small number of points [START_REF] Bronstein | Generalized multidimensional scaling: a framework for isometry-invariant partial surface matching[END_REF]. As a result, many follow-up techniques have used a relaxed version of the isometric mapping assumption, which result in more manageable optimization problems, but can often fail to guarantee a low distortion continuous mapping (e.g., [START_REF] Huang | Non-rigid registration under isometric deformations[END_REF][START_REF] Tevs | Isometric registration of ambiguous and partial data[END_REF][START_REF] Ovsjanikov | One point isometric matching with the heat kernel[END_REF][START_REF] Sahillioğlu | Coarse-to-fine combinatorial matching for dense isometric shape correspondence[END_REF][START_REF] Brunton | A lowdimensional representation for robust partial isometric correspondences computation[END_REF]). Furthermore, an additional challenge in using the isometric model assumption is that exact intrinsic isometries are extremely rare, both in theory [START_REF] Gluck | Almost all simply connected closed surfaces are rigid[END_REF] and in practice, since most deformable shapes induce some amount of distortion.

Another set of successful techniques, which are more widely applicable than those based on the isometric mapping assumption are those that assume that the mapping is conformal, and thus only preserves angles (e.g., [START_REF] Haker | Conformal surface parameterization for texture mapping[END_REF][START_REF] Wang | Conformal geometry and its applications on 3d shape matching, recognition, and stitching, Pattern Analysis and Machine Intelligence[END_REF][START_REF] Hurdal | Discrete conformal methods for cortical brain flattening[END_REF][START_REF] Lipman | Mobius voting for surface correspondence[END_REF][START_REF] Kim | Blended intrinsic maps[END_REF]). These techniques are appealing because a conformal mapping is known to exist between any pair of shapes with the same topology, but also because the set of such mappings can be parameterized relatively easily by using a canonical domain, such as a sphere for genus zero surfaces. Moreover, the resulting maps obtained by these approaches are typically continuous. At the same time, conformal mappings can often induce large area distortion, which can result in unrealistic correspondences between non-rigid shapes, which limits their use significantly.

A recent set of approaches that overcome the above-mentioned challenges to some extent is based on the functional map representation, introduced in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF]. As mentioned above, this framework is based on representing maps as linear operators acting on real-valued functions, and which can be encoded compactly by small-sized matrices in the discrete setting by using a multi-scale basis. Although the original approach and the follow-up works, including [START_REF] Kovnatsky | Coupled quasi-harmonic bases[END_REF][START_REF] Pokrass | Sparse modeling of intrinsic correspondences[END_REF], all implicitly use the isometric deformation assumption, they have been shown to be very robust to small non-isometric distortions, by extensive use of strong geometric and linear-algebraic regularization techniques. Moreover, several recent works have shown how this framework can be used in the supervised learning setting, where functional maps between unseen shapes can be obtained by exploiting information present in a small set of example maps including [START_REF] Rodolà | Dense non-rigid shape correspondence using random forests[END_REF] and Chapter 5.

Despite its practical appeal, one of the limitations of the functional map framework, is that a post-processing step is necessary to convert a functional map to a point-to-point one. The method used in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] is based on mapping Dirac delta functions. However as the points are considered independently the continuity of the resulting map is not ensured. This problem can be particularly prominent in shapes that contain intrinsic symmetries, which contain at least two equally good solutions for the optimal functional map, and the computed one is at best a linear blending of the two.

Note that, closely related to our technique, especially in the use of flows for computing continuous maps (diffeomorphisms) is the LDDMM framework [START_REF] Beg | Computing large deformation metric mappings via geodesic flows of diffeomorphisms[END_REF][START_REF] Miller | Geodesic shooting for computational anatomy[END_REF], widely used in the medical imaging community. Unlike these methods, however, our approach is purely intrinsic and operates directly on the surface of the target model, rather than deforming a template in space.

Contributions

In this chapter we propose a novel method for converting a functional map into a point-to-point one, which combines the strengths of the functional map framework that allows to compute low-distortion functional maps, with those of the conformal mapping approaches, which produce continuous correspondences. Namely, starting from a map computed using the state-of-the-art conformal-based Blended Intrinsic Map approach [START_REF] Kim | Blended intrinsic maps[END_REF], we modify it by computing the optimal vector field, whose flow, composed with the original map, would result in a functional map as close as possible to the given one. By using the recently proposed operator representation of vector fields [START_REF] Azencot | An operator approach to tangent vector field processing[END_REF] and the connection between advection and matrix exponentiation, we propose an efficient optimization approach for computing the optimal vector field entirely within the functional map framework. Moreover, we show theoretically that this approach is guaranteed to produce the correct continuous map when the input functional map represents a blending of the orientation preserving and reversing maps under certain assumptions, and demonstrate this projection step in practice.

Functional Maps Conversion

Once the functional map C is computed following the pipeline proposed in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] or by Chapter 5, the goal of the second step is to convert it to a point-to-point map. The method proposed in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF] and reused in most of the follow-up work consists in finding the nearest neighbors of the images of Dirac-delta functions on M by C among the Dirac functions on N . Namely, for each point p ∈ M , the map: ϕ(p) is computed as via ϕ(p) = arg min q ||δ q -Cδ p ||, where δ p is an indicator function on point p, written in the appropriate basis.

Main Challenges

While both steps described above are very efficient in practice, the second stage has a very serious limitation, in that it processes each point independently, meaning that the final map ϕ may not be (and often is not) continuous. The first two columns of Figure 6.6 provide examples of discontinuous maps resulting from this conversion. To illustrate this phenomenon, let us assume that the target shape N has an orientationreversing (reflectional) intrinsic symmetry S : N → N. In this case, there exist at least two equally good potential solutions for Eq. (5.1) and similarly, each point x may have several candidate correspondences. In practice the functional constraints are often not sufficient to resolve symmetric ambiguities, in large part because most robust descriptors are invariant under intrinsic isometries. The best we can hope for when approximating C ϕ is an exact functional map for symmetric functions (i.e. f , s.t. f • S = f ) and a noisy or zero functional map for antisymmetric functions (i.e. f • S = -f ). Since our approximations are obtained by solving a linear system, most likely a solution of the least squares problem will be a linear blending between the orientation preserving and reversing functional map:

C α ϕ = (1 -α)C ϕ + αC ϕ•S (6.1) N M M ϕ C ϕ 0 φ t V L 2 (N ) L 2 (M ) L 2 (M ) C C ϕ 0 exp(tD V )
Figure 6.1 -Left: the unknown continuous map ϕ C is a composition of the input ϕ 0 and the flow φ t V of a vector field V . Right: dual representation as functional maps.

Note that α = 0.5 implies that all antisymmetric functions are mapped to zero.

The conversion of C to a point-to-point map in itself gives no guaranty of continuity in the resulting map. Since each Dirac function of a point x is treated independently it can be mapped indifferently to its image ϕ -1 (x) or to its symmetric alternative S(ϕ -1 (x)). Moreover this process is not designed to be stable under the blending noise α, as in (Eq. 6.1).

In this context, the key idea developed in this chapter is to construct a point-to-point map from the functional map C by following a procedure that guarantees continuity, while being robust to blending noise. In particular, starting from an arbitrary continuous map between M and N , we find an optimal vector field, whose flow makes the final correspondence operator as close as possible to the given (e.g., computed) functional map.

Since the flow of a vector field provides a continuous, and orientation-preserving map, the final correspondence is both continuous and has the orientation of the initial map. As we show below, this can significantly improve the quality of the resulting point-to-point map, while remaining computationally tractable and avoiding expensive second-order pairwise constraints.

Algorithm Overview

The algorithm proposed in this chapter takes as input a functional map C : L 2 (M ) → L 2 (N ) and an arbitrary continuous map ϕ 0 : N → M . It then outputs a continuous point-to-point map ϕ C : N → M .

As mentioned above, the main idea of our algorithm is to construct the map ϕ C by composing ϕ 0 with the flow φ t V of a well-chosen vector field V (see Figure 6.1). We will choose the vector field V such that φ t V • ϕ 0 represented as a functional map is as close as possible to the input C. This can be done efficiently by representing φ t V as an operator (Section 6.4) and then solving a small-scale optimization problem as explained in Section 6.5. To find the map ϕ C we solve a system of ODEs with a simple solver (Section 6.6).

The main steps of the proposed algorithm are described in Algorithm 1.

Algorithm 1: Functional Map Conversion

Input : C : L 2 (M ) → L 2 (N ) functional map ϕ 0 : N → M initial continuous map Output : ϕ C : C converted into a continuous map 1 Find Optimal Vector field (Section 6.5); 2 Convert ϕ 0 to a functional map C ϕ 0 ; 3 Solve: a ⋆ ∈ arg min

a∈R n C ϕ 0 exp ( n i=1 a i D V i ) -C φ ; 4 Set: V := n i=1 a ⋆ i D V i ; 5 Compute ϕ C (Section 6.6); 6 Solve: d dt φ t V (p) = V φ t V p)), φ 0 V (p) = p ∈ N ; 7 return ϕ C := φ 1 V • ϕ 0 ;

Family of Diffeomorphisms

In this section we construct a family of diffeomorphisms which map N onto M and derive their representation as functional maps. The point-to-point map which converts the given functional map C will be chosen among this family.

Vector field flow Given a family of tangent vector fields {V i } 1≤i≤n on M , we let V be the space spanned by the linear combinations of the V i . Any vector field V ∈ V, defines a one-parameter family of maps φ t V : M → M called the flow of V . The flow is formally defined as the unique solution to the differential equation:

d dt φ t V (p) = V φ t V p)), φ 0 V (p) = p ∈ N. (6.2)
Given an arbitrary diffeomorphism ϕ : N → M we construct a family of diffeomorphisms T parametrized by t ∈ R and a ∈ R n :

ϕ t a (p) = φ t Va • ϕ(p), V a = n i=0 a i V i (6.3)
Remark that the orientation of a map ϕ t a ∈ T is given by the orientation of ϕ since the flow of a vector field is orientation preserving.

Functional Representation of the family The family of mappings T has an easy representation in the functional map framework as explained in [START_REF] Azencot | An operator approach to tangent vector field processing[END_REF]. This is because, a vector field V on a smooth manifold can be represented as an operator D V acting on a function f :

D V (f )(p) = V p , ∇f (p) p . (6.4)
Since the action of D V is linear, the operator is conveniently represented as a matrix in the discrete setting.

It is well known that g t = f • φ t V is the unique solution of the PDE:

∂g ∂t (t, p) = D V (g)(t, p), g(0, p) = f (p).
A key property of the operator representation of vector fields, introduced in [START_REF] Azencot | An operator approach to tangent vector field processing[END_REF] is that for analytic functions the functional map C φ t V is represented by the exponential of the operator D V since one has:

C φ t V f := f • φ t V = exp(tD V )(f )
Since map composition is achieved via matrix multiplication in the functional representation, this yields a simple way of describing our family of diffeomorphisms T . Let ϕ t a ∈ T then

C ϕ t a = C ϕ exp t n i=1 a i D V i . (6.5)
6.5 Optimal vector field

Optimization Problem

Our main idea, developed in the section, is to project the input functional map C onto the appropriate set of diffeomorphisms T . Namely our goal is to find a vector field V ∈ V such that the operator representation (6.5) of ϕ t a is as-close-as possible to C. This projection is easily written thanks to the operator representation, and computationally it reduces to solving the optimization problem:

min a∈R n C ϕ exp n i=1 a i D V i -C φ , (6.6)
for an appropriate choice norm . φ . Here we note briefly that the norm is chosen to be differentiable as indicated in Section 5.3.2.

In practice, the problem (6.6) can be solved using a first order method such as the L-BFGS algorithm. The main difficulty in finding the gradient of the objective function lies in the computation of derivative of exp ( n i=1 a i D V i ) in the direction V j . While there exists a vast literature on approximating the exponential of a matrix (for a survey see [START_REF] Moler | Nineteen dubious ways to compute the exponential of a matrix, twenty-five years later[END_REF]), to the best of our knowledge few methods address the problem of computing the directional derivative of the matrix exponential, which is conceptually non-trivial. As we show in Lemma 6.5.1, however, the directional derivative can be obtained as a block of the matrix exponential of a bigger operator. Note that if there are n vectors in the family of vector fields we have to compute n matrix exponentials. Lemma 6.5.1 The directional derivative d H e tA defined as

d H e tA = lim h→0 1 h (exp (t(A + hH)) -exp (tA)) ,
is one block in the matrix exponential of a bigger operator:

e tA B t d H e tA C t = exp tA 0 tH tA
Proof Let x 0 be an arbitrary vector in R n and x(t) = exp (tA) x 0 , it is well known that x(t) satisfies the ODE:

x ′ (t) = Ax(t), x(0) = x 0 . Moreover x h (t) = exp (t(A + hH)) x 0 is solution of x ′ h (t) = (A + hH)x h (t), x h (0) = x 0 .
We denote y(t) the directional derivative in the direction H:

y(t) := lim h→0 1 h (exp (t(A + hH)) -exp (tA)) x 0 := lim h→0 1 h (x h (t) -x(t)) .
After computing y ′ (t) we conclude that y(t) is the unique solution of the ODE

x ′ (t) = Ax(t), x(0) = x 0 y ′ (t) = Ay(t) + Hx(t), y(0) = 0
Expressing the solution as an exponential of a matrix leads to the results.

Properties

One of the advantages of the formulation of the problem of finding the optimal point-topoint map from a functional map via Eq. (6.6) is that it makes no assumptions on the input map C. This is particularly important since, as mentioned above, in the presence of intrinsic symmetries the functional map C can, even in the best case, be a linear blending of the functional representation of an orientation-preserving and orientation-reversing map. However, one potential problem is that the presence of the "noisy" part in the functional map can adversely affect the final output map ϕ obtained by optimizing Eq. (6.6). Fortunately, both in theory and in practice this is not the case. Namely, under some suitable assumptions, the orientation-preserving ground-truth functional map, must be a local minimum of the problem (6.6) even when the functional map C is given by the symmetry blending defined at (6.1). In particular, as we show in the Lemma 6.5.2, C ϕ must be a local minimum of Eq. (6.6) under some assumptions. Lemma 6.5.2 If the norm . φ = . 2 F is the squared Frobenius norm, the set of vector fields considered V is divergence-free and the initial transformation ϕ approximately isometric, then C ϕ must be a local minimum of Eq. (6.6).

Proof The first order necessary condition for C ϕ to be a local minimum reads:

∀dC ϕ ∈ T (C ϕ ), dC ϕ , C ϕ -C α ϕ F = 0,
where T (C ϕ ) is the tangent space of the set of all functional map at the point C ϕ . This tangent space has a simple expression. Let's consider a small perturbation of C ϕ by the flow φ t V induced by the vector field V ∈ V applied to an arbitrary function f :

lim t→0 1 t C φ t V •ϕ f -C ϕ f = lim t→0 1 t f • φ t V • ϕ -f • ϕ = d dt (f • φ t V • ϕ) | t=0 = V, ∇f • ϕ = C ϕ (D V (f )).
Additionally the deformation is nearly isometric so C ⊤ ϕ ≈ C ϕ -1 see Theorem 4.1.4. The necessary condition becomes:

∀V ∈ V, (1 -α) D V , I -C ϕ•S•ϕ -1 F ≈ 0.
Remark that the mapping ϕ • S • ϕ -1 : M → M is an intrinsic symmetry on M . Now suppose that the basis function is composed only by even and odd functions with respect to the symmetry S. Therefore the functional map associated to an internal symmetry is a diagonal matrix with 1 and -1 on the diagonal corresponding to the symmetric and antisymmetric eigenfunctions. So I -C ⊤ ϕ C S C ϕ is approximatively a diagonal matrix. By assumption the vector field V is divergence free therefore represented by a skewsymmetric operator D V as explained in [START_REF] Azencot | An operator approach to tangent vector field processing[END_REF]. Since the result does not depend on the basis C ϕ is a critical point of (6.6).

Practical Choice of the Norm

As stated before C is not reliable for antisymmetric functions. Therefore there is some function subspace on which C and C ϕ exp (t n i=1 a i D V i ) cannot agree. The choice of the norm . φ in the problem (6.6) is of critical importance. Similarly to setting in Chapter 5 above, the naive choice of the squared Frobenius norm is not well-suited for this problem since it is the sum of the squared singular values. As such, it will give a large weight on badly matched function subspace and a small weight on well matched function subspace. However, since typically we have almost no information about antisymmetric functions so the optimization problem based on this problem will put a lot of effort matching functions that we cannot hope to match and few matching interesting subspace. A better choice for . φ is a regularization of the nuclear norm. We choose

A φ = AY p ǫ,⋆ (6.7)
where Y p is a basis of p functions that we want to focus on obtained using the approach described in Chapter 5 and based on Chapter 5. In the unsupervised setting we chose Y p to be the identity. The norm . ǫ,⋆ is a defined by n i=1 σ 2 i + ǫ the σ i are the singular values of the matrix. With this norm, we give smaller weight to the subspaces that are difficult to align and focus on task we are able to complete. The parameter ǫ makes the function . ǫ,⋆ differentiable and should be small and is taken at 10 -3 . Note that the Jacobian matrix of the singular values is easily computable as explained in [START_REF] Papadopoulo | Estimating the jacobian of the singular value decomposition: Theory and applications[END_REF].

Vector Field Flow on Manifold

Once the optimal vector field V is found using the procedure described above, we obtain the final point-to-point map by composing the initial map ϕ with the flow of V . To compute this flow, we need to solve the system of equations (6.2) on the given triangle mesh. In principle any advection solver will work with our method. However since computing the flow is known to be potentially difficult, we implemented our own solution. The implementation we use gives a coarse approximation of the flow and might not be accurate for very large deformations. For more accurate solution of this problem we refer to [START_REF] Ray | Robust polylines tracing for n-symmetry direction field on triangulated surfaces[END_REF][START_REF] Myles | Robust field-aligned global parametrization[END_REF] which provide more guaranties of continuity of the flow and faster convergence. In all of our applications we assume that the shapes are given as triangulated meshes and the vector field is given as a single vector per face. Given this representation, we assume that that the vector field is constant per face and is interpolated at the edges. Three main situations can occur: the current point could be at inside face, an edge or a vertex.

At a Face Since inside a face the vector field is assumed to be constant, we follow it until we reach an edge or a vertex (Figure 6.2 from p 0 to p 1 ).

At an Edge When a point is at an edge, we try to cross the face we did not come from (Figure 6.2 from p 1 to p 2 ). If the point did not move we follow the edge by interpolating the vector field from the two neighboring faces and end up at a vertex (Figure 6.2 from p 2 to p 3 ).

At a Vertex When the point is at a vertex, we try to follow the vector of each of the neighboring faces and choose the one that goes the furthest. If the point cannot move, we try to follow the neighboring edges, using interpolated directions and to potentially end up at another vertex (Figure 6.2 from p 3 to p 4 ).

Results

For all the experiments we express all functions in the basis given by the first 150 eigenfunctions of the Laplace-Beltrami operator. We choose a family of 50 tangent vector fields for the V i given by the first eigenfunctions of the 1-form Laplace-de Rham operator, constructed following the procedure described in [START_REF] Fisher | Design of tangent vector fields[END_REF]. We have evaluated our method for computing point-to-point correspondences on the shapes on the benchmarks of Anguelov and al. [START_REF] Anguelov | Scape: shape completion and animation of people[END_REF] and of Bronstein and al. [START_REF] Bronstein | Numerical geometry of non-rigid shapes[END_REF]. In all of the cases, the input continuous map ϕ 0 is the result of the BIM algorithm [START_REF] Kim | Blended intrinsic maps[END_REF]. This map is most of the time continuous but can be very distorted in some areas. We will show that our method is able to detect the distorted areas and correct them.

Symmetry Blending

As stated above a plausible perturbation for the input functional map C is given by equation (6.1). We test our method when C is the linear blending of the ground-truth functional map and the ground-truth orientation reversing functional map for various values of α. In this experiment Y p , in Equation 6.7, is the identity matrix. For this experiment we choose a pair of shapes from the SCAPE dataset. The graph shown in Figure 6.3 shows the percentage of correspondences with a geodesic error smaller than a threshold. Of course the closer C α ϕ is to the ground-truth map the better are the correspondences. However our results are robust even when the target functional map is an exact blending of the direct and symmetric map and are always better than the map coming from BIM. Thus, even when the assumptions of our theoretical observation are not fulfilled, our method can successfully retrieve meaningful information from noisy data.

Error using a computed functional map

In a more realistic scenario, rather than using a ground-truth functional map, we compute it via the inference pipeline described in Chapter 5. In this section the experiments are conducted on several pairs of shapes: 7 human pairs (SCAPE) and 5 animal pairs (TOSCA). The functional map C is computed using the least squares problem (5.1), where each functional constraints is weighted. The weights are learned by solving problem (5.3) using the algorithm described in Chapter 5 which also outputs a matrix Y p corresponding to the p best mapped functions, where we let p equal to 70. The training set is composed of 8 randomly chosen meshes for the SCAPE example and 4 meshes for the TOSCA centaur example. We compute 310 functional constraints equally distributed among these categories:

• Heat Kernel Signature [START_REF] Sun | A concise and provably informative multi-scale signature based on heat diffusion[END_REF],

• Wave Kernel Signature [START_REF] Aubry | The wave kernel signature: A quantum mechanical approach to shape analysis[END_REF] at three different variances,

• Gaussian and Mean Curvature,

• Logarithm of the absolute value of Gaussian and Mean Curvature,

• Mesh Saliency [START_REF] Lee | Mesh saliency[END_REF].

We compare our approach with BIM, that serves as ϕ 0 , and with the functional map C converted to point-to-point map using the method proposed in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF]. The graph in Figure 6.4 shows the percent of correspondences which have geodesic error smaller than a threshold in average for SCAPE and TOSCA. In this case, we only accept direct correspondences as correct, and consider symmetric points as wrong. Note that our method shows quality improvement over Blended Intrinsic Maps. The direct conversion of C have some point with very large geodesic error due to points mapped to their symmetric counterparts. We evaluate the continuity of our map with two measures of distortion. First the maximum radius corresponding to a geodesic ball of given size. For a map ϕ this is formally given by the function: where d N is the geodesic distance on N . If the map is nearly isometric r should be close to identity. We compare this measure for different mapping in Figure 6.5a for one example from TOSCA. Our method is comparable to BIM and to the ground truth in terms of continuity while the direct conversion of C show some very large distortions. Second we compare the ratio between the triangle's area before and after deformation. Since the deformations in our examples are almost isometric this ratio should be close to one. The graph in Figure 6.5b shows the percent of triangles which have an area ration smaller than a threshold. We show only the ratio greater than one since most of the discontinuous behavior is due to large jumps. The area ratio of the exact mapping are concentrated around one which is consistent with the fact that the deformation is nearly isometric. Again the direct conversion of C show some very large area distortions compare to BIM and our method. This lack of continuity is confirmed by Figure 6.6 which provides on two examples a visualization of the point-to-point mapping using color correspondence. The direct conversion of the functional map shows some artifacts due to the blending between orientation preserving and orientation reversing maps. Our method successfully repairs the areas distorted by BIM as shown on Figure 6.7 for two different matching problems. In this example the BIM maps transfer poorly functions from the source meshes to the target meshes while our method corrects these incorrect matches by providing a more accurate transfer. A visualization of the optimal vector field is provided on Figure 6.8 for the human example. The vector field on Figure 6.8b corresponds to the displacement needed to repair the BIM map, the action of this correction can be seen on the upper row of Figure 6.7.

r(t) → max d N (p,q)≤t d M (ϕ(p), ϕ(q)), (6.8)

Parameters Dependence

In practice we consider only a small family of vector fields based on the first eigenfunctions of 1-form Laplace-de Rham operator. Therefore in this setting our method will be more efficient in repairing low frequency distortion rather than recovering a high frequency deformation that cannot be represented by the flow of low frequency vector field. Of course the bigger is the vector field basis the better will be the repairs, and the slower The texture of the first column (6.6a, 6.6d) are transferred to the second using the direct conversion of a functional map (6.6b, 6.6e) and to the third using our method (6.6c, 6.6f).

will be the method. In the experiments we presented the dimension of the vector field family can be reduced to 40 without influencing too much the point-to-point map.

Another critical parameter is the number of eigenfunctions we choose to represent the functional map C and C ϕ 0 . If the deformation is nearly isometric a small number is sufficient as the functional map C is almost diagonal. These considerations also apply to the initial map ϕ 0 : a very distorted map is badly approximated by a small number of eigenfunctions and can severely influence our method. We found that lowering the size a the function basis under 150 degrades rapidly the quality of the results. In principle our method should work for non-isometric deformations provided we are given high-quality functional map as input. To obtain such a map, the choice of the functional basis would have to be modified in order to successfully encode the functional map in a reduced basis. This direction is left as an interesting future work.

Performance

For performance evaluation the computation times are given in the Table 6.1 in various cases. All the experiments have been performed on laptop with a 1.4 GHz processor and 4Go memory without parallelization. The timings are given for the two steps of the method: solving the problem (6.6) and tracing the flow lines. The time spent solving the optimization problem is almost independent of the number of vertices. The size of

(a) (b) (c) (d) (e) (f)
Figure 6.7 -Transfer of a function on the source meshes 6.7a and 6.7d to the target meshes using BIM 6.7b and 6.7e compared to our method 6.7c and 6.7f. this problem depends only on the number of computed eigenfunctions of the Laplace-Beltrami operator and on the dimension of the vector field family, which are constant in all experiments. Note that the computation of the flow does not scale linearly with the number of the vertices. This is explained by the fact we compute a composition with the BIM map which may map many vertices to a single point.

Mesh

Conclusion, Limitations and Future Work

In this chapter we presented a method for non-rigid shape matching that is designed to output continuous maps. Our approach combines the strengths of conformal-based approaches, which often guarantee continuity with the functional map framework, which can enable low-distortion maps on the space of functions. Key to our method is enforcing continuity via the flow of a vector field, which allows our method to remain efficient by

.8 -Visualization of the direction of the optimal vector field corresponding to the experiment 6.7c: complete shape 6.8a and close-up on the face 6.8b.

avoiding expensive pairwise vertex constraints. One of the limitations of our method is that we only approximate the flow of a single vector field, whereas in practice, for complex motions, a combination of flows may be necessary. Extending our method to such cases is possible, while taking care of the robustness and non-accumulation of numerical errors. We are also planning to wider arrays of initial maps and ways to incorporate the continuity directly in the optimization of the functional maps.

Part II

Deformation to Shape

In the previous two chapters we showed how one can compute a mapping between deformable shapes. In this part, we take a more careful look at the link between the deformation of shapes and the induced metric distortion. In particular we aim at creating a characterization of surface embeddings within the operator representation framework fit for synthesis and exploration of intrinsic and extrinsic deformations. Namely, in Chapter 7 we propose an algorithm to convert shape differences to embedded surfaces. In Chapter 8 we take a different point of view and introduce a characterization of deformation fields up to isometric deformation.

Chapter 7

Functional Characterization of Intrinsic and Extrinsic Geometry

We propose a novel way to capture and characterize distortion between pairs of shapes by extending the recently proposed framework of shape differences built on functional maps. We modify the original definition of shape differences slightly and prove that, after this change, the discrete metric is fully encoded in two shape difference operators and can be recovered by solving two linear systems of equations. Then, we introduce an extension of the shape difference operators using offset surfaces to capture extrinsic or embedding-dependent distortion, complementing the purely intrinsic nature of the original shape differences. Finally, we demonstrate that a set of four operators is complete, capturing intrinsic and extrinsic structure and fully encoding a shape up to rigid motion in both discrete and continuous settings. We highlight the usefulness of our constructions by showing the complementary nature of our extrinsic shape differences in capturing distortion ignored by previous approaches. We additionally provide examples where we recover local shape structure from the shape difference operators, suggesting shape editing and analysis tools based on manipulating shape differences.

Introduction

One classic approach to comparing surfaces separates metrics of similarity into intrinsic and extrinsic measurements. Intrinsic quantities are those that can be expressed exclusively in terms of distances along the surface, whereas extrinsic quantities are those that must be defined using an embedding into space and/or surface normals. A crowning result of classical differential geometry describes local geometry in terms of two quantities: the first and second fundamental forms, which capture the intrinsic Gaussian and extrinsic mean curvatures, respectively [START_REF] Bonnet | Mémoire sur la theorie des surfaces applicables sur une surface donnée[END_REF]. Considerable research in geometry processing has been dedicated to measuring intrinsic and extrinsic curvature in an attempt to replicate this attractive characterization of shape. From a practical standpoint, however, this task remains challenging for potentially noisy or irregular meshes considered in geometry processing. After all, surface curvature is a second-derivative quantity whose approximation on a piecewise-linear mesh requires considerable discretized adaptation and mollification to deal with noise. Measurement of curvature aside, algorithms for recovering geometry from discrete curvatures remain difficult to formulate, leading to potentially non-invertible discretizations. In this chapter, we formulate an alternative characterization of surface geometry better suited for analysis, comparison, and synthesis tasks in the discrete setting. Several desiderata inform our design; a suitable framework for representing shape should

• capture and distinguish intrinsic and extrinsic geometry,

• express shape properties in a multiscale fashion to distinguish noise and fine-scale detail from large-scale structure,

• come from a smooth theory of shape to provide insensitivity to tessellation,

• be naturally expressible on continuous surfaces and on triangle mesh discretizations, and

• admit an inverse operator for reconstructing the embedded shape.

In short, we wish to pass from pointwise embeddings to a "dual" space featuring a more democratic treatment of intrinsic and extrinsic shape properties. We approach this task by extending the theory of shape differences, introduced by Rustamov et al. [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] for purely intrinsic comparisons of shape structure. Rather than defining a shape in isolation, their construction characterizes shape by considering the distortion or difference of the target shape from a fixed source shape given a functional map between them [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF]. Shape differences are couched in the language of functional analysis, indirectly measuring changes in angles and distances through the effects of these changes on inner products of functions and their gradients. They are written as linear operators whose restriction to a multiscale basis like the Laplace-Beltrami eigenfunctions distinguishes features at different levels of detail and allows for straightforward discretization via piecewise-linear finite elements (FEM).

We modify and extend this framework to derive a shape representation that is complete, encoding both the intrinsic and the extrinsic distortion without loss of information in both the continuous and discrete cases. To this end, we begin by reexamining the discretization of shape differences on triangle meshes. We modify the original definition of discrete area-based shape difference and prove an analog of a continuous property mentioned in [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] that shape differences fully capture intrinsic structure. Inspired by this fully-discrete result, we proceed to ask whether shape differences also can capture extrinsic structure. Towards this goal, we define an additional pair of shape differences on a thickened surface that captures extrinsic geometry. We then show that our full set of differences is sufficient to reconstruct a shape up to rigid motion in the discrete setting, under mild assumptions. Our discussion concludes by closing the loop between shape differences and embeddings.

In particular, we provide algorithms for recovering a shape embedding given a combination of intrinsic and extrinsic shape differences. Our algorithms are guaranteed to succeed in the presence of complete information under weak genericity conditions; when shape differences have been truncated to a smaller basis, we suggest regularizers to recover reasonable estimates of the missing data via convex optimization.

To summarize, our main contributions are:

• Theoretical discussion establishing that properly modified shape difference operators from [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] fully encode the intrinsic metric of a triangle mesh. Unlike more direct representations of meshed edge lengths, these operators enjoy connection to smooth theory-providing some degree of tessellation invariance-as well as multiscale approximation in the Laplace-Beltrami basis.

• A novel set of shape differences aimed at capturing and characterizing extrinsic or embedding-dependent information, with an associated observation that generically this set of shape differences is complete and encodes shapes up to rigid motion.

• A set of approaches for recovering geometric structure and an embedding from the shape differences that reduce to a sequence of linear solves with theoretical guarantees of recovery in the presence of complete information, and which we apply to shape editing operations based on manipulating shape differences.

We demonstrate the usefulness of these contributions on a variety of tasks, ranging from the analysis of cloth simulations by using our novel shape difference operators, sensitive to changes in extrinsic structure, to the transfer of shape structure such as geodesic distances more accurately, compared to using functional maps, and finally to recovering shape embedding even in the presence of approximate functional correspondences.

Related Work

Representation and manipulation of extrinsic and intrinsic structure is a vast theme pervading the geometry processing literature. We refer to [START_REF] Botsch | Polygon Mesh Processing[END_REF] for discussion of the basic questions of representation and interaction with continuous differential geometry. Here, we highlight research linked to our particular approach.

Functional maps Our goal of using functional maps to characterize local and global geometry builds upon the machinery of shape differences [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF]; see §7.4 for a summary. Rustamov and colleagues [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] show that in the case of smooth surfaces, shape differences fully encode intrinsic geometry. They do not, however, pursue a corresponding analysis for the discrete case. Furthermore, their work focuses solely on intrinsic geometry and hence cannot characterize extrinsic bending, critical for describing differences between nearly-isometric shapes like articulated bodies and cloth.

Shape-from-Laplacian Recovering structure from intrinsic shape differences is closely linked to recovering structure from Laplacian operators. Both in the continuous [START_REF] Rosenberg | The Laplacian on a Riemannian Manifold[END_REF] and discrete [START_REF] Zeng | Discrete heat kernel determines discrete Riemannian metric[END_REF] cases, the Laplace-Beltrami operator fully encodes intrinsic surface geometry, namely the Riemannian metric for smooth manifolds and edge lengths for discrete meshes. For triangle meshes, de Goes and colleagues [START_REF] De Goes | Weighted triangulations for geometry processing[END_REF] provide convex machinery for recovering the intrinsic structure of the mesh; their encoding of intrinsic structure using only Laplacian matrices is more compact than our pair of area and conformal shape differences, at the cost of a nonlinear objective sensitive to incomplete information. The theoretical and practical contributions proposed in this chapter provide considerable insight beyond the fundamental mathematical contributions in these other works. Specifically, the convex optimizations in [START_REF] Zeng | Discrete heat kernel determines discrete Riemannian metric[END_REF][START_REF] De Goes | Weighted triangulations for geometry processing[END_REF] operate in the case of complete, noise-free information. They cannot be used for projection-style problems, e.g. finding the closest set of edge lengths to a noisy input Laplacian approximation or to finding an intrinsic structure consistent with a truncated spectral approximation of the full operator. Additionally, we show how to use related machinery to encode extrinsic bending rather than only edge lengths.

Encoding extrinsic geometry A natural question is whether intrinsic structure can be used to reconstruct a surface embedding up to a global rigid transformation. Numerous examples of isometric smooth surface pairs disprove this notion in the continuous case [START_REF] Kreyszig | Differential Geometry[END_REF].

While exact isometries of triangle meshes are rare with the exception of inward/outward "popping" of valence-three vertices, near-isometries can often arrise and have signficant differences in the embedding, making shape recovery from intrinsic data like edge lengths a numerically ill-conditioned problem; these near-isometries appear because small variations in the input edge lengths can lead to large changes in the resulting embedding. Nevertheless, Boscaini and colleagues [START_REF] Boscaini | Shape-fromoperator: Recovering shapes from intrinsic operators[END_REF] provide an algorithm for recovering a surface embedding in R 3 from shape differences or equivalent structures. They apply the SMA-COF algorithm [START_REF] Leeuw | Applications of convex analysis to multidimensional scaling[END_REF] for multidimensional scaling to generate an extrinsic embedding that replicates shape differences in a least-squares sense. As an alternative, [START_REF] Panozzo | Frame fields: Anisotropic and non-orthogonal cross fields[END_REF] propose an algorithm for embedding from local approximations of the metric tensor; we will use an extension of this algorithm in §7.7.3. Both of these methods, however, operate using only intrinsic information and are subject to the ambiguity and instability caused by isometry invariance.

Adding extrinsic information to a shape representation allows it to be embedded in R 3 up to rigid motion. In theory, the Gauss-Codazzi equations fully characterize surfaces from the first and second fundamental forms [START_REF] Bonnet | Mémoire sur la theorie des surfaces applicables sur une surface donnée[END_REF] (see [24, pg. 236]). In geometry processing, [START_REF] Eigensatz | Curvature-domain shape processing[END_REF] reconstructs surfaces from prescribed principal curvatures, while [START_REF] Fröhlich | Example-driven deformations based on discrete shells[END_REF] use nonlinear optimization methods to recover shape from dihedral edge lengths. These methods and many subsequent techniques employ nonlinear least-squares fits with few guarantees or characterization of their behavior. Wang, Liu, and Tong [START_REF] Wang | Linear surface reconstruction from discrete fundamental forms on triangle meshes[END_REF] propose a linear technique for embedding meshes from their edge lengths, dihedral angles, and axes of rotation across mesh faces.

In this chapter, we make use of offset surfaces to introduce extrinsic information to the shape difference representation. Offset surfaces have appeared in geometry processing for some related tasks, including cage generation [START_REF] Ben-Chen | Spatial deformation transfer[END_REF] and shape optimization for printing [START_REF] Musialski | Reduced-order shape optimization using offset surfaces[END_REF]. While techniques like [START_REF] Jung | Self-intersection removal in triangular mesh offsetting[END_REF] are needed to generate "clean" offset surfaces for geometry editing purposes, in our case self-intersection and related artifacts are acceptable since the offset surface is not used for display but rather for geometric computation. [START_REF] Cohen-Steiner | Restricted delaunay triangulations and normal cycle[END_REF][START_REF] Hoffmann | A discrete parametrized surface theory in R 3[END_REF] provide curvature theories for discrete surfaces using offset geometry.

Overview

Our two main goals are to modify and extend the definition of the shape difference operators of Rustamov [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] so as to capture extrinsic distortion and to facilitate shape inference, i.e. to recover the metric and potentially the embedding of a target shape, given a base shape and a collection of shape differences.

We achieve these goals in several stages. The main ingredient for constructing both smooth and discretized shape differences is the computation of inner products between functions. So, rather than working directly with shape differences, we largely focus on matrices of functional inner products, which can be constructed on a single shape rather than a pair. Hence, after reviewing the smooth construction of shape differences ( §4.2), we reexamine the discretization of intrinsic inner products on triangle meshes and show how a simple modification of the area-based inner product fully encodes intrinsic geometry in an easily-inverted fashion ( §7.4.1).

We then capture extrinsic shape structure by introducing two operators built from intrinsic inner products on offset surfaces of a base shape ( §7.5). We accompany our construction with theoretical characterization of the new information provided by extrinsic products ( §7.5.3) and conclude by making explicit how our constructions involving inner product matrices apply to the construction of differences between shapes ( §7.6). In this section, we also consider how truncating shape differences written in the Laplace-Beltrami basis affects the linear systems we pose.

With our new definitions and analysis in place, we propose optimization procedures for recovering intrinsic and extrinsic shape structure from the shape difference operators, potentially expressed in a reduced basis ( §7.7). While the basic machinery for recovering metric information from shape differences is purely linear, we propose the use of more general convex optimization tools that add resilience to noise and incomplete information by explicitly enforcing the triangle inequality and/or smoothness. We conclude by demonstrating the ability of our constructions to capture and characterize extrinsic distortion ignored by previous approaches ( §7.8.1). We furthermore apply our methods to recovering the metric and shape embedding and to facilitating novel shape editing operations via manipulating shape difference operators ( §7.8.4).

Structure of Discrete Inner Products

By examining the derivation of formulas for computing shape differences, we can reveal how they are related to local surface geometry. This analysis not only elucidates the information encoded in a given shape difference but also will inform our design of algorithms for recovering shape embeddings from shape differences.

Discrete Inner Products

Each quantity above is straightforward to discretize in the language of finite elements over a triangle mesh; see [START_REF] Brenner | The mathematical theory of finite element methods[END_REF][START_REF] Sayas | A gentle introduction to the finite element method[END_REF][START_REF] Strang | An analysis of the finite element method[END_REF] for general introductions to this approach. To this end, suppose M is represented using a connected, orientable, and manifold triangle mesh with vertices set X and triangles set F. We model functions as vectors f ∈ R |X | interpolated to triangle interiors in piecewise-linear fashion. We will begin our fine-grained examination of shape differences by posing functional inner products on these meshes in terms of discrete geometry. Our ultimate goal is to show that before truncation in a basis, the area-based and conformal inner product matrices completely encode the intrinsic structure of meshed geometry. This property is also stated in [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] in the continuous case; their discretization, however, does not admit such completeness due to the use of lumped area weights, as explained below. Consider a single triangle T ∈ F on M , and suppose f and g are affine functions on T ; in other words, f and g are evaluated in the interior of T via barycentric interpolation of the three scalar values f 1 , f 2 , f 3 and g 1 , g 2 , g 3 defined on the vertices of the triangle. Multiplying these functions and integrating reveals that the inner product of f and g on T is given by

f, g T L 2 = µ(T ) 12 f 1 f 2 f 3   2 1 1 1 2 1 1 1 2     g 1 g 2 g 3   ,
where µ(T ) is the area of T and f i , g i denote the values of f, g on vertex v i . As a sanity check, taking f i = g i = 1 ∀i recovers the area of T . This is the exact L 2 inner product of f and g defined over the meshed surface using piecewise-linear interpolation, without mass lumping commonly introduced in finite element discretizations; this distinction is critical for our construction.

Taking inner products over all of M requires summing over triangles

T . If f, g ∈ R |X | , then f, g M L 2
is given by f ⊤ Ag, where

A vw = 1 12 •    2 T ∼v µ(T ) when v = w T ∼e µ(T ) when e = (v, w) 0 otherwise, (7.1)
where T ∼ v denotes iteration over triangles adjacent to v and T ∼ (v, w) denotes iteration over triangles adjacent to edge (v, w). This |X | × |X | "Galerkin mass matrix" A is nondiagonal but positive definite, integrating products of piecewise linear functions exactly. See e.g. [START_REF] Strang | An analysis of the finite element method[END_REF]Chapter 10,[START_REF]Continuous matching via vector field flow[END_REF]] for an example of its appearance in finite elements. We can think of A as a linear operator A(µ) : R |F | → R |X |×|X | that constructs the area-based functional inner product matrix A given a vector µ ∈ R |F | of triangle areas. We can show that A(•) is invertible in the following sense:

Proposition 7.4.1 Suppose M has a boundary or at least one interior vertex with odd valence. Then, A(µ) uniquely determines µ, recoverable via a linear solve.

Proof Equation (7.1) gives A as a linear function of µ(•). Hence, we must show that this formula is invertible. First we show how to recover the area of a single triangle on M . By the second row of (7.1), given A we have the sum of triangle areas adjacent to any edge of M . If M has a boundary, we then know the areas of the boundary triangles. Otherwise, take v with odd valence, and enumerate its adjacent triangles as T 1 , . . . , T k for odd k. Since we know the sums of adjacent areas, we have a linear system to recover µ(T 1 ), . . . , µ(T k ):

       1 1 1 1 . . . . . . 1 1 1 1             µ(T 1 ) µ(T 2 ) . . . µ(T k )      =      A 12 A 23 . . . A k1      . w v T T ′ ℓ vw ℓ v ℓ w ℓ ′ v ℓ ′ w Figure 7.1 -Notation for the conformal product W .
Consider carrying out forward substitution on the matrix. In each iteration, only the bottom row changes, from (1, 0, . . . , 0, 1) to (0, -1, 0, . . . , 0, 1), then to (0, 0, 1, 0, . . . , 0, 1) and so on with alternating sign. When k is odd, in the last step the 1 is augmented to a 2, making the final row (0, . . . , 0, 2). In other words, the matrix reduces to an upper triangular matrix with nonzero diagonal, which is invertible. Hence, in either case we can recover µ(T ) for at least one T . The remaining areas can be computed by flood filling outward from T ; given the area on one side of an edge and the sum of the adjacent areas, the adjacent area is recovered by subtraction.

The proof of this proposition and others below is in the appendix. A proposition of this nature does not hold if masses are lumped down the diagonal of A. This observation is intuitive in that a triangle mesh has approximately two times the number of triangles as vertices. If f and g are piecewise-linear functions on M , then their gradients are piecewise-constant and expressible using one vector per triangle. Taking dot products of these gradients and integrating over M shows that f, g M H 1 0 = f ⊤ W g, where

W vw = 1 8 •        µ(T ) -1 (ℓ 2 vw -ℓ 2 v -ℓ 2 w ) +µ(T ′ ) -1 (ℓ 2 vw -ℓ ′2 v -ℓ ′2 w ) when v ∼ w -u∼v W uv when v = w 0 otherwise. (7.2)
Notation for the v = w case is shown in Figure 7.1; e ∼ v denotes an edge e adjacent to vertex v, and ℓ uv is the length of the corresponding edge. This matrix is the familiar cotangent Laplacian matrix cast in terms of edge lengths and triangle areas; this form also appears e.g. in [START_REF] Boscaini | Shape-fromoperator: Recovering shapes from intrinsic operators[END_REF]. Comparing (7.1) and (7.2), scaling the edge lengths of a mesh by some factor α will correspondingly scale A by α 2 while W will be left unchanged; unless otherwise noted, we scale meshes in our experiments to have unit surface area to remove dependence on global scaling.

A crucial observation that we make here is that if the triangle areas encoded in µ are fixed then the mapping W (ℓ 2 ; µ) : R |E| → R |X |×|X | taking squared edge lengths ℓ 2 ∈ R |E| to a conformal inner product matrix W is linear. Note also that W is fully determined by its values W vw for v ∼ w. Thus, if we represent the list of inner products W vw as a vector c in R |E| then for a fixed set of area weights µ, there exists a matrix B µ such that c = B µ ℓ 2 . The entries of B µ are, of course, given in Eq. 7.2. In the pipeline that we propose below, we will first recover the triangle areas and then use those to recover edge lengths from the corresponding inner products. The following proposition shows that "generically" the matrix B µ is invertible. I.e., the set of weights µ for which B µ is singular has measure 0.

Proposition 7.4.2 Assume that the mesh M is manifold without boundary. Then, for almost all choices of areas µ, the map W (ℓ 2 ; µ) uniquely determines ℓ, which is recoverable via a linear solve.

Proof By construction W (ℓ; µ) takes squared edge lengths ℓ and outputs the matrix W . Extracting elements of W corresponding to edges on M yields a linear operator B : R |E| → R |E| with matrix

B ij = 1 8    µ(T i ) -1 + µ(T ′ i ) -1 if i = j -µ(T ) -1
if i, j are edges of T 0 otherwise.

Here, indices i, j refer to edges on M ; for a given edge i, we label its adjacent triangles T i and T ′ i . Remark that B can be written as a weighted sum:

B = k 1 8 µ(T k ) -1 B k ,
where each B k is a matrix such that:

B k ij =    1
when i = j, and i belongs to triangle k. -1 when i, j are edges of triangle k. 0 otherwise.

It is easy to see that the intersection of the kernels of all B k is empty, since B k is nonsingular when restricted to the values on edges of triangle k. Moreover, by considering the determinant of B as a multivariate polynomial with real coefficients, we conclude that B is either singular for any choice of values of µ(T k ) -1 , or for a finite set of coefficients, which thus have measure zero.

To complete the proof we note that if B is singular for any choice of values of µ(T k ) -1 , then the matrix pencil B = k a k B k is singular (i.e., B is singular for any choice of coefficients a k ). Using Lemma 3.4 from [START_REF] Muhic | On the singular two-parameter eigenvalue problem[END_REF] and the fact that B k are symmetric, we see that in that case for every choice of a k , there must exist a vector x such that x T B k x = 0 for every B k , and Bx = 0. Now, given the values of x on some triangle, this means that its values on the adjacent triangle are either uniquely determined by the corresponding two equations (one linear, one quadratic), or these equations cannot be satisfied. By inspecting the resulting equalities, it is easy to see that at least two of the values on every triangle must be equal, and by considering any closed loop of triangles, these equations cannot be consistent for every choice of weights a k . Thus, B cannot be a singular matrix pencil, and therefore B is invertible for almost any choice of values µ(T k ) -1 .

This proposition implies that the linear map W (ℓ 2 ; µ) is invertible for a small (possibly zero) perturbation of any set of area weights µ. Nevertheless, there exist cases in which Figure 7.2 -A mesh for which W (ℓ 2 ; µ) is not invertible when µ = 1.

the squared edge lengths are not recoverable via inversion of the linear map W (ℓ 2 ; µ) for a fixed set of area weights. One example of such a shape is shown in Figure 7.2, consisting of two tetrahedra glued at their bases. In this case, all the triangles have equal area weights, and it can be seen that the resulting linear system is singular. We also remark that the condition of no boundary is necessary in the Prop. 7.4.2 above, as it is possible to construct meshes for which the map W (ℓ 2 ; µ) is singular for all choices of µ (e.g., a pair of triangles glued along a shared edge). For all the meshes that we tried in practice ( §7.8), we have observed that the resulting system is both invertible and typically well-conditioned. We leave the formulation of the necessary and sufficient conditions on the mesh and the weights µ for the invertibility of W (ℓ 2 ; µ) as a question for future work.

Encoding Extrinsic Structure

Intrinsic inner products capture the metric tensor (first fundamental form) of a surface, so to complete our representation we show how a related structure can be used to encode its second fundamental form. In keeping with previous discussion, we will use additional inner product matrices to derive a multiscale representation of this missing information. While there exist many possible ways to measure extrinsic distortion, this "functional" language facilitates a connection between continuous and discrete characterizations and unifies our treatment of intrinsic and extrinsic distortion. These added structures complement the area-based and conformal products by making our representation of a shape unique up to rigid motion. In addition to providing a lossless representation of surface geometry in the presence of complete information, we demonstrate how the new products can capture and encode geometric relationships that are not captured by purely intrinsic analysis.

Extrinsic Alternatives

In discrete language, the inner product matrices A(µ) and W (ℓ 2 ; µ) determine the edge lengths of a triangle mesh but not its dihedral angles, illustrated in Figure 7.3(a). Additionally providing dihedral angles is sufficient to recover a mesh up to rigid motion. There are many expressions of extrinsic shape that potentially encode these angles; before presenting our final solution, we mention a few straightforward alternatives to explain why they are less desirable.

At the most basic level, any technique encoding one value per edge of a triangle mesh could be used to represent dihedral angles. For instance, since the angles are in a vector θ ∈ R |E| , we could use an analog of Proposition 7.4.2 to store them in the matrix W (θ; µ). This matrix roughly corresponds to taking products of functional gradients under the second fundamental form h as h(∇f, ∇g) dA. The resulting matrix is not positive semidefinite, however, which prevents the definition of a smooth analog via the Riesz Representation Theorem (which applies only to positive definite inner products) and causes numerical issues due to departure from the cone of semidefinite matrices. Dihedral angles also are known only up to a period of 2π, providing potential for ambiguity in the expression of the vector θ. Lastly, if the mesh is flat, a definition based on dihedral angles will lead to a degenerate inner product.

In an attempt to bring back the positive definiteness enjoyed by the intrinsic formulation, we might attempt to encode the edge lengths of a dual mesh, shown in Figure 7.3(b). These lengths indirectly encode dihedral angles up to sign but are unable to distinguish between inward and outward folding directions, as shown in the figure. Obvious techniques for disambiguating the inward and outward folds generally accompany edge lengths with signs, reintroducing the problems discussed in the previous paragraph. An alternative construction might define extrinsic shape differences via the Gauss map, or map from a surface into the unit sphere based on normal direction; see [START_REF] Meyer | Discrete differentialgeometry operators for triangulated 2-manifolds[END_REF] for an example in geometry processing. While the Gauss map is used in classical differential geometry to derive extrinsic properties of surfaces, we find it to be unstable within the shape difference framework. In particular, the image of the Gauss map is composed of many overlapping spherical triangles that change rapidly from vertex to vertex. Projection of this information into low-frequency Laplace-Beltrami bases tends to remove the majority of the meaningful geometric signal. In a sense, however, we can view the offset surface construction proposed below as a means of smoothing out this construction. Before proceeding, we should remark that strictly speaking it may not be necessary to provide extrinsic information at all. According to a classical result by Gluck [START_REF] Gluck | Almost all simply connected closed surfaces are rigid[END_REF], almost all triangulated simply connected closed surfaces are rigid. Although this result might imply that triangle edge lengths are, in general, sufficient to reconstruct the mesh up to rigid motion, this is only true if the metric is known exactly; moreover, it is highly ℓ decreases ℓ increases Figure 7.4 -Edge lengths change according to curvature of the offset surface. nontrivial to recover an embedding even if one is known to exist. When the edge-lengths are perturbed or are approximated, the corresponding embedding might either not exist or be very far from the desired shape. As we show below, the presence of explicit information about the extrinsic distortion can greatly help in both direct and inverse problems and is largely complementary to the intrinsic distortion addressed by prior methods in practice.

Offset Surfaces

Our construction of an extrinsic shape representation is an extension of the dual mesh idea from §7.5.1 that does not suffer from sign ambiguity. Instead, we are able to rely upon the positive definiteness of inner product matrices directly to encode both intrinsic and extrinsic information. In short, rather than encoding a metric and its derivative, we encode a metric and a slightly deformed metric, both of which admit natural positive definite representations. The intuition for our construction is illustrated in Figure 7.4. Suppose we wish to recover the embedding of the blue torus. As discussed in the previous section, it may be difficult to reconstruct the torus purely from its list of edge lengths. Instead, suppose we generate an offset surface by displacing each vertex and face along its outward normal a fixed distance t. The operation is extrinsic, since the mesh moves through the surrounding space, modulating edge lengths ℓ based on the curvature of the surface. The edge lengths in the interior of the torus shrink while the edge lengths on the exterior expand, effectively distinguishing the bend direction. In the continuous case, we can formalize the effect of offsetting a surface as follows: Proposition 7.5.1 Suppose M is a compact orientable Riemannian 2-manifold without boundary. Consider a family of immersions

F t : U ⊂ R 2 → M t ⊂ R 3 satisfying ∂F t ∂t (p) = n t (p), ∀(p, t) ∈ U × R + ,
where n t denotes the outward unit normal of M t := F t (U ). At all time t the normal n t remains equals to the normal n 0 of M 0 := M .

Moreover, if g ij (t) := ∂Ft ∂x i , ∂Ft ∂x j is the metric of the embedded surface and h ij (t) := ∂Ft ∂x i , ∂n ∂x j is its second fundamental form, then

g ij (t) = g ij + 2th ij + t 2 h il g lm h mj and µ t = (1 + tH + t 2 K)µ,
where H := (g ij h ij ) t=0 is the mean curvature, 2K := H 2 -(g ij h jk g kl h il ) t=0 is the Gaussian curvature at t = 0 and µ t := det g(t).

Proof First let's show that the normal vector is constant during the flow. The vector n is unit norm therefore ∂ t n t , n t = 2 ∂ t n t , n t = 0 meaning that ∂ t n t belongs to the tangent of M t . Since any tangent vector V can be written

V = V, ∂ i F g ij ∂ j F , we have ∂n t ∂t = ∂n t ∂t , ∂F t ∂x i g ij (t) ∂F t ∂x j = -n t , ∂n t ∂x i g ij (t) ∂F t ∂x j = 0.
As consequence the immersion has the following closed form expression at all time: F t = tn 0 + F 0 . Therefore the metric is quadratic with respect to time:

g ij (t) = ∂F t ∂x i , ∂F t ∂x j = t ∂n 0 ∂x i + ∂F 0 ∂x i , t ∂n 0 ∂x j + ∂F 0 ∂x j = g ij + 2th ij + t 2 ∂n 0 ∂x i , ∂n 0 ∂x j .
The proposition statement is obtained using the Gauss-Weingarten relations linking ∂ i n with the second fundamental form. Because the normal is unit norm ∂ i n is a tangent vector whose expression is given by:

∂n 0 ∂x i = ∂n 0 ∂x i , ∂F 0 ∂x j g jk ∂F 0 ∂x k = h ij g jk ∂F 0 ∂x k .
The evolution of the local area is obtained by a direct computation of the determinant of metric tensor.

Results of this nature are fairly well-known for offset surfaces and sometimes referred to as Steiner formula; see e.g. [START_REF] Patrikalakis | Shape Interrogation for Computer Aided Design and Manufacturing[END_REF] for related discussion. Informally, the proposition shows that the second fundamental form of M is encoded through the change in metric while the surface is being offset along its normal directions. In case of surfaces, a closed form expression for h can be derived. preciselly, there is only one square root of 2-matrices (out of four) satisfying ∂ t g ij = 2h ij at time zero:

h ij (t) = 1 τ t t g ij (t) + 1 t δ t τ t -1 g ij , -ǫ < t < ǫ
where δ t = µ(t)/µ, and

τ t = g ij g ij (t) + 2δ t .
When M is an oriented triangle mesh, there are many potential constructions of discrete offset surfaces, and several likely would suffice for the proofs in this chapter. For Bottom layer Top layer Figure 7.5 -Topology of offset mesh.

mathematical simplicity, we choose the construction in Figure 7.5. On the left we show triangles of the original mesh M in blue. On the right, we define the topology of the offset mesh in red, which contains a vertex for every vertex of M and every triangle of M . For a fixed constant t > 0, we place the vertices distance t above M along its face/vertex normals; any reasonable definition of a unit-length vertex normal suffices. Offset vertices associated with triangles are placed directly above the barycenter of the triangle.

Recovery of Embedding

In the end, we encode the geometry accompanying a fixed triangle mesh topology using four structures: the intrinsic area-based and conformal inner product operators and the same operators for the offset surfaces with fixed normal offset distance t > 0. We denote the offset surface of M as M t . In this section, we show-at least before truncation-that these four difference matrices are sufficient for fully reconstructing a shape. The challenge of reconstructing a triangle mesh from its edge lengths arguably comes from the fact that there are many ways to glue together two adjacent triangles by fixing different dihedral angles. Rigidity may imply that only one such embedding exists, but it is not obvious from local relationships. In contrast, an oriented tetrahedral mesh is easy to reconstruct from its list of edge lengths simply by gluing individual tetrahedra face by face. Hence, our intuition for why a mesh plus its offset are enough to reconstruct the mesh comes from a volumetric perspective. This intuition is confirmed by Proposition 7.5.1 and the Gauss-Codazzi equations, since offset geometry provides extrinsic curvature information, which in turn determines an embedding. We chose the topology of the offset mesh (Figure 7.5) specifically to allow for a canonical "thickening" of the offset slice into a tetrahedral mesh, shown in Figure 7.6. Using mesh-based shape differences, we can recover the edge lengths of the bottom and top layers of the thickening. By construction of the offset mesh, we are able to recover the lengths of the interior edges of the thickening, effectively proving the following proposition: Given the topology of M , the area-based and conformal product matrices A(µ) and W (ℓ 2 ; µ) of M , and the area-based and conformal product matrices A t (µ t ) and W (ℓ 2 t ; µ t ) of M t , the geometry of M can (almost always) be reconstructed up to rigid motion.

Proof The previous propositions show that µ, ν, µ t , and ν t are (almost always) sufficient to recover the edge lengths of the base and offset surfaces. The remaining edges of the canonical thickening are between the inner and outer layers and are recoverable essentially by convention. Specifically:

• The edges along surface normals are length t by definition.

• The bottom edge lengths of the "upward-facing" tetrahedra (Figure 7.6) are known because they are on the base surfaces. The remaining edges of these tetrahedra can be computed because the upward-facing tetrahedron is generated via normal offset from the barycenter of the base triangle by a distance t.

• "Mesh edge" tetrahedra are adjacent to "upward-facing" tetrahedra and outer faces of the thickening and hence have edge lengths fixed by their neighbors' construction.

• Similarly, "downward-facing" tetrahedra have one normal edge of length t, and the remaining edges are on the outer surface or adjacent to an "upward-facing" tetrahedron.

The embedding of a single oriented tetrahedron is fixed up to rigid motion given its edge lengths, so the proposition follows by gluing the tetrahedra of the canonical thickening according to the topology of the construction.

Discussion

We pause to summarize the theoretical development in the previous sections. We began by reconsidering the construction of inner products and shape differences from first-order finite elements. When area elements are not lumped, we showed that inner product matrices fully determine the edge lengths of a mesh and that they can be recovered by solving two linear systems of equations: one for recovering the triangle areas, and the other for recovering the edge lengths. Moreover, generically, both systems are non-singular.

In both the continuous and discrete cases, these intrinsic measurements are not enough to distinguish isometric shapes. Even worse, the space of near-isometric shapes can be very large. Hence, we propose generating an offset surface M t from a mesh or surface M . In the continuous case, the geometry of M t determines the extrinsic structure of M by encoding its second fundamental form. In the discrete case, combining edge lengths of M with edge lengths of M t fully determines M up to rigid motion. The main development is that we can completely determine a shape using functional inner products via the constructions above.

Our theoretical contributions deal with the noise-free, non-truncated case. Roughly, they show that if intrinsic/extrinsic shape differences were computed from an embedded mesh M with fixed topology, then the embedding of M almost always can be covered from those differences up to rigid motion. We evaluate sensitivity to noise and the possibility of recovering geometry from truncated shape differences empirically in §7.8.

From Inner Products to Shape Differences

With the goal of working with quantities that exist when meshes are not in vertex-forvertex correspondence, we shift from working with matrices of inner products to shape differences. This shift is needed to propose algorithms in §7.7 for estimating the dense structure of a target mesh given a source mesh and an approximate relationship between the source and the target, represented as a functional map computed e.g. using assorted correspondence techniques.

Discrete Shape Differences

We begin by considering two meshes M and N in vertex-for-vertex correspondence, with areas µ M , µ N ∈ R |F | and squared edge lengths ℓ 2 M , ℓ 2 N ∈ R |E| . Based on the continuous definitions in §4.2, the "full" area-based and conformal shape difference between meshes M and N are [106, §5, "option 1"]

D A = A(µ M ) -1 A(µ N ) D C = W (ℓ 2 M ; µ M ) -1 W (ℓ 2 N ; µ N ). (7.3)
A straightforward corollary of the discussion in §7.4.1 is that these two differences completely determine the edge lengths and triangle areas of N given the geometry of M . Notice the first relationship is still linear in µ N and the second in ℓ 2 N , preserving the proposed system of equations for reconstruction. Similarly, the extrinsic differences are simply the shape differences between the offset surfaces:

D E A = A(µ Mt ) -1 A(µ Nt ) D E C = W (ℓ 2 Mt ; µ Mt ) -1 W (ℓ 2
Nt ; µ Nt ).

(7.4)

The discussion in §7.5.3 implies that the tuple

(D A , D C , D E A , D E C
) is sufficient to reconstruct N up to rigid motion given M .

Source-Truncated Correspondence

More commonly, suppose Φ ∈ R |X |×k contains the orthonormal Laplace-Beltrami basis of M , truncated to k functions. Assuming M and N are still in vertex-for-vertex correspondence, we can write "reduced" shape differences as

D Φ A = Φ ⊤ A(µ N )Φ D Φ C = diag(-{λ M i }) + Φ ⊤ W (ℓ 2 N ; µ N )Φ, (7.5)
where the eigenvalues of the Laplacian on M are λ M i . These differences no longer determine angles and edge lengths exactly but still encode a multiscale notion of geometry that is valuable for understanding the relationships between M and N ; extrinsic differences can be defined analogously from the offset surface. We can still define linear systems for computing µ N and ℓ 2 N from D A , D C , µ M , and ℓ 2 N using these relationships, although they are unlikely to be full-rank for small k; we provide regularizers in the next section. These truncated differences essentially correspond to removing rows and/or columns from the full shape differences after writing them in the Laplace-Beltrami eigenbasis. Such a computation can be useful for multiscale analysis of surface deformations, in which vertex-for-vertex correspondence is known but high-frequency changes may not be useful to analyze. What remains, however, is to consider the case when M and N are not in vertex-for-vertex correspondence and both have incomplete bases.

Source-and Target-Truncated Correspondence

Suppose we are given truncated bases Φ M ∈ R |X M |×k M and Φ N ∈ R |X N |×k N for the eigenspaces of M and N , respectively, and a functional map matrix C ∈ R |X N |×|X M | taking functions written in the Φ M basis on M to functions in the Φ N basis on N . Following [106, §5], we define shape differences in this case as

D Φ M ,Φ N A = C ⊤ C D Φ M ,Φ N C = diag(-{λ M i }) + C ⊤ diag(-{λ N i })C. (7.6)
Whereas the truncated shape differences in (7.5) contain a limited window of values from the full shape difference matrix, in this final case the non-truncated entries of the shape difference matrices also undergo some change. This is because even if a function on M is in the column space of Φ M , it will not be transported fully to N by the functional map C due to removal of high frequencies. These shape differences are discretizations of analogous linear operators in the smooth setting. For this reason, even though the differences in (7.6) no longer satisfy exact equality relationships like those in (7.5) for recovering areas µ N and squared edge lengths ℓ 2 N from shape differences and the geometry of M , we will pose approximate relationships

D Φ M ,Φ N A ≈ (Φ M ) ⊤ A(µ M ←N )Φ M D Φ M ,Φ N C ≈ diag(-{λ M i }) + Φ M ⊤ W (ℓ 2 M ←N ; µ M ←N )Φ M . (7.7)
The unknown variables µ M ←N and ℓ 2 M ←N can be thought of as "pullbacks" of the metric of N to that of M , in the sense that they attempt to assign areas and edge lengths to the topology of M to mimic inner products on N . The first condition is linear in

µ M ←N ∈ R |F M | and the second in ℓ 2 M ←N ∈ R |E M | .
Since the shape differences in (7.6) are the most realistic test cases, we will assume in our experiments that truncated shape differences are computed in this fashion unless noted otherwise. That is, we will assume that we are given a source-and target-truncated shape difference. The experiments in §7.8.2 verify that this approximation is reasonable as long as k M and k N are sufficiently large.

Recovery of Intrinsic and Extrinsic Structure

Having established theoretical aspects of intrinsic and extrinsic shape differences, we now provide algorithms for recovering a shape N given a base shape M and shape differences to N and its offset. First, we recover triangle areas from the base and offset surfaces from corresponding area-based shape differences. With these areas fixed, we then recover edge lengths, which were shown in §7.5.3 to completely determine the surface. Both steps can be carried out using linear solves when shape differences are not truncated. When dealing with truncated or inexact functional maps, we augment the optimization with constraints ruling out unreasonable structures. We also show how to apply existing techniques for recovery of an embedding from edge lengths of the surface and its offset.

Triangle Area Computation

We first show how to recover areas of triangles given an area-based shape difference. Our approach is an extension of the basic linear technique outlined in the proof of Proposition 7.4.1, extended to deal with truncation and noise. Following §7.6, suppose D A is the area-based shape difference between M and N in the Laplace-Beltrami basis Φ M . Recall that our goal is to pull the geometry of N back to the mesh of M . Hence, the the area-based difference from M to the reconstructed target shape N * should satisfy

D * A = Φ ⊤ M A(µ N )Φ M .
If the reduced basis Φ M on M has k functions, this linear system for µ N has k 2 equations and |F| unknowns. So, we need at least k ∼ |F| to have a well-posed system. The quality of the solution found by solving this system without regularization depends on two factors: the quality of D A and the conditioning of the resulting linear problem. We find that both limitations are improved considerably by introducing a nonnegativity constraint, leading to the following optimization problem for µ N :

min µ N Φ ⊤ M A(µ N )Φ M -D A 2 Fro s.t. µ N (T ) ≥ 0 ∀ triangles T ∈ F. (7.8)
We solve this and other convex programs using the Mosek toolbox [START_REF] Mosek Aps | The MOSEK optimization toolbox for MATLAB manual[END_REF]. When the system is highly underdetermined, we additionally add a regularizing viscosity term

ε µ N -µ M 2 2
for small ε > 0, under the assumption that triangle areas should change minimally unless there is evidence to do otherwise; we set ε = 10 -4 in all the experiments in this chapter.

Edge Length Computation

Now that we can compute triangle areas, we can recover edge lengths. As in the last section, we start from Proposition 7.4.2 to propose a basic linear system for squared edge lengths and then provide regularization techniques for dealing with inexact or truncated differences.

The conformal shape difference encodes the transformation of the cotangent Laplacian through the deformation. Again borrowing from §7.6, the geometry of N can be pulled back to M via the following linear condition on squared edge lengths ℓ 2 N given fixed areas

µ N : diag(-λ M i )D C = Φ ⊤ M W (ℓ 2 N ; µ N )Φ M . (7.9)
Solving this linear system of equations for ℓ 2 N depends critically on the approximated areas µ N ; numerical or discretization error from the method in §7.7.1 invalidates this step, regardless of the quality of D C . To provide resilience to this issue and to noise in D C , we add constraints to this system ruling out unrealistic edge lengths ℓ 2 N . To define a triangulation, the squared edge lengths ℓ 2 N must be nonnegative; furthermore, ℓ 2 N (T ) must respect the triangle inequality in each mesh triangle T . We enforce the latter constraint via the following proposition: Proposition 7.7.1 The symmetric matrix E defined by

E = 1 2   2x 1 x 3 -x 1 -x 2 x 2 -x 1 -x 3 x 3 -x 1 -x 2 2x 2 x 1 -x 2 -x 3 x 2 -x 1 -x 3 x 1 -x 2 -x 3 2x 3  
is positive semidefinite if and only if x 1 , x 2 , x 3 are nonnegative and their square roots satisfy the triangle inequality.

Proof We denote (e 1 , e 2 , e 3 ) the canonical basis and the indices {i, j, k} ∈ {1, 2, 3}. If √ x k ≤ √ x i + √ x j , then there exist three points (v 1 , v 2 , v 3 ) which define an embedding of a triangle. Let E be the matrix with columns

v 3 -v 2 , v 1 -v 3 and v 2 -v 1 then E = E ⊤ E.
The matrix E is therefore positive semidefinite. Since E is symmetric positive semidefinite, the Cauchy-Schwartz inequality holds. Expanding the expression (e i + e j ) ⊤ E(e i + e j ) yields

(e i + e j ) ⊤ E(e i + e j ) = x i + x j + 2e ⊤ i Ee j ≤ x i + x j + 2 √ x i x j ≤ ( √ x i + √ x j ) 2 .
At the same time, a direct computation shows (e i + e j ) ⊤ E(e i + e j ) = x k which implies that

√ x k ≤ √ x i + √ x j .
We also can link squared edge lengths to the computed triangle areas µ N (T ). This link is provided by the submatrices E k defined as

E k = 1 2 2x i x k -x i -x j x k -x i -x j 2x j ,
where {i, j, k} = {1, 2, 3}. With this definition in place, we leverage the following proposition:

Proposition 7.7.2 E is positive semidefinite if and only if x k ≥ 0 for all k ∈ {1, 2, 3} and det(E 3 ) ≥ 0. Moreover, if E 0, then det(E k ) = 4µ N (T ) 2 .
Proof Starting with the second statement, direct computation of the determinant shows

4 det(E k ) = x 1 x 2 x 3   -1 1 1 1 -1 1 1 1 -1     x 1 x 2 x 3   .
When the x i 's are squared triangle edge lengths, this is a formulation of Heron's area formula.

The first statement is proved using a well-known theorem on positive block matrices (property of Schur complements) [START_REF] Boyd | Convex Optimization[END_REF]: E 0 if and only if E 3 0 and

2x 3 - x 2 -x 1 -x 3 x 1 -x 2 -x 3 ⊤ E -1 3 x 2 -x 1 -x 3 x 1 -x 2 -x 3 ≥ 0.
Notice that

E k 1 1 = - 1 2 x j -x i -x k x i -x j -x k ,
and hence the first condition is met whenever x 3 ≥ 0. Moreover, E 3 is a 2 × 2 matrix and therefore is positive semidefinite if and only if x 1 ≥ 0, x 2 ≥ 0 and det(E 3 ) ≥ 0.

Enforcing constraints derived from these relationships in the computation of edge lengths from a shape difference leads to the following optimization problem:

min ℓ 2 N Φ ⊤ M W (ℓ 2 N ; µ N )Φ M -diag(-λ M i )D C 2 Fro s.t. ℓ 2 N ≥ 0 det(E 3 (T )) = 4µ N (T ) 2 ∀ triangles T ∈ F.
This problem, however, is large and non-convex due to the determinant constraint. A convex relaxation is possible by noticing that the cone of symmetric positive semidefinite matrices with determinant ≥ 1 is convex; this observation derives from the convexity of the function A → -log(det A)) [START_REF] Boyd | Convex Optimization[END_REF]. So, the former problem can be relaxed to a convex problem:

min ℓ 2 N Φ ⊤ M W (ℓ 2 N ; µ N )Φ M -diag(-λ M i )D C 2 Fro s.t. ℓ 2 N ≥ 0 det(E 3 (T )) ≥ 4µ N (T ) 2 ∀ triangles T ∈ F. (7.10)
The determinant constraint is handled using the rotated quadratic cone optimization in the Mosek toolbox [START_REF] Mosek Aps | The MOSEK optimization toolbox for MATLAB manual[END_REF]. While (7.10) contains a relaxation of the full set of constraints, we find empirically that this relaxation generally is tight; we leave it to future work to prove conditions for "exact recovery" akin to those in [START_REF] Dym | Exact recovery with symmetries for Procrustes matching[END_REF] for mesh alignment problems. As in §7.7.1, we can additionally regularize by adding

ε ′ ℓ 2 N -ℓ 2 M 2
2 to the objective; our experiments use ε ′ = 10 -4 ..

Global Extrinsic Reconstruction

At this point, we have presented algorithms for recovering edge lengths for the entire canonical thickening defined in §7.5.3. As suggested in the proof of Proposition 7.5.2, if these edge lengths are computed without error, the thickening can be reconstructed greedily; then, the embedding of N from M is the inner envelope of this thickening. In reality, the squared edge lengths in ℓ 2 N likely exhibit numerical error. For this reason, we employ the algorithm in [START_REF] Panozzo | Frame fields: Anisotropic and non-orthogonal cross fields[END_REF] for reconstructing a triangle mesh given its edge lengths. We adapt their approach to take into account the tetrahedra defined by the offset surface, by using the same ARAP-style deformation energy, defined on each triangle facet of each tetrahedron, and using the same alternating optimization strategy. We note, in particular, that this approach does not require embedded surfaces to be manifold, and can easily incorporate edges shared by more than two triangles, which only changes the computation of the gradient of the energy. Hence this allows us to reconstruct the entire set of triangles in the canonical thickening rather than the inner or outer surfaces only. We provide the thickening of M as a starting point for their alternating optimization algorithm. Whereas their method is subject to isometric ambiguity when embedding manifold meshes, reconstructing the entire thickened structure reduces ambiguity and more reliably provides an extrinsically correct embedding.

Experiments

In this section we illustrate the utility of the constructions presented above in a variety of practical application scenarios. We start by showing how the extrinsic shape differences can be useful for shape exploration and analysis, by complementing the information provided by the intrinsic differences of Rustamov et al. [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF]. We then show how our metric and shape recovery methods can be used to both infer shape structure and ultimately recover the embedding from approximate, truncated shape differences.

Shape Space

An example application of shape differences that does not rely on exact reconstruction of local geometry involves the extraction of variability within a collection of related shapes. Suppose we choose an arbitrary base shape and compute its shape difference matrices with the remaining shapes in a collection. Then, a simple low-dimensional description of shape variability is to do PCA on the collection of matrices, resulting in the embedding of each shape as a point in PCA space. (no offset) (with offset) (with offset)

Figure 7.7 -PCA on collections of shape differences reveals the axes of variability within a collection; each shape on the left is colored the same as its corresponding points in the plots. The area-based and conformal differences are unable to distinguish the inward and outward bumps in the top example, leading to clusters of four points.

Figure 7.8 -Human models from Figure 7.7 sorted by the first PCA dimension for area-based shape differences (top) and area-based differences including an offset surface (bottom). The differences without offsets distinguish body type, while the differences with offsets distinguish pose.

We use PCA experiments to illustrate the power of our proposed extrinsic differences. For instance, Figure 7.7 illustrates embeddings of two-parameter shape collections into the plane using the procedure above (k M = 50, k N = 100). The top row illustrates the need for extrinsic differences most clearly. Here, we generate cubes with smooth bumps, smoothly varying from an inward bump to an outward bump. Intrinsic shape differences are identical for inward and outward bumps, leading to PCA embeddings that cluster sets of four shapes together. Adding extrinsic information disambiguates the embedding problem, separating the clustered points. Similarly, the extrinsic area-based shape difference best separates the parametric human models evenly among the two axes (57.2% variability along the principal axis, 38.3% along the secondary axis); interestingly, conformal shape differences among offset surfaces do not exhibit much variability for this particular class of surfaces. Figure 7.8 highlights how intrinsic and extrinsic shape differences can measure different properties of shape. We sort the collection of human models by the one-dimensional embeddings (x-axis) of intrinsic (top) and extrinsic (bottom) area-based shape differences (k M = 50, k N = 100). The intrinsic shape differences distinguish the body type of the model and are invariant to the pose of the arms and legs; this ordering reflects the property that articulated deformations of humans are nearly isometric. Complementing this embedding, the extrinsic differences distinguish pose and are less sensitive to body type. This property is also visible in Figure 7.7 since the area-based embeddings without and with the offset are transposed from one another. Figure 7.9 shows a similar experiment applied to shapes from individual frames of animation sequences. Both intrinsic and extrinsic shape differences are able to recover the cyclic structure of a galloping horse animation; this indicates that the galloping motion contains both intrinsic and extrinsic deformation modes. Contrastingly, the intrinsic differences severely underperform in recovering an animated sequence of deforming cloth. The physics of cloth naturally avoids intrinsic stretching and shearing, maintaining the initial developable structure. Thus, intrinisic shape differences provide little-to-no information, while the extrinsic differences capture the evolution of the animation. From a wider perspective, the experiments in this section reveal the value of explicitly representing both intrinsic and extrinsic deformation in navigating datasets of 3D surfaces. A sizable fraction of geometry processing algorithms, including the original work on shape differences, focuses on shape exploration based exclusively on intrinsic structure. Yet, motions like the deformation of a piece of cloth cannot be captured by this representation. While cloth deformation may be an extreme example, based on these results we advocate inclusion of both intrinsic and extrinsic structures in shape analysis rather than discarding the extrinsic information.

Effects of Truncation

The propositions in this chapter show that discrete shape differences completely encode geometric structure when they are written in a full basis. For many applications, however, we approximate shape differences in a truncated low-frequency spectrum. While the effects of this truncation are difficult to characterize mathematically, in this section we evaluate the effects of this approximation numerically.

There are two potential sources of truncation error in the twice-truncated differences discussed in §7.6.3: The choice of k M and the choice of k N . As mentioned in §7.6.3 these two parameters have slightly different effects; decreasing k M corresponds to removing rows or columns of the shape difference matrices, while decreasing k N can affect the values of the entries. Figure 7.10 illustrates the results of an experiment varying k M and k N for intrinsic shape differences and using the pipeline described in §7.7 to recover areas and edge lengths; recall that this technique extracts areas and edge lengths on N using calculations on M . We choose a pair of meshes with a ground-truth map to avoid additional error due to map approximation

(|X M | = |X N | = 1000).
Each color plot shows the relative error of assorted quantities extrapolated from the truncated shape differences: face areas (µ), edge lengths (ℓ), truncated eigenvalues ({λ} k N i=1 ), full eigenvalues ({λ}

|X | i=1
), and entries of the Laplacian (∆). We assume k N ≥ k M , providing the upper-triangular structure of the plots; the vertical axis represents k M (range: k M ∈ [60, 500]) and the horizontal axis represents k N (range: k N ∈ [START_REF] Khavari | Regge Calculus as a Numerical Approach to General Relativity[END_REF]500]). We choose ε so that the viscosity regularizer contributes < 10% of the optimal objective. These plots show that even truncated shape differences can be used to extract per-face and per-edge information about the mesh using our pipeline. Even with 15% of the Laplacian eigenvectors, we can relatively reliably extract the face areas and edge lengths of the target mesh. Even on challenging tasks like recovering the full spectrum of the target mesh-beyond the eigenvalues used to compute the shape difference-our algorithm has some success. The choice of k N is particularly important. Intuitively, this phenomenon might be explained by the fact that modulating k N changes the values in the shape difference matrices rather than just their size. The top row of the matrix also exemplifies a pattern we observed across our experiments; below a certain value for k M , there is not enough information to get a meaningful indication of local geometry from the shape difference matrices.

Intrinsic Recovery

The experiments in §7.8.2 illustrate a remarkable observation, that we are able to recover local information about the target of a shape difference from a truncated shape difference. That is, the nonnegativity and semidefinite constraints proposed in §7.7 paired with regularization are sufficient to avoid the null space of the truncated linear systems for recovering areas and edge lengths.

As the cotangent Laplacian of a triangle mesh (with or without area weights) is an intrinsic structure,we can use our computed vectors µ and ℓ to pull back the Laplacian operator from N to M . This technique is illustrated in Figure 7.11. In this experiment, we compute shape differences from a 100 × 200 functional map to compute µ and ℓ on the source surface; we then use (7.2) to construct a new Laplacian operator on M using µ and ℓ pulled back from N and show eigenfunctions of the resulting operator. Not only do the eigenvalues of the pulled-back Laplacian better approximate the Laplace-Beltrami eigenvalues of N , but qualitatively the eigenfunctions of the pulled-back Laplacian exhibit more structure in common with the eigenfunctions of N . Boxed examples in Figure 7.11 show particularly striking differences between the source and reconstructed eigenfunctions. Figure 7.12 illustrates an application of recovering edge lengths from truncated shape differences. Without constructing an embedding, we use pulled-back edge lengths to compute two commonly-used intrinsic functions: single-source geodesic distances and the wave kernel signature [START_REF] Aubry | The wave kernel signature: A quantum mechanical approach to shape analysis[END_REF]. Our edge lengths enable computation of these functions on the source mesh using the metric of the target, given a functional map between them. As a baseline, computing these functions on the target and pulling them back to the source using the functional map (right column) is less accurate; this is due to truncation of the functional map, which removes high frequencies e.g. at the center point of the geodesic function. .12 -Our technique can be used to recover the pullback metric and therefore compute geodesic distances without direct access to the target mesh. We compare three geodesic pullbacks on the source mesh: the exact pullback using point-to-point correspondence, the geodesic computed by reconstruction of the metric from the shape differences, and the pullback of the geodesic function using a functional map. For each pullback we compute the L 2 distance d to the exact version. Our method achieves better reconstructions than the direct usage of a functional map.

Reconstruction

Figures 7.13, 7.14 and 7.15 illustrate experiments in which geometry is reconstructed after estimating local structure from shape differences. To highlight our method's effectiveness on extrinsic motion, we show behavior on human shapes and cloth simulation data. Figure 7.13 applies our method to reconstructing models of humans from shape differences. From a coarse human base mesh (|X M | = |X N | = 502,; M = 100, k N = 200 in truncated experiments), we recover various poses. We compare reconstructions using only the intrinsic shape difference (right of each pair) to reconstructions using intrinsic and extrinsic differences together (left of each pair); we also compare using a truncated basis for shape differences (second column) to using a full basis (third column). As a baseline, we compare to [START_REF] Boscaini | Shape-fromoperator: Recovering shapes from intrinsic operators[END_REF], which uses only intrinsic geometry (rightmost column). Reconstruction from intrinsic information shows considerable artifacts due to the non-uniqueness of the solution of the embedding problem. Our provably complete intrinsic/extrinsic description is much more stable and close to the solution. The truncation of the basis, discussed in §7.8.2, tends to smooth out the sharp creases as they are represented as high frequency features. In 1089, k M = 100, k N = 200). After running a cloth simulation with coarse time steps, we compute the shape difference between subsequent frames (t ∈ [0, 1]). We then use the method in §7.7.3 to construct plausible motion between the frames by interpolating linearly between the computed shape differences (t = 0.5). We further extrapolate the motion beyond the t ∈ [0, 1] range to t = 1.5, effectively exaggerating the deformation between the frames. As expected, the extrinsic shape differences allow for reconstruction of largely isometric cloth motion. In this case, we reconstruct the same cloth simulation sequence but vary the topology of the source and target meshes. Now it is impossible to pull back the deformation exactly to the new mesh topology, but we still reconstruct plausible motion, with the notable exception of artifacts near the boundary of the patch.

Timings

Figure 7.16 shows timings by stage for our pipeline, applied to meshes of various sizes and topologies. We employ a simplistic single-threaded implementation in Matlab, using the Mosek toolbox [START_REF] Mosek Aps | The MOSEK optimization toolbox for MATLAB manual[END_REF] in the CVX library for convex optimization [START_REF] Grant | CVX: Matlab software for disciplined convex programming[END_REF]; for this 0. reason, the timings should be viewed as relatively pessimistic upper bounds. Even so, our implementation-including semidefinite constraints, regularization, and the like-can handle meshes of up to several thousand vertices.

Discussion & Conclusion

In this chapter, we introduced a new way to express intrinsic and extrinsic shape information through functional shape differences. Not only do we prove that discrete shape differences can be used to recover shape, but we also extend to characterizing shapes up to rigid motion rather than isometry. Our four shape differences together-two intrinsic and two extrinsic-comprise a powerful description of shape that applies to a wide range of variability, including not only non-isometric shapes but also models obtained from physical simulation and animation. We also show that the inverse problem of recovering shape structure from shape differences can be meaningful even in the under-determined truncated case. While this work offers the possibility of direct application in pipelines for shape search, embeddings of shape space, and approximate reconstruction, it also suggests myriad avenues for future research. On the theoretical side, a better understanding of the effect of Laplace-Beltrami eigenfunction truncation may provide better guidance for the minimal-sized shape differences needed to reconstruct a shape; spectral truncation is a common part of the geometry processing pipeline, so any relevant theory would have the potential to affect understanding of many existing algorithms. On the practical side, the primary limitation of our proposed reconstruction methods is the introduction of semidefinite constraints in computing the squared edge lengths ℓ; multi-scale or lighter-weight optimization methods would enable application to larger-scale meshes. Furthermore, the regularization proposed for recovery of µ and ℓ in §7.7 is very generic and can be ineffective for noisy or highly truncated shape differences. Application of machine learning techniques may allow for the characterization of edge length and triangle area distributions specific to a given class of shapes, considerably reducing the search space for our recovery algorithms. Figures 7.17 18 show examples illustrating these potential avenues for improving our pipeline. Figure 7.17 shows how reconstruction can fail when shape differences are over-truncated; stronger regularizers might fill in missing information when truncated shape differences are insufficient to recover edge lengths to high precision. Figure 7.18 shows results of shape reconstruction in the presence of noise. Here, we add noise directly to the shape difference matrix so that it no longer corresponds to an embedded surface. At some point, increasing k N does not improve the reconstruction result, because noise in the entries dominates added high-frequency shape information.

Chapter 8

Functional Characterization of Deformation Fields

In this chapter we present a novel representation for infinitesimal deformations of 3D shapes, by considering the induced changes in the underlying metric. Unlike traditional representations given, for example, by the infinitesimal displacement of each point, our framework allows to consider shape deformations as operators acting on real-valued functions defined on the shapes. This enables a wide variety of applications such as deformation design through precise control of metric distortion, deformation analysis and transfer and even improved shape matching by considering the composition of the deformation with other functional operators. Fundamentally, our approach helps to establish a direct connection between extrinsic deformation fields and changes in intrinsic metric quantities, which can be useful in many deformation processing scenarios.

Introduction

Designing and analyzing shape deformations is a central problem in computer graphics and geometry processing, with applications in scenarios such as shape manipulation [START_REF] Yu | Mesh editing with poisson-based gradient field manipulation[END_REF][START_REF] Sorkine | As-rigid-as-possible surface modeling[END_REF], animation and deformation transfer [START_REF] Sumner | Deformation transfer for triangle meshes[END_REF], shape interpolation [START_REF] Kilian | Geometric modeling in shape space[END_REF][START_REF] Von-Tycowicz | Real-time nonlinear shape interpolation[END_REF], and even anisotropic meshing [START_REF] Panozzo | Frame fields: Anisotropic and non-orthogonal cross fields[END_REF] among myriad others. Traditionally, shape deformation has been motivated by interactive applications in which the main goal is to design a deformation that satisfies some user-prescribed handle constraints while preserving the main structural properties of the shape. In other applications, such as shape interpolation and deformation transfer, that lack handle constraints, the goal is to design a global deformation field that would satisfy some structural properties as well as possible.

In both types of applications, most deformation methods are based on specifying a deformation energy and providing a method to optimize it. On the other hand, another very productive line of work has demonstrated that by choosing an appropriate representation for shape deformations, many tasks can become significantly easier, and in particular can help to enforce certain properties of the deformation field, which are otherwise very difficult to access and optimize for. In addition to the classical per-vertex displacement vectors, such representations have included gradient-based deformations [START_REF] Yu | Mesh editing with poisson-based gradient field manipulation[END_REF][START_REF] Zayer | Harmonic guidance for surface deformation[END_REF], Laplacian-based approaches [START_REF] Lipman | Differential coordinates for interactive mesh editing[END_REF][START_REF] Sorkine | Laplacian surface editing[END_REF] and Möbius transformations in the context of conformal deformations [START_REF] Crane | Spin transformations of discrete surfaces[END_REF][START_REF] Vaxman | Conformal mesh deformations with möbius transformations[END_REF] among others.

At the same time, several recent works have demonstrated that several basic operations in geometry processing can be represented as linear operators acting on real-valued functions defined on the shape. This includes the functional representation of mappings described in Chapter 4 or correspondences [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF][START_REF] Pokrass | Sparse modeling of intrinsic correspondences[END_REF], representations of vector fields as derivations [START_REF] Pavlov | Structure-preserving discretization of incompressible fluids[END_REF][START_REF] Azencot | An operator approach to tangent vector field processing[END_REF] and formulation of shape distortion via shape difference operators introduced in [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] and described in Chapter 4. One advantage of these representations is that linear operators can be naturally composed, which makes it easy to define, for example, the push-forward of a vector field with respect to a mapping, if both are represented as linear operators, or to solve for Killing vector fields, by composition between a derivation and the Laplacian operator.

While tangent vector fields are classically understood as operators (derivations) in differential geometry, extrinsic vector fields do not enjoy a similar property. Our main goal is to provide a coordinate-free representation of extrinsic vector fields (infinitesimal deformations) as functional operators, which will prove useful for analysis and design of shape deformations. As we demonstrate below our representation greatly simplifies certain tasks such as deformation transfer, deformation symmetrisation and even the computation of shape correspondences by composition with other operators. Moreover, it provides an explicit link between infinitesimal deformations and the changes in intrinsic metric quantities, which can be useful in a variety of analysis and deformation processing tasks.

For example, consider two shapes shown in Fig. 8.1 (left). By using our framework, it is possible to combine local extrinsic deformation constraints with intrinsic objectives such as constructing a deformation field that is as-isometric-as-possible (center). Moreover, our coordinate-free representation allows to relate deformations on multiple shapes, enabling deformation transfer and joint design even when only soft (functional) correspondences are known.

Related Work

Shape deformation is one of the oldest and best-researched topics in computer graphics and geometry processing. We therefore only mention works most directly related to ours and refer the interested reader to surveys including [START_REF] Nealen | Physically based deformable models in computer graphics[END_REF][START_REF] Botsch | On linear variational surface deformation methods[END_REF] and Chapter 9 in [START_REF] Botsch | Polygon Mesh Processing[END_REF]. A multitude of methods exists for surface deformation starting with the seminal work of [START_REF] Terzopoulos | Elastically deformable models[END_REF], its early follow-ups including [START_REF] Celniker | Deformable curve and surface finite-elements for free-form shape design[END_REF][START_REF] Welch | Variational surface modeling[END_REF] and the multi-scale variants, such as [START_REF] Zorin | Interactive multiresolution mesh editing[END_REF][START_REF] Kobbelt | Interactive multiresolution modeling on arbitrary meshes[END_REF][START_REF] Guskov | Multiresolution signal processing for meshes[END_REF] among many others. Similarly to our approach, many of these techniques are based on optimizing the so-called elastic thin shell energy that measures stretching and bending, and which is often linearized for efficiency. In the majority of cases, the deformation is represented explicitly as an extrinsic vector field defined on a surface, making deformation transfer difficult in the absence of precise pointwise correspondences. A number of methods have proposed alternative representations for deformation fields, which greatly simplify certain tasks in design and analysis. This includes gradient-based techniques [START_REF] Yu | Mesh editing with poisson-based gradient field manipulation[END_REF][START_REF] Zayer | Harmonic guidance for surface deformation[END_REF] which consider the deformation field by aligning its gradient with a set of local per-triangle transformations. By working in gradient space, constraints can be posed independently on the triangles and then optimized globally by solving the Poisson equation. Similarly, Laplacian-based techniques [START_REF] Sorkine | Laplacian surface editing[END_REF][START_REF] Lipman | Differential coordinates for interactive mesh editing[END_REF][START_REF] Nealen | A sketch-based interface for detail-preserving mesh editing[END_REF] are based on defining shape deformations by manipulating per-vertex differential coordinates (Laplacians) in order to match some target Laplacian coordinates. Such differential coordinates enable direct editing of the local shape properties, which can be especially beneficial for preserving and manipulating the high-frequency details of the surface. However, these coordinates are typically not rotationally invariant and additional steps are necessary to introduce invariance [START_REF] Sorkine | Laplacian surface editing[END_REF][START_REF] Lipman | Differential coordinates for interactive mesh editing[END_REF][START_REF] Paries | Simple and efficient mesh editing with consistent local frames[END_REF]. More recently, a number of methods have introduced representations for mesh deformations specifically geared towards particular shape manipulations, such as computing conformal transformations by designing special maps into the space of quaternions [START_REF] Crane | Spin transformations of discrete surfaces[END_REF] or by using face-based compatible Möbius transformations [START_REF] Vaxman | Conformal mesh deformations with möbius transformations[END_REF]. These techniques are rotationally invariant and coordinate-free, while being restricted to special types of manipulations. Another technique, closely related to ours, designs shape deformations by constructing a continuous divergence-free vector field [START_REF] Funck | Vector field based shape deformations[END_REF], and applying path line integration to obtain a deformed shape. We also consider the effect of the deformation on the metric, but both analyze the distortion of arbitrary extrinsic vector fields and show how they can be represented in coordinate-free way as linear functional operators. Our approach of considering the deformation via its induced metric distortion is also related to the work of [START_REF] Eigensatz | Positional, metric, and curvature control for constraint-based surface deformation[END_REF] and [START_REF] Sela | Computational caricaturization of surfaces[END_REF] who manipulate shapes by explicitly editing their curvature properties. Moreover, our use of the strain tensor in characterizing metric distortion is closely related to the applications in various physically based deformation scenarios including [START_REF] Thomaszewski | Continuum-based strain limiting[END_REF][START_REF] Müller | Strain based dynamics[END_REF] among many others (see also the surveys on physically based elastic deformable models [START_REF] Nealen | Physically based deformable models in computer graphics[END_REF][START_REF] Rumpf | Geometry processing from an elastic perspective[END_REF]). Our approach is also related to the works that aim to design as-isometric-as-possible shape deformations [START_REF] Zhang | Fast as-isometricas-possible shape interpolation[END_REF][START_REF] Solomon | As-killing-aspossible vector fields for planar deformation[END_REF][START_REF] Martinez Esturo | Generalized metric energies for continuous shape deformation[END_REF]. Similarly to the latter work, however, our framework is general and allows an arbitrary prescribed distortion, although our method works directly on surface representations and moreover enables applications such as joint deformation design.

Finally, our joint deformation design is related to the deformation transfer and interpolation techniques such as [START_REF] Sumner | Deformation transfer for triangle meshes[END_REF][START_REF] Baran | Semantic deformation transfer[END_REF] and [START_REF] Kilian | Geometric modeling in shape space[END_REF] to name a few. However, we place special emphasis on relating infinitesimal deformations between shapes with only soft (or functional) correspondences, which are often much easier to obtain than detailed point matches.

In contrast to the majority of existing techniques our goal is to devise a coordinate-free representation of infinitesimal deformations as linear functional operators, by making an explicit connection between the extrinsic deformations and the change in intrinsic metric quantities. Thus, although we build on classical constructions such as the infinitesimal strain tensor, we show how they can be exploited to create a functional representation of shape deformation, which can be used in conjunction with other operators. As we demonstrate below, our representation is particularly useful for relating deformations that exist on multiple shapes, with only soft correspondences between them. In particular, it enables joint deformation design and transfer without triangle or point-to-point matches, and allows to introduce extrinsic information in the computation of functional maps, greatly increasing the applicability of these techniques in a wide range of application scenarios.

Overview

Our main strategy for devising a functional representation of infinitesimal deformation is to consider the effect of the deformation on the intrinsic shape metric. Namely, we will characterize the extrinsic vector fields by the isometric distortion that they induce. For this purpose we will start by considering the shape difference operators described above, which provide a way to capture the isometric distortion induced by a map between two shapes as a pair of linear operators acting on functions on one of the shapes. Thus, in this chapter our motivation is to use the information contained in those operators to design and synthesize shape deformations. The first challenge presented by this approach is the fact that the shape difference operators only contain intrinsic information which is typically not enough to specify an embedding, and thus a deformation. To tackle this problem we restrict our attention to infinitesimal deformations which can be represented by extrinsic vector fields, and which, as we show under certain genericity conditions, are fully encoded by the isometric distortion they induce. Second, the definition of the shape difference operators introduced in [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] is divided between an area and conformal-based distortions. While this distinction might be desirable for deformation analysis, it is cumbersome for deformation synthesis since it would require us to represent each extrinsic vector field via a pair of functional operators. For this reason we introduce a unified shape difference operator, that characterizes intrinsic distortion fully. Finally, although the shape difference operators characterise intrinsic distortion, designing a deformation that would reproduce a given shape difference leads to challenging optimization problems as shown in [START_REF] Boscaini | Shape-fromoperator: Recovering shapes from intrinsic operators[END_REF]. As we demonstrate below, our approach of considering infinitesimal deformations greatly simplifies the reconstruction procedure, by allowing us to phrase it with a single linear system of equations. To summarize, we proceed in three stages: first we introduce a unified shape difference operator that fully characterizes isometric distortion. Second, we show how extrinsic vector fields can be represented by infinitesimal shape difference operators, which act linearly on the real-valued functions, in the same way as shape difference operators. Finally, we demonstrate that extrinsic vector fields can be recovered from infinitesimal shape difference operators by solving a linear system of equations. Our main contribution therefore include:

• Defining a unified shape difference operator that fully characterizes isometric distortion.

• Introducing infinitesimal shape difference operators as a way to represent extrinsic vector fields as operators acting on functions, represented as matrices in the discrete setting.

• Showing how infinitesimal shape difference operators can be used to naturally add extrinsic information (second fundamental form) into the computation and analysis of functional maps.

• Presenting various applications that demonstrate the usefulness of our representation, in particular by showing how the infinitesimal deformations can be analysed, designed, transfered and used in map computation by exploiting the operator representation.

Extrinsic Vector Fields as Operators

In this section we provide a coordinate-free representation of extrinsic vector fields by considering their action on the underlying shape metric. As mentioned in the previous section, our main strategy will be to represent extrinsic vector fields via infinitesimal shape difference operators. Therefore, we first define a unified isometric operator in Section 8.4.1. We then define infinitesimal shape difference operators in Section 8.4.2 and list some of their key properties in Section 8.4.3. Throughout this section we assume that we are dealing with smooth surfaces without boundary embedded in R 3 . The appropriate discretization of all the concepts introduced in this section will be given in Section 8.5.

Isometric Shape Difference Operator

While the two shape difference operators defined in Chapter 4 are very convenient to separate measure (area) and conformal distortion, they also imply that two operators are needed to fully characterize a deformation, which will be cumbersome for representing extrinsic vector fields. Moreover, these two operators are not commensurable since D A is defined through products of functions, unlike D C which is a differential operator defined through inner products of gradients. Thus, we introduce a single unified shape difference operator below.

Unified shape difference The main reason for which D C is only sensitive to conformal changes is that both the inner product and the integration are taken on the target shape.

To see this, let us re-write the conformal shape difference D C by integrating on a fixed domain using the pullback metric tensor as follows:

f, D C (g) M H 1 0 = M C ϕ -1 ( ∇C ϕ (f ), ∇C ϕ (g) ) d(ϕ ⋆ µ N ).
This definition is equivalent to the one given in Eq. (4.3) in Chapter 4, but here we simply highlight the integration with respect to the pushforward measure d(ϕ ⋆ µ N ) rather than the measure on M itself. Here, the linear operator C ϕ -1 is simply the functional map with respect to the inverse diffeomorphism ϕ -1 . Since we are dealing with surfaces, the change in the area will cancel out on the right, leaving D C only sensitive to conformal changes.

To define a unified shape difference taking into account all intrinsic changes one should compare the pullback metric to the metric on M while keeping the integrating measure fixed. We thus propose a unified shape difference operator D I that fully characterizes isometric distortion.

Definition 8.4.1 Assuming that ϕ : N → M is a diffeomorphism, the unified shape difference D I : H 1 0 (M ) → H 1 0 (M ) is defined implicitly by:

f, D I (g) M H 1 0 := M C ϕ -1 ( ∇C ϕ (f ), ∇C ϕ (g) ) dµ M .
The existence of D I is once again guaranteed by the Lax-Milgram theorem in H 1 0 . Moreover, as we claimed above, the following proposition shows that the unified shape difference fully characterizes isometric deformation. Proposition 8.4.2 The following equivalence holds:

• D I (f ) = f for all f ∈ H 1 0 (M ) if
M ∇f, ∇D I (g) dµ = M C -1 ( ∇C(f ), ∇C(g) ) dµ = M ∇f, A∇g dµ
where A is a symmetric (1,1)-tensor defined as A ij = ((ϕ -1 ) ⋆ g N ) ik g M kj . This definition allows us to prove both equivalences.

• If ϕ is an isometry then (ϕ -1 ) ⋆ g N = g M so M ∇f, ∇D I (g) dµ M = M ∇f, ∇g dµ Therefore D I (g) = g for all function g in H 1 0 . If D I (f ) = f then M ∇f, ∇(Id -A)g dµ = 0 for all f, g ∈ H 1 0 (M ). Lemma 4.2.2 tells us that ϕ is an isometry.

• If ϕ is an area-preserving then det((ϕ -1 ) ⋆ g N ) = det(g M ) so the tensor field A is equal to the tensor field B as defined in Equation 4.6.

If D I (f ) = D C (f ) then, using Lemma 4.2.2, the equality holds between the tensor A and B. We conclude that det((ϕ -1 ) ⋆ g N ) = det(g M ).

This new shape difference is defined up to isometric deformation and a pullback operation can be defied to transport a shape difference on another manifold. More precisely, Proposition 4.2.4 can be restated for the unified shape difference.

I : H 1 0 (M ) → H 1 0 (M ) associated to ϕ • φ : N → M is given by D ϕ•φ I = D ϕ I • C -1 ϕ • D φ I • C ϕ .
Proof The proof relies only on Definition 8.4.1:

P C ϕ ∇f, ∇D ϕ•φ I (g) dµ = P C -1 φ ( ∇(f • ϕ • φ), ∇(g • ϕ • φ) ) dµ = P ∇(f • ϕ), ∇D φ I (g • ϕ) dµ = P C ϕ ∇f, ∇D ϕ I D φ I (g • ϕ) • ϕ -1 dµ.
This yields the equality

D ϕ•φ I (g) = D ϕ I D φ I (g • ϕ) • ϕ -1 for all g ∈ H 1 0 (M ).
Example Computing D I is challenging when studying meshes with different connectivity as it requires the map and its inverse. However, when using it to define an operator representation for infinitesimal deformations, we will only be dealing with shapes with a fixed connectivity, for which we are able to obtain an appropriate discretization. Moreover, even in the general case, we are able to approximate this quantity in order to compute it using the functional map framework, as shown in Section 8.5.

To illustrate the properties of the three shape differences we use a simple low-dimensional description of a shape collection. Here we choose a fixed base shape and compute the shape difference matrices with respect to the remaining shapes in a collection. Then, we represent each shape by its shape difference matrix and plot them as points in PCA space. Figure 8.2 represents the conformal deformation of a bunny into a sphere as viewed by the three shape differences. As expected the conformal shape difference is almost identity while the area and isometric shape differences both capture the distortion. In the second experiment, shown in Figure 8.2, we explore another collection obtained by the shearing of a plane patch. As this deformation is area preserving, the area-based shape difference provides no information, unlike the other two shape differences which are equal in this case.

Infinitesimal Deformations as Operators

One aspect of the shape difference operators, highlighted in [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF], that makes them particularly useful for capturing and representing shape deformations is that both the The PCA applied to shape differences confirms the presence of large area and isometric distortion in contrast to small conformal distortion. Bottom row: Area preserving deformation of a plane (top row right). The area shape difference is almost constant while conformal and isometric differences agree.

domain and the range of these operators are functions on a fixed base shape. In other words, although they are obtained by considering the distortion with respect to some target shape, the final functional operators both act on and produce real-valued functions on the base shape. Our main goal in this section is to consider a one-parameter family of shapes M t , given by displacing the points of a base shape along some fixed deformation field. This family of shapes will result in a one-parameter family of shape difference operators (D I ) t on the base shape, that we will use to define the infinitesimal shape difference operator ∂(D I ) t /∂t at t = 0, as a coordinate-free representation of the deformation field.

Deformations via Extrinsic Vector Fields Specifically, given a surface M embedded in R 3 we consider a family M t of deformations of M , parameterized by a scalar t and isometrically immersed by the local mappings

F t : U ⊂ R 2 → M t ⊂ R 3 .
This family of manifolds is generated by the displacement of the points along the smooth vector field

V (p) ∈ T p M × T p M ⊥ ≃ R 3 : ∂F t ∂t (p) = V (p), (p, t) ∈ M × R + , (8.1)
where F 0 is the local immersion of M . We call V a displacement field, or, alternatively, an extrinsic vector field to highlight the distinction with tangent vector fields.

Our main goal below is to characterize extrinsic vector fields as linear functional operators. For this, as mentioned above, we consider the family of diffeomorphisms ϕ t : M t → M , given trivially via ϕ t (p) = F t (p) -tV (p), and the associated functional maps C ϕt mapping functions from M 0 to M t . This creates a one parameter family of shape difference operators D V t (which can be taken either to be the area or conformal-based operators or the unified shape difference operator defined above). We then introduce the infinitesimal shape difference, as a functional representation of the deformation field V as follows: Definition 8.4.4 The functional representation of an extrinsic vector field V on a surface M is given by the functional operator E V :

E V := ∂D V t ∂t t=0 (8.2)
We will call E V the infinitesimal shape difference operator associated with the deformation field V .

We will use E V A , E V C and E V I to denote the infinitesimal shape difference operators arising from the area-based, conformal and unified shape differences respectively. Note that although we characterize the properties of all three infinitesimal shape difference operators, in all practical applications we will only use E V I as the functional representation of an extrinsic vector field V . This is because this operator allows us to represent an extrinsic vector field using a single linear operator, and moreover, as we show below, under certain genericity conditions, in the discrete case, E V fully characterizes V up to rigid motion. Remark that since E V is defined as a derivative of a one-parameter family of linear operators acting on real-valued functions on a surface, both the range and the domain of E V are also real-valued functions on M . As such it can naturally be represented in a multi-scale basis, composed with other functional operators, and analyzed thus enabling applications including design, processing and transfer of extrinsic vector fields, as well as their use in other tasks such as improved map (correspondence) computation. Moreover, as we demonstrate below E V I is very closely related to classical constructions in computer graphics, geometry processing and mechanics, based on the infinitesimal strain tensor.

Main Properties of Infinitesimal Shape Differences

Below we provide the defining characteristics of the infinitesimal shape difference operators E V , which will be useful in both deriving the appropriate discretization in the case of shapes represented as triangle meshes, and in the applications shown below.

The Levi-Civita Covariant Derivative We first need to introduce some fundamental notions from differential geometry. In particular, we will use the Levi-Cevita connection to define derivatives on a surface. For this, consider a tangent vector u at some point p ∈ M , and an extrinsic vector field V on M . Then, given an arbitrary curve γ(t) on M such that γ(0) = p and γ ′ (0) = u, we obtain ∇u V = ∂V (γ(t)) ∂t t=0

. Here ∇u V is the covariant derivative of the ambient space. Note that at a fixed point p ∈ M , ∇u V is a vector in R 3 . We can project the covariant derivative onto the tangent plane at p to obtain a vector in the tangent plane, which we denote simply by ∇ u V where ∇ is the Levi-Cevita connection on M extended naturally to extrinsic vector fields, ( [START_REF] Carmo | Riemannian Geometry, Mathematics: Theory & Applications[END_REF] p. 126), using the local coordinates ∇ i V j = ∂ i V, ∂ j F . We also remark that for any vector x in the tangent space, ∇ u V, x = ∇u V, x , which we will use in our discretization. This definition allows us to express the divergence of the extrinsic vector field V as the trace of the Levi-Civita connection, given, at a fixed point p by div(V ) = g ij ∇ i V j .

The fundamental object that we consider below is the infinitesimal strain tensor, which can be understood as a bilinear form, acting on pairs of vectors x, y in the tangent plane of a point p ∈ M . Namely, given an extrinsic vector field V , the infinitesimal strain tensor L V g(x, y) is defined as:

L V g(x, y) = x, ∇ y V + ∇ x V, y (8.3)
This quantity has the advantage of being linear in the vector field V which makes it easy to handle for deformation and vector field design and therefore has been used in a wide variety of works in computer graphics [START_REF] Nealen | Physically based deformable models in computer graphics[END_REF]. With these definitions in hand we can express the fundamental properties of the infinitesimal shape difference operators described above.

Proposition 8.4.5 Given a one parameter family of surfaces described in Eq. (8.1), for a fixed point p, the first-order change in the metric tensor g and in the local area element µ = det(g) are given as:

∂g(t) ∂t t=0 = L V g (8.4) ∂µ(t) ∂t t=0 = div(V )µ. (8.5)
Proof Those properties are easily proven when using local coordinates. Given a family of diffeomorphisms ϕ t the Lie derivative of the metric tensor with respect to the vector field V denoted L V g is by definition:

L V g := ∂ ∂t (ϕ -1 t ) ⋆ g(t) t=0 .
Since the local immersion F t use a common chart system, the coordinates of the pullback metric

((ϕ -1 t ) ⋆ g t ) ij are equal to the metric on M t in local coordinates g ij (t) = ∂ i F t , ∂ j F t .
The computation of derivative is then straigthforward:

∂ ∂t ((ϕ -1 t ) ⋆ g(t)) ij t=0 = ∂ ∂t ( ∂ i F t , ∂ j F t ) t=0 = ∇ i V j + ∇ j V i
From here Eq. 8.5 is obtained easily:

∂µ(t) ∂t t=0 = ∂ ∂t det(g(t)) t=0 = 1 2µ det(g)g ij (∇ i V j + ∇ j V i ) = div(V )µ
This property makes apparent the fact that the strain tensor is equal to the Lie derivative of the metric tensor with respect to the extrinsic vector field V , and that the first-order change in the local area is controlled by the divergence of the vector field, as is similarly the case for tangent vector fields.

Infinitesimal Shape Differences As mentioned in Section 8.4.1, shape differences operators can be expressed using the pullback metric and pushforward measure. Using Proposition 8.4.5 we have all the tools to describe the infinitesimal shape difference operators. The following proposition characterizes these new operators.

Proposition 8.4.6 Let V be a smooth deformation field on M , the derivatives of D A , D C and D I at time zero satisfy for all smooth functions f, g:

f, E V A (g) M L 2 = M div(V )f g dµ, f, E V C (g) M H 1 0 = M div(V ) ∇f, ∇g -L V g(∇f, ∇g) dµ, f, E V I (g) M H 1 0 = - M L V g(∇f, ∇g) dµ.
Proof The first statement is obtained by using (8.5):

f, ∂ t D A (g) M L 2 = ∂ ∂t Mt C t (f )C t (g) dµ(t) t=0 = ∂ ∂t M f g det(g(t)) det(g) dµ t=0 = M div(V )f g dµ.
The second statement is a consequence of Equation (4.6):

f, ∂ t D C (g) M H 1 0 = ∂ ∂t Mt ∇C t (f ), ∇C t (g) dµ(t) t=0 = ∂ ∂t M ∂ i f g ij (t)∂ j g det(g(t)) det(g) dµ 0 = M div(V ) ∇f, ∇g dµ - M ( ∇f, ∇ ∇g V + ∇ ∇f V, ∇g ) dµ
Starting from Definition 8.4.1:

f, ∂ t D I (g) M H 1 0 = ∂ ∂t M C -1 t ( ∇C t (f ), ∇C t (g) ) dµ t=0 = ∂ ∂t M ∂ i f g ij (t)∂ j g dµ t=0 = - M ∇f, ∇ ∇g V + ∇ ∇f V, ∇g dµ.
As can be seen, the infinitesimal shape differences inherit the properties of the original operators. Namely, E V A vanishes if and only if div(V ) is equal to zero, i.e., whenever V infinitesimally preserves the volume form. On the conformal side, finding an extrinsic vector field V such that E V C = 0 is equivalent to solving the conformal Killing equation: L V g = div(V )g characteristic of infinitesimal conformal vector field. Both properties combined lead to an infinitesimal isometric deformation induced by the vector field V captured by Eq. (8.4). Moreover Prop. 8.4.6 reveals a clear link between shape differences:

f, E V I (g) M H 1 0 = f, E V C (g) M H 1 0 -1, E V A ( ∇f, ∇g ) M L 2 . (8.6)
Thus, intuitively, the operator E I , representing isometric distortion, can be decomposed into an area and conformal parts.

The operator representation of extrinsic vector fields, allows us to derive a some interesting properties specific to this representation and useful for applications, while making a connection between previously proposed shape difference operators and the classical strain tensor.

Linearity The first key property of E I is its linearity with respect to both vector fields and functions:

E αU +βV I (f ) = αE U I (f ) + βE V I (f ), E V I (αf + βg) = αE V I (f ) + βE V I (g). (8.7)
This makes it easily usable for vector field design, by representing the target deformation field in some fixed basis.

Lie derivative representation From Proposition 8.4.2 we deduce that E V I (f ) = 0 for all f if and only if L V g = 0, i.e., if the extrinsic vector field V preserves the intrinsic shape metric to first-order. The operator E V I clearly quantifies how the Lie derivative of the metric affects gradients of functions. So the linear dependence between shape operators shown in Eq. (8.6) can be understood as the decomposition of the matrix L V g into a trace free part, linked to the conformal Killing equation, and a divergence part, related to the change in area.

Vector field representation In general the operator E V I does not uniquely define an extrinsic vector field. From Prop. 8.4.2 the kernel of V → E V I coincides with the vector fields satisfying L V g = 0. In case of a volumetric manifold (i.e. M ⊂ R 3 ) the infinitesimal change in the metric ∂ t g defines the vector field V up to rigid motion (see [27] and all the Korn's inequalities). In the case of a surface embedded in R 3 the kernel of E V I includes infinitesimal isometries such as Killing vector fields but also local normal fields in planar areas. However such fields seem either rare or nonexistent generically. No rigidity result seems to be known for smooth surfaces. However, as we demonstrate below, in the discrete case it is can be shown that for almost all surfaces the kernel of V → E V I consists only of rigid deformations (Prop. 8.5.6).

Second-fundamental form representation One interesting special case to consider is the interpretation of E I when the deformation field is the normal field n. By using Eq. (8.3) it is possible to see that the covariant derivative of the normal yields the second fundamental form denoted by h p : T p M × T p M → R, more precisely L n g = 2h. Therefore the operator E n I captures the action of curvature on functions, since:

M ∇f, ∇E n I (g) dµ = -2 M h(∇f, ∇g)dµ.
From a theoretical point of view the knowledge of the Laplace-Beltrami operator gives access to the first fundamental form and E n I yields information about the second. Thus these two operators jointly provide a coordinate-free representation of the embedding. The operator E n I can be used to obtain a multi-scale representation of curvature information on the triangle mesh, as shown in Figure 8.3. In particular, the eigenfunctions corresponding to the largest eigenvalues of E n I , are those that align the best with the maximal principal curvature direction, and can be obtained even if E n I is represented in a reduced functional basis, making the computation less sensitive to noise in the triangulation. Moreover, as we demonstrate in Section 8.7.2, the operator E n I can be used to inject extrinsic information into the computation of functional maps.

Composition with mappings In many applications, we are interested in the relation between deformation on multiple surfaces related by a mapping. In particular given a deformation field U N of shape N and a diffeomorphism ϕ : N → M with the associated functional map (pullback) C ϕ of functions from M to N , one can define an infinitesimal deformation V M of shape M that produces the same metric distortion. Instead of looking directly at the deformation of the metric, which might require a mapping between individual triangles [START_REF] Sumner | Deformation transfer for triangle meshes[END_REF], we account for the action of the metric on functions:

E U N I C ϕ (f ) = C ϕ E V M I (f ) ∀f ∈ H 1 0 (M )
In other words, V M can be obtained by considering an extrinsic vector field, whose operator representation has the same effect on functions when composed with the functional Figure 8.4 -Example of deformation fields that commute with the underlying diffeomorphism. Note that both the direction and the magnitude of the vector field have to adapt to the underlying geometry to produce the same metric distortion.

map C ϕ as E U N I . This property allows us to transfer deformation fields without requiring point-to-point correspondences between shapes, by simply considering the commutativity of the operators C ϕ and E I . We illustrate it in Figure 8.4 and use it in Section 8.7.1 for deformation transfer and deformation symmetrization on meshes with different connectivity with only a functional map known between them. Furthermore, this approach is applicable to design deformations jointly on two shapes, such that they are consistent with the functional map C ϕ and even as a regularizer in map computation.

Discrete Setting

In this section we provide the discretization of the concepts introduced above. For this, we first define a discretization of the unified shape difference by considering Definition 8.4.1.

We then obtain a discrete version of the infinitesimal shape difference, by considering a particular discretization of the Levi-Civita connection and the Lie derivative of the metric on the triangle mesh. Finally, we demonstrate that the property stated in Proposition 8.4.6 which links the Lie derivative of the metric and the infinitesimal shape difference is satisfied exactly in the discrete case. Throughout this section, we assume that we are given a manifold triangle mesh. We denote by (X , E, F) respectively the set of vertices, edges and faces. To compute the shape differences we start from the discretization of the inner product ., . H 1 0 using standard first order finite elements. We will denote by L the classical cotangent Laplacian matrix, W the inner product of H 1 0 and A the lumped mass matrix such that L = A -1 W . As before µ is a measure and µ(T ) denotes the area of triangle T .

Discrete unified shape difference

The discretization of the unified shape difference is straightforward when N and M are triangle meshes and share the same connectivity. In Definition 8.4.1 given above, the gradients and the point-wise scalar products are taken on N while the measure dµ M comes from M . Therefore the right hand side can be discretized by a modified cotangent weight formula:

W M D I = W M N ,
where

(W M N ) i,j = 1 2 µ M (T α ) µ N (T α ) cot α N ij + µ M (T β ) µ N (T β ) cot β N ij . (8.8)
Here T α , T β are the two triangles adjacent to edge i, j, which is opposite to angles α and β, while µ M and µ N are the triangle areas on shapes M and N respectively. Note that W M N differs from the standard cotangent weight matrix W N only by the ratio of weights per triangle. Moreover, one can remark that if the transformation is area preserving for all triangles then D I reduces to the conformal shape difference defined in [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] (Option 1 in Section 5). From the expression above it follows that D I = W -1 M W M N . For comparison the conformal shape difference has a close formulation D C = W -1 M W N . Remarkably, Theorem 4.2.1 admits a discrete parallel as the operator D I entirely encodes the discrete metric. Proof In [START_REF] Zeng | Discrete heat kernel determines discrete Riemannian metric[END_REF] the authors prove that W uniquely determines the edge length up to global scaling. Our proof simplifies and adapts their arguments to our problems. We will first build a convex energy whose gradient corresponds to the modified cotangent weights then show that this energy have a unique minimizer. Let C be the set of squared edge length respecting the triangle inequality, namely the set:

C = x ∈ R |E| : x ij ≥ 0, ∀(i, j) ∈ E ; √ x ij ≤ √ x jk + √ x ki , ∀(i, j, k) ∈ F .
As shown by Proposition 7.7.1, C is convex. Thanks to the Heron's formula the triangle areas can be expressed as a function of the squared edge lengths. More precisely we have:

µ(T ) 2 = 1 16   ℓ 2 ij ℓ 2 jk ℓ 2 ki   ⊤   -1 1 1 1 -1 1 1 1 -1     ℓ 2 ij ℓ 2 jk ℓ 2 ki   , ℓ 2 ∈ C.
Therefore the gradient of µ(T ) with respect to the squared edge lengths is equal to the vector of cotangent:

∇ ℓ 2 µ T (ℓ 2 ) = 1 4   cot(α ij ) cot(α jk ) cot(α ki )   , ℓ 2 ∈ C,
We recognize here the cotangent weight Laplacian L N with lumped area weights, namely A M L N . In the case of meshes with different connectivity, this remark suggests the following approximation of the isometric shape difference, valid only in a discrete sense, for an arbitrary linear functional map C between M and N :

f ⊤ A M L M D I g ≈ f ⊤ A M C -1 L N Cg.
In the reduced basis of the Laplacian eigenvectors, the approximation of the shape difference becomes D I ≈ Λ -1 M C -1 Λ N C, which preserves the principal property of the operator: D I is identity if and only if the deformation is an isometry since the Laplacian on N has to be equal to the Laplacian on M . We used this discretization in Figure 8.2 and observed that the two expressions given above typically produce similar results.

Infinitesimal shape difference

Throughout this section we will consider the deformation field, which we also call an extrinsic vector field, to be a given as a three-dimensional vector per vertex. As we mentioned above, there are two possible ways of discretizing the infinitesimal shape difference operator. First, one can use Definition 8.4.1 and take the derivative at time zero of a discrete shape difference operator. The second possibility is to start from Proposition 8.4.6 and find the discrete equivalent of the Lie derivative of the metric. Remarkably, those two paths lead to the same discretization in our case.

Discrete connection

To build the discrete operator E V I we need a consistent discretization of the Levi-Civita connection. While several discrete connections have been proposed (e.g. [START_REF] Azencot | Discrete derivatives of vector fields on surfaces-an operator approach[END_REF][START_REF] Liu | Discrete connection and covariant derivative for vector field analysis and design[END_REF]), because of the special nature of our problem, we choose to build our own. This is because, in applications such as parallel transport it is fundamental that the vectors u, v and ∇ u v are expressed in the same space (at vertex or face or edge) so often an averaging step has to be introduced to transfer, for example, a face-based representation of a vector to an edge based representation. In our setting such a requirement is not needed and it is easier to distinguish tangent vector fields that will be expressed by one vector per face and extrinsic vector fields expressed at vertices. Thus, our goal is to obtain a connection of the ambient space ∇u V where u is a tangent vector and V is an extrinsic vector field :

∇ : R 3|F | × R 3|V| → R 3|F | (u, V ) → ∇u V
We build the connection ∇ using finite differences as follows. Since extrinsic vector fields are defined at vertices the differences are taken along the edges.

Definition 8.5.2 In a given triangle T ∈ F the ambient covariant derivative along the edge e ij is defined by

∇ e ij e ij V T = V i -V j e ij .
Thus the ambient connection in the directions E = (e ij , e jk ) can be stored in a matrix

( ∇E V ) T = V i -V j V j -V k .
Then, given any tangent vector x = Eα, the covariant derivative in its direction can be computed as ∇x V = ( ∇E V )α.

Given the expression above, the discrete Lie derivative of the metric inside triangle T follows immediately, using Eq. (8.3). Namely for any pair of tangent vectors x = Eα, y = Eβ in the triangle T , we have:

L V g(x, y) T = x, ( ∇E V )β + ( ∇E V )α, y . (8.9)
After integration we obtain the discrete infinitesimal shape difference:

f ⊤ W M E V I g = - T ∈F L V g(∇f, ∇g) T µ(T ).

Shape difference derivative

The above discrete formulation can also be considered from another point of view. Recall that we first introduced the infinitesimal shape difference as a derivative of a one-parameter family of shape differences. Interestingly, it is possible to reproduce the equivalence shown in Proposition 8.4.6 in the discrete setting using our discretization of the unified shape difference operator. Suppose that each vertex p i of the mesh is displaced by the vector V i by p t i = p i + tV i . This produces a family of triangle meshes (X t , E, F) with identical connectivity. Using this setup properties similar to the continuous case can be derived. Proposition 8.5.3 Given a one parameter family of meshes, the first-order change in the metric tensor g T = E ⊤ E and in the area at a triangle T ∈ F, is given as:

∂ t g T ∂t t=0 = E ⊤ (∇ E V ) T + (∇ E V ) ⊤ T E, ∂ t µ(T ) ∂t t=0 = Tr g -1 T E ⊤ (∇ E V ) µ(T ).
Proof First let's remark that the derivative of the discrete metric can be expressed with respect to discrete connection:

∂g t T ∂t t=0 = E ⊤ (∇ E V ) T + (∇ E V ) ⊤ T E.
Since the metric is linked to the triangle area by µ t (T ) = 1 2 det(g t T ) the statement obtained by a direct computation of the derivative.

In the continuous case the divergence was defined using orthonormal tangent vectors. The formula above is the analog for a pair non orthonormal vectors, therefore we will denote div(V ) T = Tr g -1 T E ⊤ (∇ E V ) . Taking the derivative of the discrete unified shape difference in Eq. 8.8 might be challenging. However, using the formulation of Eq. (3.7) in Chapter 3 leads to an equivalent formulation of Eq. (8.8) is:

f ⊤ W D t I g := T ∈F ∇ t f, ∇ t g t T µ(T ). (8.10)
By taking the derivative of his expression at time t = 0, we obtain an alternative discretization of the infinitesimal shape difference E I .

Proposition 8.5.4 The discrete infinitesimal shape difference reads

E V I (u) = W -1 M H
, where H is a Laplacian matrix whose weights depend on the extrinsic vector field:

(H) ij = 1 2 j∼i (c(T α ij ) + c(T β ij )), c(T ) = ( e jk , V j -V i + e ij , V j -V k ) 1 4µ(T )
-div(V ) T e jk , e ki µ(T ) .

Proof Using the FEM gradient, e.g. Equation (3.6), to discretize the unified shape difference written in Eq. 8.10 leads to:

f ⊤ W D t I g = T ∈F 1 4 f j -f k f j -f i ⊤ g t T µ t (T ) 2 g j -g k g j -g i µ(T ). (8.11)
We can now compute the derivative with respect to time by using Proposition 8.5.3:

f ⊤ W E V I g = T ∈F 1 4µ(T ) f j -f k f j -f i ⊤ L T g j -g k g j -g i , L T = E ⊤ (∇ E V ) T + (∇ E V ) ⊤ T E -2Tr g -1 T E ⊤ (∇ E V ) g T .
The matrices L T can be written in a form similar to the discrete metric (see Prop. 3.2.2):

L T = 1 2 2a ij a ki -a jk -a ij a ki -a jk -a ij 2a jk , where a ij = e ij , V i -V j -2div(V ) T l 2 ij , (8.12)
leading to the point-wise formulation:

(W M E V I ) ij = 1 2 j∼i (c(T α ij ) + c(T β ij )), c(T ) = -a ki + a jk + a ij 4µ (

T ) .

Remarkably these two discritizations are strictly identical.

Proposition 8.5.5 The discretization of E I based on the discrete Levi-Civita connection is equivalent to the one obtained by differentiating the unified shape difference operator.

Proof In Equation (8.9) the tangent vectors in a given triangle have be to expressed in the basis form by two edges of the triangle. Following the discussion in Section 3.2.3, the FEM gradient at a face T can be written in two equivalent ways:

∇f T = 1 2µ(T ) R 90 • E f j -f k f j -f i = Eg -1 T f j -f i f k -f j .
Therefore the discrete strain tensor at triangle T follows immediately:

L V g(∇f, ∇g) = f i -f j f j -f k ⊤ g -1 T ( ∇E V ) + ( ∇E V ) ⊤ g -1 T g i -g j g j -g k ,
From Proposition 8.5.3, we recognize the term the derivative of the inverse metric. Using Equation (3.10), one can further modified the expression to:

L V g(∇f, ∇g) T = - 1 4 f j -f k f j -f i ⊤ ∂ ∂t g t T µ t (T ) 2 t=0 g j -g k g j -g i . (8.13)
The right had side term appears in the discrete isometric shape difference as written in Equation (8.11). This leads to the equality between the different discretization:

∂ ∂t f ⊤ W D t I g t=0 = - T ∈F L V g(∇f, ∇g) T µ(T ).
Interestingly, the discrete metric tensor can be defined for higher dimension simplicial complexes in accordance with the FEM discretization. Therefore Proposition 8.5.5 holds in more general cases.

Properties

Interestingly, many of the properties of the continuous operators are satisfied exactly by their discrete counterparts.

Linearity

The discretization E V I (f ) naturally preserves the linearity with respect to both V and f . Shape difference decomposition Since the discretization using a discrete connection and through a the time derivative agree, the decomposition described by Eq. (8.6) is also satisfied exactly. With the family of meshes conformal shape difference is:

f ⊤ W D t C g := f ⊤ W t g = T ∈F ∇ t f, ∇ t g t T µ t (T ),
where W t is the cotangent weight matrix associated to the deformed mesh at time t. As noted above the derivative of ∇ t f, ∇ t g t T will lead to E V I . The derivative of the area will produce a term close to M div(u) ∇f, ∇g dµ thanks to Proposition 8.5.3. Namely:

f ⊤ W E V C g = f ⊤ W E V I g + T ∈F div(V ) T ∇f, ∇g T µ(T ).
Thus, the decomposition of E V I , representing isometric distortion, into area and conformal parts given in Eq. (8.6) in the continuous case holds exactly in the discrete case as well.

Vector Fields representation In the continuous setting the kernel of V → E V I is the set of infinitesimal isometries. However, to the best of our knowledge, there is no characterization of how often this set is reduced to rigid motion. In the particular setting of our discretization some standard results can be applied, however. Proposition 8.5.6 For almost all triangle meshes M without boundary, the operator E V I uniquely defines the extrinsic vector field V up to rigid motion. Proof The proof is organized as follows: we show that we can recover the matrices L T from the infinitesimal shape difference then we use a standard results in combinatorix to prove that L T = 0 if and only if the extrinsic vector field is a rigid motion.

• Kernel of L T → E V I . Information about the extrinsic vector field are solely contained by the matrices L T . Those matrices agree on edges so they can be reduced to the vector a ∈ R |E| as defined in Eq. (8.12). The mapping a → E V I is linear and almost always invertible as proven by Proposition 7.4.2.

• Rigidity Theorem. As shown previously the kernel of a → E V I is almost always reduced to the zero element. Going back to the matrices L T , the extrinsic vector field in the kernel should satisfy:

g -1 T L T = g -1 T E ⊤ (∇ E V ) T + g -1 T (∇ E V ) ⊤ T E -2Tr g -1 T E ⊤ (∇ E V ) Id = 0.
Taking the trace in both sides implies that div(V ) T = Tr(g -1 T ∂ t g t T | 0 ) should vanish so it is equivalent to have all matrices E ⊤ (∇ E V ) T + (∇ E V ) ⊤ T E equal zero. It follows that the extrinsic vector field satisfies e ij , V i -V j = 0 at all edges. It has been proved in [START_REF] Gluck | Almost all simply connected closed surfaces are rigid[END_REF] that almost all simply connected closed surfaces only admit rigid deformation as solution of this equation. Therefore, unlike continuous setting, we can guarantee that E V I is a coordinate-free representation of almost any extrinsic vector field V .

Representation in basis

Basis for the Function Space

Our infinitesimal shape difference is a linear operator acting on smooth functions defined on the surface. In practice, it is handy to use a basis for the function space, so any function can be represented as a linear combination of some basis functions φ i . In this basis, the operator E V I can be seen as a (possibly infinite) matrix. For a function f = i α i φ i we use the linearity property E V I (f ) = E V I ( i α i φ i ) = i α i E V I (φ i ). The choice of basis depends on the application. Since we want to represent smooth deformations of a surface, we take a subset of the smoothest functions given by the first eigenfunctions of the Laplace-Beltrami operator.

Vector field basis

As mentioned in Section 8.4.3 the infinitesimal shape difference operator E V I is linear with respect to V . Therefore, if the vector field is given in some basis V = i β i X i then the operator reads E V I = i β i E X i I . This means that when designing infinitesimal deformation V we can consider an objective as a function of the coefficients β, and moreover if the function is quadratic in β then recovering V can be done by solving a single linear system of equations. Note that this allows us to combine constraints both on the vector field itself (e.g., handle constraints) and on its operator representation, such as enforcing commutativity with other operators, while remaining easy to optimize for.

Choice of basis

In practice we use two alternatives for the basis of extrinsic vector fields. The simplest option is to take the eigenfunctions of the Laplace-Beltrami operator as the basis for each component of the deformation field. While simple, this basis might not preserve rotation invariance. Thus, alternatively we construct a basis via modal analysis of a deformation energy. In particular we consider an energy of the form V → M ∇V 2 dµ. This corresponds to the Dirichlet energy on a particular discretization of the Bochner Laplacian of extrinsic vector fields. To obtain the basis we take the eigenvectors of the Hessian of the energy function, which correspond to smooth deformation fields.

Experiments

In this section we apply our constructions to various tasks in deformation design and analysis. As our framework relies on manipulating and inverting moderately-sized matrices, all of the applications are very efficient and even when combining multiple objectives our method remains near interactive.

Deformation transfer

Pose transfer Given frames from an animation sequence and a functional map we can use our method to transfer the deformation to an arbitrary mesh. At each step we take a given deformation field U and find the transferred deformation V by minimizing Source (3k vertices) Target (5k vertices) Target (5k vertices)

Style transfer

In a similar experiment we use our approach to transfer style across the poses of different face shapes, shown in Figure 8.7, as we all shapes in the FAUST dataset [START_REF] Bogo | FAUST: Dataset and evaluation for 3D mesh registration[END_REF], shown in We then use our framework to transfer U to another shape in a different pose and with different mesh structure. In Figure 8.8 our method consistently preserves the global structure, although some high frequency details of the deformation are lost because of the basis reduction. The deformation transfer for faces in Figure 8.7 is successful because the deformation is very smooth and therefore well represented in the basis.

Symmetry transfer One interesting feature of the functional representation is that it is "shape aware." For example in Figure 8.9 we transfer the shrinking of the right leg to the left leg by looking for the operator which commutes with the operator representation of the symmetry map. Since both legs are in different positions this transfer is not easy by a simple point-to-point transfer of the vector field or even by transferring it using local coordinates. As shown by the vector field representation in Figure 8.9 bottom row, the transferred vector field adapts to the geometry. field. In Figure 8.10 we require that at a point p the deformation field matches a given vector u: V (p) = u in addition to other global constraints. We find the most isometric deformation by minimizing V → E V I 2

Deformations design

F . Given a self-map S, we design a symmetric vector field imposing a constraint of the form E V I C S -C S E V I = 0. We observe a similar deformation on each tentacles. We can also impose an anti-symmetry constraint with E V I C S + C S E V I = 0. Alternatively, regularization technique for deformation can be tested by imposing the commutativity with the Laplace-Beltrami operator, interestingly the deformation tends to spread to the entire shape. Figure 8.11 presents an example of joint deformation design. Namely, we impose a set of directional constraints on two different shapes and want a deformation on each shapes that is aware in some way of the deformation of the other shape. To do so, we find U and V minimizing (U, V ) → E U I M C ϕ -C ϕ E V I N 2 and respecting the local constraint on their respective shape. As a result, the constraints on one shape are transferred to the other. Moreover the area that could lead to contradictory deformation remains still.

Functional map inference

The infinitesimal shape difference is not only useful for design and analysis of deformations, it can also be used as a regularization in a shape matching problem. Below we show how this representation can be used to add extrinsic information to the computation of functional maps [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF]. In this framework it is challenging to use extrinsic information in another way than just using extrinsic based descriptors and while remaining rotationally invariant.

In this experiment we propose to solve the following problem: given two shapes and a Source Target In collection Figure 8.8 -The deformation field defined by the blue shapes (first column) is transferred to the same shape in different poses (top white shapes). While the style is consistent across the poses (red shapes) some details of the deformation are lost due to the basis representation. The style transfer are compared to the corresponding shape in the collection (bottom white shapes). sparse set of correspondences recover a dense map. The shapes come from the Faust dataset [START_REF] Bogo | FAUST: Dataset and evaluation for 3D mesh registration[END_REF] and we are given five corresponding landmarks at the hands, feet and head. The baseline method following the logic of the original paper is to represent the landmark points as delta functions δ M and δ N and look for the most isometric functional map C : L 2 (M ) → L 2 (N ), by enforcing commutativity with the Laplace-Beltrami operator. Thus, the straightforward approach would be to solve the optimization problem:

min C C∆ M -∆ N C 2 F s.t. Cδ M = δ N .
In order to add extrinsic information to this problem, we would like C to commute with the operator E n I representing the information about the second fundamental form. Note Initial shape Deformation Symmetrization Once the functional maps are obtained they are converted to a point-to-point map using the knn-algorithm as described in [START_REF] Ovsjanikov | Functional maps: a flexible representation of maps between shapes[END_REF]. The results are shown for two shape matching problems a nearly-isometric and a non-isometric. Figure 8.12 shows the percentage of correspondence within a given geodesic distance, and it can be seen clearly that the additional extrinsic information provides valuable information since the correspondence exhibits less error. Figure 8.13 provides a visualization of the point-topoint correspondences by transferring the coordinates functions encoded as RGB channels.

Interestingly our new constraint make the map smoother in both cases. We preform a similar experiment in Figures 8. 14 and 8.15, where instead of using the normal field as regularization we deform by hand the shapes to be matched with an equivalent deformation. This extrinsic vector field is used to regularized the shape matching problem. In both case, the maps obtained using commutativity constraint (right column) tend to be smoother than the naive approach (middle column).

Conclusion and Future Work

In this chapter we presented a method for representing extrinsic vector fields as linear operators acting on functions on the shapes, by considering the metric distortion induced by the deformation. For this we adapted the previously proposed shape difference operators by first introducing a unified shape difference, which fully encodes isometric distortion.

We then defined infinitesimal shape difference operators, which provide a compact functional representation of deformation fields. We showed how this representation can be used to analyze, transfer, and design deformations and to introduce extrinsic information into the computation of functional correspondences.

In the future, we are planning to use the newly introduced functional representation for shape animation, by exploiting its simplicity and ability to relate deformations on multiple shapes even without pointwise correspondences.

Chapter 9

Conclusion

In this thesis we have studied functional characterization of diffeomorphisms and embeded surfaces. The field of operator representation is relatively new in computer graphics and many questions remain open. Here we summarize some possible improvements and future works.

Shape to Deformation

This thesis provides a whole pipeline for shape matching. The first stage is the computation of an unknown functional map between shapes. The descriptor selection and stable subspace learning prove to be resilient to noise and robust to non-isometric matching.

The learning procedure is guided by given correspondences in a collection. However, this type of data might be difficult to obtain, whereas non-supervised optimization seems an achievable goals with functional maps. For example spectral descriptors are often related to geodesic distances. So tracking inconsistencies among descriptors might be a way of learning common intrinsic changes and moving parts in models.

The second stage, which consists of converting the functional map into a continuous point-to-point map attempts to improve currently available methods by adding continuity.

The method relies on the functional representation of vector field flow. However some diffeomorphisms are badly or not at all represented due to the use of a reduced function basis. Therefore, only a subset of vector fields generates representable maps. Currently we restrict vector fields to smooth bases but a more careful analysis of representability is needed to adapt the size of the vector field space to the size of the function space. Besides the hierarchical nature of eigenfunctions may be exploited into a "spectral multi-grid" algorithm.

A more long term challenge is to rethink this pipeline to solve non-isometric matching.

The conversion to a point-to-point map do not make any assumptions on the nature of the deformation as-long-as a valid functional map is given as input. The computation of the functional map however depends on function correspondences which are most of the time unreliable for large deformations. Moreover, the function basis, based on the Laplacian eigenfunctions, can be quite different on each mesh, meaning that large functional maps are needed. It seems that improvements can be obtained by considering adaptive basis function and localized descriptors.

Deformation to Shape

A functional characterization of deformable surfaces offers many possible future works in particular it defines an alternative shape space (e.g. [START_REF] Kilian | Geometric modeling in shape space[END_REF][START_REF] Heeren | Splines in the space of shells[END_REF]). Infinitesimal shape differences constitute a tangent space of a shape difference. This is a first step toward building a complete representation of a manifold of Riemannian metrics. We can then explore the shape space by following geodesics between surfaces or finding barycentric shapes. The most straightforward method to do so is by linear interpolation of the edge length and dihedral angle [START_REF] Winkler | Multi-scale geometry interpolation[END_REF] -discrete analogs of the metric and mean curvature. However, it is challenging to extent this view to other types of data whereas operator based methods rely on well-studied discretizations. To reach this goal, however, we are lacking a notion of metric between infinitesimal shape difference operators. The naive Frobenius norm does not have meaning outside the discrete setting and may lead to ill-conditioned problems. An analysis of the space of shape differences would help defining a suitable metric and geodesic distances as different metric and curvature may lead to different notion of interpolation of meshes. Defining a functional curvature representation stands as an unsolved problem. The solution provided by Chapter 7 is not entirely satisfactory since it requires to transform a surface to a tetrahedral mesh increasing the size of the problem. Moreover, global constraints linking the shape differences of the bottom layer to the top layer so it represents an embedded surface, are lacking and may not be easy to use in practice. Chapter 8 also suggests an answer but without theoretical proof of complete representation of triangle meshes. Nevertheless, this issue is fundamental in order to perform deformation based on curvature analysis (parametrization is an extreme example) or curvature interpolation. Two properties are desirable: the completeness of information, in the sense that we are able to recover the embedding from the shape differences and the curvature representation, and the insensitivity to noise and tessellation. Le deuxième développement concerne le problème inverse : étant donné une déformation représentée comme un opérateur, comment déformer une surface en conséquence ? Dans une première approche, nous analysons un encodage de la structure intrinsèque et extrinsèque d'une forme en tant qu'opérateur fonctionnel. Dans ce cadre, l'objet déformé peut être obtenu, à rotations et translations près, en résolvant une série de problèmes d'optimisation convexe. Deuxièmement, nous considérons une version linéarisée de la méthode précédente qui nous permet d'appréhender les champs de déformation comme agissant sur la métrique induite. En conséquence la résolution de problèmes difficiles, tel que le transfert de déformation, sont effectués à l'aide de simple systèmes linéaires d'équations.

Title : Functional representation of deformable surfaces for geometry processing Keywords : Geometry processing, Shape matching, Geometry algorithm Abstract :

Creating and understanding deformations of surfaces are recurring themes in geometry processing. As smooth surfaces can be represented in many ways from point clouds to triangle meshes, one of the challenges is being able to compare or deform consistently discrete shapes independently of their representation. A possible answer is choosing a flexible representation of deformable surfaces that can easily be transported from one structure to another. Toward this goal, the functional map framework proposes to represent maps between surfaces and, to further extents, deformation of surfaces as operators acting on functions. This approach has been recently introduced in geometry processing but has been extensively used in other fields such as differential geometry, operator theory and dynamical systems, to name just a few. The major advantage of such point of view is to deflect challenging problems, such as shape matching and deformation transfer, toward functional analysis whose discretization has been well studied in various cases. This thesis investigates further analysis and novel applications in this framework. Two aspects of the functional representation fra-mework are discussed. First, given two surfaces, we analyze the underlying deformation. One way to do so is by finding correspondences that minimize the global distortion. To complete the analysis we identify the least and most reliable parts of the mapping by a learning procedure. Once spotted, the flaws in the map can be repaired in a smooth way using a consistent representation of tangent vector fields.

The second development concerns the reverse problem : given a deformation represented as an operator how to deform a surface accordingly ? In a first approach, we analyse a coordinate-free encoding of the intrinsic and extrinsic structure of a surface as functional operator. In this framework a deformed shape can be recovered up to rigid motion by solving a set of convex optimization problems. Second, we consider a linearized version of the previous method enabling us to understand deformation fields as acting on the underlying metric. This allows us to solve challenging problems such as deformation transfer are solved using simple linear systems of equations.

Figure 3 . 1 -

 31 Figure 3.1 -Two charts and a transition map on a manifold M .

  The geodesic distance d : M × M → [0, +∞) is the minimum over all curves on M :d(p, q) := min c:[0,1]→M L(c)subject to c(0) = p, c(1) = q.

Figure 3 .

 3 Figure 3.2 shows an example of geodesic distance field on a surface.

Figure 3 . 2 -

 32 Figure 3.2 -Isolines of the geodesic distance, computed with the spectral method [34].

Figure 3 . 3 -

 33 Figure 3.3 -Three categories of mapping according to the intrinsic property they preserve.

Figure 3 . 4 -

 34 Figure 3.4 -Example of Laplace-Beltrami eigenfunctions computed on a surface. The first eigenfunction, not represented here, is constant.

Figure 3 . 5 -

 35 Figure 3.5 -Notations used for the Finite Element Method.

Figure 3 . 6 -

 36 Figure 3.6 -Local basis in a tetrahedra.

Theorem 4 . 1 . 1

 411 Let C : L 2 (M ) → L 2 (N ) be a non-zero linear operator. Then, there exists a map ϕ : N → M such that C = C ϕ if and only if:

Theorem 4 . 1 . 4

 414 Let C ϕ : L 2 (M ) → L 2 (N ) be a composition operator then,

1 ., 2 .

 12 Compute a basis (or an orthonormal family) of piecewise linear functions on M and N Use Equation (4.2) to obtain C M N .

  4 is still valid in a discrete sense if the meshes have same connectivity. Proposition 7.4.1 in Chapter 7 states that the mass matrix uniquely determines the triangle areas. So C ⊤ M N C M N is the identity matrix if and only if M and N have same triangle areas.

Figure 4 .

 4 [START_REF] Aigerman | Orbifold tutte embeddings[END_REF] shows examples of functional maps in three cases.

Figure 4 . 1 -

 41 Figure 4.1 -The first two rows represent three maps via color correspondence and the third their corresponding functional maps computed with 40 eigenfunctions. Three test cases are organized in three columns: a nearly isometric map is represented by an almost diagonal functional map, a partial matching becomes a slanted-diagonal matrix and a non-isometric (but almost conformal) map which do not exhibit a very sparse structure.

Figure 4 . 2 -

 42 Figure 4.2 -Left: The function f on M and its exact transfer on N . Middle: Representation of the function C M N f on the shape N when the functional map C M N is computed for varying k N and k M . Right: The relative error f • ϕ -C M N f N L 2 / f • ϕ N L 2 for various values of k N and k M . This experiment highlights the effect of basis sizes k N and k M on the transfer of a function by a functional map. When k M < k N the function is badly represented on M but the projection is well transferred to N . When k N < k M the function is well represented on M but the projection is badly transferred to N .

( 4 . 5 )

 45 allows us to rewrite Definition (4.3): M ∇f, ∇D C (g) dµ = M ∇f, B∇g dµ,

Figure 4 . 3 -Proposition 4 . 2 . 4

 43424 Figure 4.3 -Left: the point-to-point maps φ, ϕ linking manifolds M, N, P . Right: dual representation as functional maps.

Figure 4 . 4 -

 44 Figure 4.4 -Most distorted functions corresponding to the area and conformal-based shape differences for a nearly isometric deformation of the armadillo.

Figure 4 . 5 -

 45 Figure 4.5 -Left: The function f on M and its image by the conformal shape difference defined between the source and target meshes. Middle: Representation of the function D C (f ) when the functional map is computed for varying k N and k M . Right: The relative error between D C (f ) computed with a functional map for various values of k N and k M and the exact conformal shape difference. This experiment highlights the effect of basis sizes k N and k M on the conformal shape difference operator. When k M < k N the shape difference cannot represent high frequency distortion. When k N < k M the shape difference may represent parasite distortion due to the troncation of the basis.

  (a) Smoothed Gaussian curvature. (b) Logarithm of the absolute value of Gaussian Curvature.

Figure 5 . 1 -

 51 Figure 5.1 -Probe functions computed independently on two shapes. One carries meaningful information (a) and the other is misleading (b).

Figure 5 . 2 -

 52 Figure 5.2 -Pipeline of the proposed algorithm with the corresponding section. Top: learning from a given collection. Bottom: processing a new shape.

Figure 5 . 3 -

 53 Figure 5.3 -Left: Optimal weights for different strategies after training with 9 cats. Right: Effect of the different weights on the distribution of the singular value of the difference X n+1 (D) -C n+1 .

Figure 5 . 4 -

 54 Figure 5.4 -Visualization of the first component of the extracted basis. First row: the reference shape. Second row: transfer using different functional maps: (5.4a) naive map, (5.4b) unique weight, (5.4c) one weight per category of function, (5.4d) one weight per functions.

Figure 5 . 5 -

 55 Figure 5.5 -Left: Accuracy of the extracted function basis measured with Eq. (5.5). Right: Comparison between full map (plain lines) and reduced map (dash lines). Symmetric correspondences are considered correct.

Figure 5 . 6 -

 56 Figure 5.6 -Visualization of the first two components of the extracted basis on the reference shape (man) then mapped to the unseen shape (woman) using the functional map without optimization (left) and with the learned weight (right).

Figure 5 . 7 -

 57 Figure 5.7 -Visualization of the first two components of the extracted basis on the reference shape (gorilla) then mapped to the unseen shape (woman) using the functional map without optimization (left) and with the learned weight (right).

Figure 6 . 2 -Figure 6 . 3 -

 6263 Figure 6.2 -Example of a path trace starting point at p 0 .

Figure 6 . 4 -

 64 Figure 6.4 -Improvement of the BIM map using our method.

Figure 6 . 5 -

 65 Figure 6.5 -Comparison of the distortion induced by various for a pair of centaur (TOSCA).

Figure 6 . 6 -

 66 Figure 6.6 -Visualization of the point-to-point mapping through color correspondence.The texture of the first column (6.6a, 6.6d) are transferred to the second using the direct conversion of a functional map (6.6b, 6.6e) and to the third using our method (6.6c, 6.6f).

Figure 7 . 3 -

 73 Figure 7.3 -Two potential ways to encode extrinsic mesh structure.

MeshFigure 7 . 6 -

 76 Figure 7.6 -Canonical thickening (top right) of a triangle mesh (top left); types of tetrahedra in the thickening (bottom).

Figure 7 . 9 -

 79 Figure 7.9 -PCA on shape differences applied to recovering the sequence of animated frames for a simulated piece of cloth (top) and galloping horse (bottom).

Figure 7 . 10 -

 710 Figure 7.10 -Effects of truncation on computation of mesh structure. See §7.8.2 for discussion.

Source 2 .

 2 [START_REF] Gluck | Almost all simply connected closed surfaces are rigid[END_REF] 3.20 5.68 7.90 10.5 20.6 22.1 33.2 34.3 50.3 51.2 58.1 67.6 75.6 84.5 91.7 99.3 127. 136. Target 2.69 3.35 8.33 9.08 12.9 19.9 23.2 36.0 44.4 49.2 53.2 59.4 76.2 77.8 83.1 97.1 106. 131. 131. Recon. 2.78 3.52 8.10 9.46 12.7 20.8 23.7 36.9 44.9 51.0 54.9 61.0 79.6 81.0 84.4 102. 109. 133. 135. λ 2 λ 3 λ 4 λ 5 λ 6 λ 7 λ 8 λ 9 λ 10 λ 11 λ 12 λ 13 λ 14 λ 15 λ 16 λ 17 λ 18 λ 19 λ 20

Figure 7 .

 7 Figure 7.11 -Our machinery can be used to pull back Laplacian operators from a target mesh to a source. Here, we use a truncated functional map (100 Laplace-Beltrami basis functions on source, 200 on target) to compute revised edge lengths on the source mesh. Eigenfunctions of the Laplacian before and after edge length adjustment are shown with eigenvalues; boxed columns provide examples where the eigenfunction changes structure significantly.

  Figure7.12 -Our technique can be used to recover the pullback metric and therefore compute geodesic distances without direct access to the target mesh. We compare three geodesic pullbacks on the source mesh: the exact pullback using point-to-point correspondence, the geodesic computed by reconstruction of the metric from the shape differences, and the pullback of the geodesic function using a functional map. For each pullback we compute the L 2 distance d to the exact version. Our method achieves better reconstructions than the direct usage of a functional map.

Figure 7 .Figure 7 . 13 -

 7713 Figures 7.13, 7.14 and 7.15 illustrate experiments in which geometry is reconstructed after estimating local structure from shape differences. To highlight our method's effectiveness on extrinsic motion, we show behavior on human shapes and cloth simulation data. Figure7.13 applies our method to reconstructing models of humans from shape differences. From a coarse human base mesh (|X M | = |X N | = 502,; M = 100, k N = 200 in truncated experiments), we recover various poses. We compare reconstructions using only the intrinsic shape difference (right of each pair) to reconstructions using intrinsic and extrinsic differences together (left of each pair); we also compare using a truncated basis for shape differences (second column) to using a full basis (third column). As a baseline, we compare to[START_REF] Boscaini | Shape-fromoperator: Recovering shapes from intrinsic operators[END_REF], which uses only intrinsic geometry (rightmost column). Reconstruction from intrinsic information shows considerable artifacts due to the non-uniqueness of the solution of the embedding problem. Our provably complete intrinsic/extrinsic description is much more stable and close to the solution. The truncation of the basis, discussed in §7.8.2, tends to smooth out the sharp creases as they are represented as high frequency features. In Figure 7.14, we interpolate between frames of an animation sequence (|X M | = |X N | =

Figure 7 .

 7 15 illustrates a more challenging experiment (|X M | = 669, |X N | = 1089, k M = 60, k N = 180).

Figure 7 . 14 -Figure 7 . 15 -

 714715 Figure 7.14 -Mesh recovery and interpolation from a source mesh and the intrinsic/extrinsic shape differences. The target meshes come from a cloth simulation sequence.

Figure 7 . 16 -= 100 Figure 7 . 17 -

 716100717 Figure 7.16 -Performance measured on a 2015 iMac 3.3GHz.

d

  H = 0.128 d H = 0.073 d H = 0.065 d H = 0.036 d H = 0.023 Source Target k N = 100 k N = 200 k N = 300 k N = 400 Full No noise d H = 0.161 d H = 0.161 d H = 0.112 d H = 0.109 d H = 0.134 Source Target k N = 100 k N = 200 k N = 300 k N = 400 Full Noisy shape diffirence

Figure 7 .

 7 Figure7.18 -Impact of the basis truncation on from a source mesh and noisy shape differences. In this experiment, we fix the number of basis functions on the source shape to k M = 100 and reconstruct the embedding for various k N . Top row: With no additional noise, the quality of the embedding increases with k N . Bottom row: With added noise, larger k N -which normally yields better transfer of high frequency deformation-does not increase the quality of the reconstruction. The noisy shape differences do not correspond to an actual embedding.

  and 7.

Figure 8 . 1 -

 81 Figure 8.1 -Example of deformation design using our framework where all objectives are easily expressed as linear constraints. Left: A set of local constraints for the deformation fields on two different shapes. Middle: The deformations are computed separately on each shape to be as-isometric-as possible. Right: Joint deformation design with a soft (functional) map and no pointwise correspondences between the shapes.

Proposition 8 . 4 . 3

 843 Assume that D ϕ I : H 1 0 (M ) → H 1 0 (M ) represents the distortion of the metric between the surfaces M and P induced by the diffeomorphism ϕ : P → M and D φ I : H 1 0 (P ) → H 1 0 (P ) the distortion between the surfaces P and N linked through φ : N → P . The distortion D ϕ•φ

Figure 8 . 2 -

 82 Figure 8.2 -Middle row: Approximately conformal deformation of a bunny into a sphere (top row left).The PCA applied to shape differences confirms the presence of large area and isometric distortion in contrast to small conformal distortion. Bottom row: Area preserving deformation of a plane (top row right). The area shape difference is almost constant while conformal and isometric differences agree.

Figure 8 . 3 -

 83 Figure 8.3 -Two eigenfunctions associated with the largest eigenvalues of the operator infinitesimal shape difference E n I for the normal field n. The gradients of these functions represent the direction of principal curvature (Bottom row).

Theorem 8 . 5 . 1

 851 The modified cotangent weight matrix W M N (ℓ N ) uniquely determines the edge length ℓ N on N . In particular the following equivalence holds true: D I = Id if and only if ℓ M = ℓ N .

Figure 8 . 5 -

 85 Figure 8.5 -Deformation transfer from a given animation sequence (top row) to: another mesh with different connectivity (middle row), and in a different pose (bottom row) in an orientation-invariant way without using point or triangle correspondences.

Figure 8 . 8 .Figure 8 . 6 -

 8886 Figure 8.6 -Transfer of expressions in a faces collection. The expressions of the reference face (top row) are transferred to a different face (middle row) and compared to the corresponding faces of the collection (third row).

Figure 8 . 7 -

 87 Figure 8.7 -The deformation between the two shapes in the first column is transferred to a face in various expressions (red faces top row). The obtained deformation is consistent across the poses (red faces bottom row) even if the meshes have different connectivity and only an approximate functional map is used for the transfer.

Figure 8 . 9 -

 89 Figure 8.9 -An initial deformation (first two columns), corresponding to the shrinking of the right leg of a human model, is transferred to the left leg by imposing the commutativity between the infinitesimal shape difference and the symmetry map. Both legs are in different position so the transfer has to adapt to the geometry.

Figure 8 . 10 -Figure 8 . 11 - 2 F 2 F

 81081122 Figure 8.10 -We design deformations respecting the directional constraints shown on the far left and minimizing various criteria (from right to left): the infinitesimal shape difference leading to the most isometric vector field, the commutativity with the a selfmap, the anti-commutativity with the same self-map and the commutativity with the Laplace-Beltrami operator.

Figure 8 . 12 -

 812 Figure 8.12 -Percentage of correspondences within a geodesic ball for a nearly-isometric and a non-isometric shape matching problems. Note that adding the commutativity constraint with E n I (red lines) improves the correspondences.

Figure 8 . 13 -Figure 8 . 14 -

 813814 Figure 8.13 -Representation of the point-to-point map evaluated in Figure8.12. The RGB channel (left column) represents the xyz-coordinates, which are transferred using the recovered point-to-point map. In both cases, the maps obtained using commutativity with E n I (right column) tend to be smoother than the naive approach (middle column).

Figure 8 . 15 -

 815 Figure 8.15 -Representation of the point-to-point map evaluated in Figure 8.14.The RGB channel (left column) represents the xyz-coordinates, they are transferred using the point-to-point map. In both case, the maps obtained using commutativity constraint (right column) tend to be smoother than the naive approach (middle column).

Titre:

  Représentation fonctionnelle des surfaces déformables pour l'analyse et la synthèse géométrique Mots clés : Imagerie 3D, Association forme, Algorithmes géométriques Résumé :La création et la compréhension des déformations de surfaces sont des thèmes récurrents pour le traitement de géométrie 3D. Comme les surfaces lisses peuvent être représentées de multiples façons allant du nuage de points au maillage polygonal, un enjeu important est de pouvoir comparer ou déformer des formes discrètes indépendamment de leur représentation. Une réponse possible est de choisir une représentation flexible des surfaces déformables qui peut facilement être transportée d'une structure de données à une autre. Dans ce but, les "functional map" proposent de représenter des applications entre les surfaces et, par extension, des déformations comme des opérateurs agissant sur des fonctions. Cette approche a été introduite récemment pour le traitement de modèle 3D, mais a été largement utilisée dans d'autres domaines tels que la géométrie différentielle, la théorie des opérateurs et les systèmes dynamiques, pour n'en citer que quelques-uns. Le principal avantage de ce point de vue est de détourner les problèmes encore non-résolus, tels que la correspondance forme et le transfert de déformations, vers l'analyse fonctionnelle dont l'étude et la discrétisation sont souvent mieux connues. Cette thèse approfondit l'analyse et fournit de nouvelles applications à ce cadre d'étude. Deux questions principales sont discutées.Premièrement, étant donné deux surfaces, nous analysons les déformations sous-jacentes. Une façon de procéder est de trouver des correspondances qui minimisent la distorsion globale. Pour compléter l'analyse, nous identifions les parties les moins fiables du difféomorphisme grâce à une méthode d'apprentissage. Une fois repérés, les défauts peuvent être éliminés de manière différentiable à l'aide d'une représentation adéquate des champs de vecteurs tangents.

  , M. Ovsjanikov and A. Chambolle, Continuous matching via vector field flow, in Proceedings of the Eurographics Symposium on Geometry Processing, vol. 34, Eurographics Association, 2015, pp. 129-139. [32] • E. Corman, S. Solomon, M. Ben-Chen, L. Guibas and M. Ovsjanikov, Functional Characterization of Intrinsic and Extrinsic Geometry, currently under minor revision for Transaction On Graphics • E. Corman and M. Ovsjanikov, Functional Characterization of Deformation Fields, to be submitted (2016)

	Chapter 3

Table 6 .

 6

		Vertices Optimization Flow
	Horse	19248	369s	29.4s
	Dog	25290	300s	20.4s
	Centaur	15768	381s	39.0s
	SCAPE	12500	231.8s	30.8s

1 -Average CPU time of each step for different mesh size.

  and only if ϕ is an isometry,• D I = D C in H 1 0 (M )if and only if ϕ is an area-preserving.

	Proof Using Equation (4.5) we can rewrite Definition 8.4.1:

Remerciements

where α ij denotes the angle opposite to the edge (i, j). We introduce the energy term E : C ⊂ R |E| → R defined as

By construction the gradient yields ∇ ℓ 2 E(ℓ 2 ) = -w M (ℓ 2 ) + w, where w M is the vector containing the modified cotangent weight of Equation (8.8) computed with ℓ instead of ℓ N and w is a target vector of weights.

The triangle areas are also accessible through the determinant of the discrete metric, i.e. µ(T ) = det(g T )/2 in Prop 3.2.2. The function A → -log(det(A)) is strictly convex, i.e. [START_REF] Boyd | Convex Optimization[END_REF], and ℓ 2 → g T is an injective linear map so E is strictly convex on its domain of definition. Thus the optimization problem:

has, when it exists, a unique solution characterized by w M (ℓ 2 ) = w. As claimed, there is a one-to-one correspondence between C and w M (C). Suppose D I = Id it implies that the equality w M (ℓ 2 N ) = w M (ℓ 2 M ) holds true, thus the only solution is ℓ N = ℓ M . Theorem 8.5.1 suggests that the unified shape difference can be used for recovering intrinsic structure as developed in Chapter 7. However, it will most likely lead to highly non-convex optimization problems as we perform a projection of the weights in the least-squares sense.

Expression in a basis Similarly to the construction given in [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] we can also express the unified shape difference when the basis on the source shape M is given by the eigenfunctions of the Laplace-Beltrami operator Φ M . In that case, using a diagonal matrix Λ M of eigenvalues, the expression for D I becomes:

This expression has the advantage of avoiding the inverse of a large sparse matrix, and can be used to analyze deformation of a shape with fixed connectivity in multi-scale basis, which can make the computations resilient to local perturbations (see Option 3 in Section 5 of [START_REF] Rustamov | Map-based exploration of intrinsic shape differences and variability[END_REF] ).

Approximation with a functional map Note that both expressions above assume that the source and target meshes share the same connectivity. When the meshes have different connectivity this discretization requires a map between triangles making it challenging to use in practice. To overcome this problem we approximate this discrete formulation by transferring the weights on triangles to lumped weights on vertices. The approximation then reduces to the usual discrete quantities: