The matching problem between functional shapes via a BV penalty term: a Γ-convergence result

The matching problem often arises in image processing and involves finding a correspondence between similar objects. In particular, variational matching models optimize suitable energies that evaluate the dissimilarity between the current shape and the relative template. A penalty term often appears in the energy to constrain the regularity of the solution. To perform numerical computation, a discrete version of the energy is defined. Then, the question of consistency between the continuous and discrete solutions arises. This paper proves a Γ-convergence result for the discrete energy to the continuous one. In particular, we highlight some geometric properties that must be guaranteed in the discretization process to ensure the convergence of minimizers. We prove the result in the framework introduced in [11], which studies the matching problem between geometric structures carrying on a signal (fshapes). The matching energy is defined for L 2 signals and evaluates the fshapes difference in terms of varifolds norm. This paper maintains a dual attachment term, but we consider a BV penalty term in place of the original L 2 norm.


Introduction
Context and previous work.This paper discusses some theoretical aspects of the matching problem, which has several applications in image processing.The matching problem seeks a bijection between a current surface (or curve) and a target one by minimizing a dissimilarity function called matching energy.Such an error function is the sum of two terms: the penalty term defining the regularity of the optimal solution and the attachment term estimating the dissimilarity between current and target surfaces.
The solution to the matching problem defines the geometric transformation that links the two surfaces, enabling their comparison or transformation to a standard template.Moreover, this variational approach is implemented in the discrete framework using the stepwise descent algorithm, which provides the evolution related to the mentioned geometric transformation.The image processing community broadly studies the matching problem to define robust methods for shape analysis or image registration [13,24,7,11,25].
We point out, in particular, the increasing use of tools from geometric measure theory (currents, varifolds) for the definition of matching energies [11,28,12,19].F. Almegren has developed the theory of varifolds [4], afterward generalized by W.K. Allard [3], to generalize the notion of manifold in the framework of measure theory.Varfiold's theory considers a shape as a rectifiable measure and enables the definition of attachment terms in the sense of measures.This theoretical framework is suitable for obtaining compactness and existence results, but the computation of weak metrics remains a challenging problem.
Nowadays, new developments in non-invasive acquisition techniques, such as Magnetic Resonance Imaging (MRI) or Optical Coherence Tomography (OCT), provide geometric-functional data for several diseases (e.g., cortical thickness in the study of Alzheimer's disease or thickness of retina layers for the evolution of glaucoma).Then, new methods are needed to deal with data containing geometrical and functional information [20,21] to describe anatomical variability and produce statistical analysis for related diseases.
In this context, [11] introduces a new framework to study the matching problem of surfaces equipped with a signal (functional shapes or fshapes).The matching energy proposed in [11] considers the L 2 norm of the signal as a penalty term and a varifold-type attachment term.To overcome the difficulties linked to the discretization of weak varifold norm, the authors develop the model in the framework of Reproducing Kernel Hilbert Spaces (RKHS), which allows them to compute a varifold norm via the dual representation theorem in Hilbert spaces.
As usual, for minimization problems on surfaces, the question of error estimates for discrete solutions arises.The definition of a suitable discrete framework is central to ensuring that discrete solutions give a good approximation of continuous ones.In particular, the approximation quality should improve when the triangulation approximating the surface is finer (for example, in the context of finite elements, when the diameter of the triangles of the triangulation is sufficiently small).
Unfortunately, this is not true in general, and a very famous example called the Schwartz lantern (see [9] Section 3.9) shows that a cylinder can be approximated (in the Hausdorff topology) by a sequence of triangulations whose areas diverge.This example shows that the discretization process of a surface can highly affect the quality of the optimal solution.We finally note that, although we focus on the fshapes matching problem, the previous remarks hold for all variational problems on surfaces.
This paper aims to determine the conditions ensuring that the discrete solution represents a consistent (in terms of geometry) and close (in terms of energy) approximation of the continuous optimal solution.In particular, we prove that when the discretization process verifies some geometrical conditions, the discrete solution is a good approximation of the continuous one.The proof is given in the framework of the Γ-convergence theory, proving in particular that the discrete minimizers converge to the continuous one as the diameter of the triangles in the triangulation goes to zero.Error estimates are well established for Euclidean finite elements.However, to our knowledge, this kind of result has not been established for variational problems on surfaces and represents the main contribution of this paper.
Contributions.We call functional shape (fshape) any couple (X, f ) composed of a smooth surface X with boundary and a signal f defined on X.
A functional shape defines a varifold by considering the measure µ (X,f ) = H 2 X ⊗ δ T X (x) ⊗ δ f (x) .We study the matching problem between two fixed surfaces X, Y to optimize the signal f on X with respect to a target signal g defined on Y .The optimal signal is found by minimizing the following matching energy where ∥.∥ denotes a dual norm, and the minimum is taken on BV (X).This formulation differs slightly from the original model presented in [11].First, we consider a BV penalty term instead of the L 2 norm, proving the result for classic non-regular metrics.The regularity of f can strongly influence the behavior of the optimal solution.Standard L 2 signals allow one to work with smooth norms, but accumulation or oscillation phenomena may appear in the optimal solution.For this reason, gradient-dependent norms are increasingly used in image processing to guarantee more realistic optimal configurations.Furthermore, we prove that our results generalize to penalty terms The other main difference concerns the definition of the dual norm.This work considers a standard dual norm instead of the RKHS-based norm proposed in [11].This choice provides the result in a more general context and allows it to be read even without expertise in RKHS theory (Remarks 3.5 and 4.3 show how to adapt our results to the RKHS-based norm).Finally, the proof is given in the case of fixed geometry (X is fixed), and the generalization to the optimization problem concerning (X, f ) will be addressed in further work.
The first problem we address concerns the definition of admissible triangulation to obtain a geometrically consistent approximation of surfaces.To ensure a complete discretization of the surface in a neighborhood of its boundary, we consider a larger triangulation whose excess part of the projection of X has a small area (Definition 5.3).This allows us to overcome the bijection problems at the boundary of X.
Moreover, as explained above, the proximity in terms of Hausdorff distance does not guarantee a consistent surface approximation.The cited example of Schwartz lantern shows that when approximating a cylinder with triangulated surfaces, the area of the triangulations can even diverge depending on the geometric properties of the triangles.In particular, [23] shows that convergence of areas can be ensured by a condition on the angle between the normal vectors at corresponding points (via the projection) of the surface and the triangulation.This additional requirement is added to the definition of admissible triangulations to guarantee a convenient approximation of the related surface (see Lemma 5.5 and Hypothesis 2).
Once the set of admissible triangulations is defined, we recall the discrete version of the problem following [11].Section 6 describes in detail how to set the problem in the setting of finite elements via the projection of the triangulations on the surface.Then, the continuous problem can be approximated by a sequence of discrete problems defined on some triangulations whose triangle diameter goes to zero.The triangle diameter is also the parameter indexing the family of discrete problems, which ensures the convergence results for fine triangulations.
The main result proves that discrete solutions approximate continuous minimizers.The result is proven using the Γ-convergence theory, which is a notion of convergence of functionals introduced by E. De Giorgi [16,8], allowing to justify the transition from discrete to continuous problems.Γconvergence guarantees, in particular, the convergence of the discrete minimizers to the continuous one as the diameter of the triangles goes to zero (see Theorems 7.3 and 7.4).
We finally point out that, beyond the specificities of the matching problem for functional shapes, the quality problem for discrete approximations concerns many numerical problems defined on surfaces.This work shows how the approximation quality depends on the discretization process used for the numerical approach.Moreover, Γ-convergence is the proper framework to establish the consistency of numerical results.Then, our approach can generalize to other problems with promising theoretical and numerical results.

Structure of the paper.
In Section 2, we present the theory of BV functions on manifolds with boundary and adapt the classical results of approximation and compactness.In Section 3, we recall the framework of functional varifolds and the link with [11].In Section 4, we set the continuous problem and prove the related existence result.Section 5 is dedicated to the definition of admissible triangulations (Definitions 5.3) and the geometric conditions to guarantee the areas convergence (Hypothesis 2 and Proposition 5.6).Section 6 defines discrete operators on triangulations in the framework of finite elements.Finally, in Section 7, we consider the discrete problem and prove the Γ-convergence result (see Theorem 7.3 and 7.4).

BV functions on manifolds
In this section, we introduce the central notions and results about BV functions on manifolds.To develop a much larger setting beyond our framework, a similar definition is given in [6].
Although this work concerns approximation results on surfaces, we present here the BV theory on a slightly more general framework where X is a general compact d dimensional manifold (that is not supposed to be a finite 2D compact submanifold on R 3 ).
The definition of BV functions depends on introducing a divergence operator (or, equivalently, a volume form) and a local notion of length.A Riemannian structure provides these two things.
Let X denote an oriented smooth (at least C 1 ) compact d-dimensional Riemannian manifold, possibly with boundary denoted ∂X, and let vol X be the associated Riemannian volume form.The boundary ∂X is supposed to be a C 1 compact (d−1)-dimensional manifold and we have vol X (∂X) = 0. Finally, let us denote which is a C 1 manifold without boundary.When X is without boundary, the previous construction gives X 0 = X.We say that f ∈ L 1 (X) is a function of bounded variation on X if where χ 1 c (X 0 ) denotes the set of C 1 vector fields u : X → T X on X compactly supported in X 0 and div X is the divergence operator on X defined by where (e 1 , • • • , e d ) is an orthonormal frame on T X.Here ∥u∥ ∞ = sup x∈X g x (u(x), u(x)) 1 2 where g is the metric tensor associated with the Riemannian structure.
We recall the integration by parts formula : where h ∈ C 1 c (X 0 ) and u(f ) denotes the derivative of f along the vector fields u that is defined by [u(f )](x) = d x f (u(x)) for any x ∈ X 0 .We retrieve the usual submanifold setting when considering the metric induced on the submanifold by the ambient space R 3 .
The functional space BV (X) endowed with the norm is a Banach space.
Definition 2.1.The space BV (X) can also be equipped with the following notions of convergence, both weaker than the norm convergence: 1. Weak- * topology.Let {f h } h ⊂ BV (X) and f ∈ BV (X).We say that the sequence {f h } h weakly- * converges in BV (X) to f if where * ⇀ denotes the weak convergence in the space of measures on X.

Strict topology.
Let {f h } h ⊂ BV (X) and f ∈ BV (X).We say that the sequence is a distance in BV (X) inducing the strict convergence.
We recall the main compactness result for functions of bounded variation extended to the manifold framework.
Theorem 2.2 (Compactness).Let X be a smooth manifold and {f h } h be a sequence of BV (X) such that ∥f h ∥ BV (X) is uniformly bounded.Then {f h } h is relatively compact in BV (X) with respect to the weakly- * convergence.
Proof.The proof consists in adapting that one given for the Euclidean case (see [5], Theorem 3.23) to a manifold with a boundary.As X is compact, we can consider a finite atlas {(U i , φ i )} i=1,...,n where {U i } i=1,...,n is a finite open cover of X (note that {U i ∩ X 0 } i=1,...,n is a finite open cover of X 0 ) and φ i : is the dimension of the manifold).Without loss of generality, we can suppose that each V i is a Lispchitz domain, and we can also consider a partition of unity {η i } i=1,...,n on X subordinate to the cover {U i } i=1,...,n .Now, by Proposition 2.2 in [6], as )) and we can apply the classical theorem (see [5], Theorem 3.23) to {(f h η i ) • φ −1 i } h for any local chart.The extension argument used in the classical compactness theorem allows us to define a subsequence converging in L 1 (V i ).
Then, via at most n extractions of subsequences, we can define a subsequence (not relabeled) such that, for every i, {f h η i } h weakly*-converges to some We now establish an approximation result of BV functions by smooth functions.Interestingly, even if there exists several density results using smoothing through the heat kernel semi-group for geodesically complete Riemannian manifold (see [22], [10] and [17]), to the best of our knowledge, an approximation result in the case of a compact manifold with boundary does not seem to be available.We give below such a result : Theorem 2.3.Let X be an orientable Riemannian compact manifold with boundary ∂X and let (2.2) Proof.The proof shares similar ideas to the proof of the classical approximation result in the Euclidean case (see Theorem 3.9 in [5]), Furthermore, it is based on the parameterization of the function by the local charts.When considering the problem in local charts, introducing a non-constant volume term introduces several new elements.The proof is detailed in the Appendix.
In the following, we return to the case of a submanifold of R 3 and show how the previous result can be improved.To achieve this goal, we need to introduce an extension operator to be able to expand the manifold and its signal in a consistent way.
Let X denote a C p (p ≥ 2) compact oriented 2 dimensional submanifold of R 3 with non empty boundary denoted ∂X.We denote n X 0 the C p−1 vectors field of positively oriented normal along X 0 = X \ ∂X.Note that n X 0 can be continuously extended on X and denoted in that case n X .For any x ∈ ∂X, we can define the unit vector ν(x) pointing outward and orthogonal to n X (x) and ∂X.For sufficiently small r > 0, we consider the open subset where π X 0 (z) is the orthonormal projection of z on X 0 and |t X 0 (z)| is the distance of z to X 0 .Now, there exists η 0 > 0 small enough, such that the mapping is well defined and is a C p−1 diffeomorphism on an open neighborhood of ∂X in R 3 (see Figure 1).Moreover, Ext maps ∂X × {0} × {0} to the boundary ∂X of X and for 0 < η < η 0 We can now prove the following approximation result: Theorem 2.4.Let X be a C p (p ≥ 2) compact oriented 2D submanifold of R 3 with non empty boundary denoted ∂X.Let f ∈ BV (X) and let ϵ > 0. Then there exists η > 0 and f Moreover, there exists 0 < η ′ < η 0 such that Let δ > 0, we can easily check that for η ′ small enough since by considering the change of variable induced by Ext we notice that (x, s) → π X (x + sν(x)) has a determinant converging to 1 when (x, s) convergence to a point on ∂X × {0} with s < 0.Moreover, we verify that so that, choosing δ = 1 and η ′ small enough, we have and we get immediately that proves the result.
Remark 2.5 (L 2 and H 1 norms).Theorems 2.3 and 2.4 can be established, by similar arguments, for the L 2 and H 1 norms.We refer to [18] for an introduction to Sobolev spaces on manifolds.

Functional varifolds
This section presents the main definitions and results for functional varifolds as introduced in [11].We recall that this framework is a generalization of the varifold theory [29] to spaces with a scalar component; this allows the representation of geometric shapes equipped with a signal in the framework of measure theory.We refer to [29,11] for a more detailed presentation of the respective theories.The functional approach considers mathematical objects, called fshapes, containing geometrical and functional components.As detailed below, the role of fshapes in the functional varifolds framework is the same as that of manifolds in the varifold theory; they allow one to define measures supported on surfaces and are the bridge between geometry and measure theory.
In the following, we work with smooth surfaces and denote by X a generic 2-manifold verifying the following properties : Hypothesis 1. X denotes a 2-submanifold of R 3 (surface) with boundary ∂X and which is smooth (C 2 at least), oriented, connected and compact.The smoothness assumption implies that its interior, denoted by X 0 , is a two-dimensional smooth manifold, and ∂X is a one-dimensional manifold of the same regularity.Definition 3.1 (fshape).We define a fshape as a couple (X, f ) where X is a surface verifying Hypothesis 1 and f ∈ BV (X) is a signal defined on X.
Of course, the penalty term in the matching problem defines the regularity of the signal.In the following, we consider fshapes endowed with L 2 or H 1 signal when studying the problem with L 2 or H 1 penalty terms, respectively.
Similarly to the classical theory of varifolds, we can define a functional varfiold via the dual of C 1 0 (R 3 ×G(3, 2)×R), which denotes the closure, with respect to the C 1 norm, of the set of C 1 functions with compact support on R 3 ×G(3, 2)×R.We denote by G (3,2) the Grassmannian of the non-oriented 2-dimensional linear subspaces of R 3 .

Definition 3.2 (fvarifolds). A 2-dimensional functional varifold (fvarifold) is any operator
In the following, we consider the following norm for fvarifolds : We can associate a fvarifold to a given fshape (X, f ) by considering the following measure that acts as a linear functional in the following way : where T X (x) denotes the tangent space to X at x and H 2 is the 2-dimensional Hausdorff (or volume) measure.
Finally, we prove some lemmas that will be useful in the following, concerning some convergence properties for fvarifolds supported on fhsapes.In particular, the convergence to the null fvarifold depends on the geometric component : Proof.For every h and for every The following lemma links the convergence of signals and the weak- * convergence of fvarifolds : which proves the results by taking the supremum and using the L 1 convergence.
Remark 3.5 (The choice of the fvarifold norm).This work uses the dual norm as a distance between fvarifolds.However, in [11], another metric for fvarifolds is defined by using the framework of Reproducing Kernel Hilbert Spaces (RKHS) [31].We refer to [11] for a general definition, and we recall that their framework allows them to define the following dual product for fvarifolds supported on fshapes: with where σ e , σ t , σ f are three positive constants and n T represents the unit normal vector to T .In particular, W ′ denotes the dual space of W , the RKHS associated with the kernel k e ⊗k t ⊗k f (see Propositions 2 and 4 in [11]).As W is continuously embedded into , its dual norm defines a metric for fvarifolds.
The choice of this metric has many advantages.It allows them to prove an existence result for the matching problem with L 2 signals (see Section 4).Moreover, in numerical applications, using the Gaussian kernel allows the localization of the matching at a given scale.
In this work, we decided to consider the matching problem defined for BV or H 1 signals and use the standard dual norm for fvarifolds.This framework is strong enough to prove the existence of optimal solutions and a Γ-convergence result.Throughout the paper, we will demonstrate how to adapt the results to the original framework introduced in [11].

Existence of optimal solutions for the matching problem
In this section, we define the matching energy between two fshapes and prove an existence result for the optimal solution.
Let be X a surface and (Y, g) a target fshape.We consider the following energy defined for a generic fshape (X, f ) and we aim to solve the minimization problem We recall that we optimize only with respect to the signal, which implies that the initial and optimal configurations have the same geometric support.However, the geometry is taken into account in the attachment term.Theorem 4.1.Let (Y, g) a given fshape and X a 2-manifold verifying Hypothesis 1.Then, there exists at least one solution to the problem Proof.Let {f h } h be a minimizing sequence belonging to BV (X).We can suppose that ∥f h ∥ BV (X) is uniformly bounded and, by Theorem 2.2, {f h } h converges (up to a subsequence) to some f ∈ BV (X) with respect to weak- * topology.The result follows by remarking that the fvarifold norm is continuous with respect to the L 1 topology (Lemma 3.4) and that the BV norm is lower semicontinuous with respect to the L 1 topology.

Remark 4.2 (
The H 1 model).We get the same result if we consider signals belonging to H 1 (X) instead of BV (X).It follows from the fact that the unity ball of H 1 (X) is compact with respect to the weak topology.Moreover, H 1 (X) is compactly embedded in L 2 (X), which implies (up to a subsequence) the L 1 convergence of the minimizing sequence.
Remark 4.3 (L 2 model in [11]).The L 2 model is defined (see [11], part I) via the following matching energy where γ f , γ W are two positive constants, and the fvarivold norm is induced by (3.3).We point out that, for every minimizing sequence {f h } h , a bound on E(f h ) guarantees only the L 2 weak compactness for the signals, which is not enough to get the semicontinuity of the fvarifold term (a result similar to Lemma 3.4 holds for the W ′ -norm with respect to the a.e.convergence of the signals).This justifies the choice of the RKHS-based metric described in Remark 3.5, whose properties play a central role in the existence result : Theorem 4.4 (Proposition 7 in [11]).Let X, Y be two finite volume bounded 2-rectifiable subsets of If the ratio γ f /γ W is large enough, then there exists at least one solution to the minimization problem inf and every minimizer belongs to L ∞ (X).Moreover, if X is a C p surface and W → C m 0 (R 3 ×G(3, 2)×R) with m ≥ max{p, 2}, then every minimizer belongs to C p−1 (X).
Finally, there exists a constant C > 0 (independent of X and Y ) such that, every minimizer The proof is based on the relaxation of the energy to the class M X defined as follows: Note that µ (X,f ) ∈ M X for every fshape (X, f ).
Then, the energy E can be relaxed to the following functional and it holds : The relaxed energy provides compactness for minimizing sequences in the space of measures.It can be shown that the minimizing measure ν * of Ẽ is associated with a fshape, so that ν * = µ (X,f * ) for some f * ∈ L 2 (X).The proof relies on the Implicit Function Theorem, which needs the hypothesis on γ f /γ W .We refer to Proposition 7 and Lemma 2 in [11] for more details.
However, the L 2 model has two main issues.Firstly, the existence result depends on the weights used to define the energy.Moreover, the L 2 penalty does not prevent some oscillating configurations for the optimal signal.For this reason, we introduced a penalty on the signal derivative that justifies the BV and H 1 models studied in this work.

Surfaces and triangulations
This section is dedicated to defining suitable triangulations of surfaces to ensure a good approximation of continuous functions by their discretization.
As explained in the Introduction, the quality of approximation and the error estimates highly depend on the geometric properties of triangulations.The famous example of the Schwartz lantern (see [9] Section 3.9) exhibits a sequence of polyhedral surfaces converging to the cylinder in the Hausdorff metric and whose areas diverge.This example points out that triangulations must verify some specific hypothesis to represent a geometrically consistent surface approximation.
Section 5.1 recalls some general facts about triangulations and defines the set of admissible triangulations in Definition 5.3.In Section 5.2, we give a sufficient condition (Hypothesis 2), ensuring the convergence of areas for admissible triangulations.Finally, in Section 5.3, we define the sampling method for signals and the error estimate for Sobolev norms.

Triangulations of a surface
As X verifies Hypothesis 1, following [1,27,26], we give the following definition of triangulation: Definition 5.1 (Triangulations).A triangulation T is a two-dimensional manifold (with boundary) consisting of a finite set ∆ T of affine triangles such that: 1. any point p ∈ T lies in at least one triangle T ∈ ∆ T ; 2. each point p ∈ T has a neighbourhood that intersects only finitely many triangles of ∆ T ; 3. the intersection of any two non-identical triangles T, T ′ ∈ ∆ T is either empty or consists of a common vertex or edge.
In the following, we denote by T a generic triangle of ∆ T and by ∂T the boundary of the manifold T .
The vertices of T are the vertices of the triangles belonging to ∆ T .The diameter of T is defined as follows : We also assume that every triangulation is regular, which means that for some C > 0, where h T is the diameter of T and ρ T is the diameter of the sphere inscribed in T .
We define the distance function to X as the following function : For every We also recall that the Hausdorff distance between two surfaces X, Y ⊂ R 3 is defined as Definition 5.2 (Tubular neighborhoods).Let X be a surface satisfying Hypothesis 1.We denote by U r (X) the subset of R 3 of the form such that every point x ∈ U r admits a unique projection π X (x) ∈ X.We refer to [15] for the proof of the existence of such a tubular neighborhood U r for some r > 0.
To guarantee the injectivity of the projection, we introduce the following tubular neighborhood where we denote by n X (x) the unit normal vector to X at x.Then, we have x ∈ N r (X) if and only if d X (x) < r and Then, we can split every triangulation into two parts to isolate the set (T in ) of points that can be bijectively projected onto the surface: Generally, a triangulation is not in bijection with a surface with a smooth boundary through the normal projection.The bijectivity can fail close to ∂X because of the curvature of ∂X as depicted in Figure 2a.Locally, the normal projection of ∂X on a hyperplane must not be a line, and it is impossible to project ∂X on the edge of a triangle.
Then, we introduce a new class of admissible triangulations to ensure the projection bijectivity between triangulation and surface except on a small part close to the boundary: Definition 5.3 (h-admissible triangulations for a surface).Let h > 0. We say that a triangulation T is h-admissible for the surface X if the following properties hold (i) T lies in N h (X +η ) for some η > 0, where X +η is an extension of X defined in (2.4); (ii) T in = T ∩ N h (X) and X are in one-to-one correspondence through π X ; Figure 2a: the triangulation T 1 is not in one-to-one correspondence with X through the projection map (for instance, the part of the smooth star in red exceeds the triangulation).Figure 2b: the subset T in 2 of T 2 is in one-to-one correspondance with X .
This definition introduces a discretization framework based on triangulations larger than the corresponding surface but converging to it with respect to the Hausdroff distance (as h → 0).Condition (i) means that T is adapted to discretize an extended surface and, because of (ii), X can be completely projected onto the triangulation.However, the part T out , corresponding to the triangles along the boundary, is assumed to be small via condition (iii).Finally, condition (iv) allows us to identify the family parameter with the triangulation diameter in order to describe finer triangulations as h → 0. Fig. 2 shows an example of admissible and inadmissible triangulations.
We point out that the vertices of triangulations are not, as usual, a set of sampled points on the surface.This approach corresponds better to the data acquisition routines used in image processing.For instance, biomedical images (OCT, functional MRI) are often modified (segmentation, deblurring, denoising) to improve their quality for numerical experiments.Then, the data do not correspond to a sampled version of the imaged objects, but instead, they represent an approximation, often noisy, of the natural structures.

Convergence of areas
As pointed out in the previous section, if T is an h-admissible triangulation of X for some h > 0 then d H (X, T ) = O(h).However, this does not guarantee a similar error estimate for the respective areas.
In this section, we introduce some geometric conditions to ensure geometrically consistent triangulations and avoid pathological cases as the already cited Schwartz lantern (see [9] Section 3.9).

Definition 5.4 (Angle between normals). Let h > 0. Assume that T is an h-admissible triangulation for X. For every x ∈ T in we define the angle α x as follows
• if x belongs to the interior of some triangle, then α x is the angle belonging to [0, π/2] between the two normals n X (π X (x)) and n T (x); • if x belongs to an edge of a triangle, then α x is the biggest angle belonging to [0, π/2] between n X (π X (x)) and the normals of the triangles containing x.
In the following, we set Lemma 5.5.Let h > 0 and T be an h-admissible triangulation for X.We have as α 2 max , d H (X, T ) → 0, where α max is introduced in Definition 5.4.
Proof.For every x ∈ X we consider on the tangent space T x X the basis B(x) = {e 1 (x), e 2 (x)} given by the two principal directions.We denote by κ 1 (x) and κ 2 (x) the principal curvatures of X at x. Similarly, for every x ∈ X we can consider the basis B(x) = {e 1 (x), e 2 (x), n X (x)} for R 3 .Consider the differential Dπ X : R 3 → T π X (x) X of π X .Note that, for every x ∈ N h (X), we have Dπ X (x)(v) = 0 for every variation v in the direction n X (π X (x)) orthogonal to T π X (x) X.Then, we should consider the tangential variations to calculate the projection's Jacobian.
In [30] it is proved that for any x ∈ U h (X) and v parallel to T π X (x) X we have where the matrix of Dπ X (x) written with respect to the basis B(π X (x)) and B(π X (x)) is be a subset of a triangle of T .This implies that the Jacobian of the projection on X restricted to A is given by We have cos α x = 1 + O(α 2 max ).Moreover, as the principal curvatures are uniformly bounded and as α 2 max , d H (X, T ) → 0. We can conclude by the usual change of variables in the area formula.
Because of Lemma 5.5, let us introduce the following set of assumptions: Hypothesis 2. Let {T h } h be a sequence of triangulations indexed by a parameter h → 0 such that (i) for any h > 0, the triangulation T h is h-admissible for X; (ii) the sequence α h max = O(h) as h → 0, where α h max is the angle defined in Definition 5.4 for T h .
We point out that, because of Definition 5.3, we get d H (X, T h ) = O(h), which implies that, under the assumption of Hypothesis 2, the area error computed in Lemma 5.5 is also O(h).Then, Hypothesis 2 implies the geometric consistency of the approximation and helps avoid pathological cases such as the Schwartz lantern.
The following result holds, generalizing Corollary 5 of [23] to surfaces with smooth boundaries : Proposition 5.6 (Convergence of the area).Let X be a surface satisfying Hypothesis 1 and {T h } h a sequence of triangulation satisfying Hypothesis 2. Then we have lim h→0 Proof.The proof follows from Lemma 5.5 and Definition 5.3 (iii).

From the triangulation to the surface
This section defines how to carry a signal from the triangulation to the surface.

Definition 5.7 (Projection).
For every function f defined on an h-admissible triangulation T for X (h > 0) we define the projection of f onto X by We point out that the function f ℓ carries on X the signal defined on T in .
Proposition 5.8.Let h > 0 and T be and an admissible h-triangulation for X.Then, for every f ∈ W 1,∞ (T , R), we have Proof.The first equality is proved by performing the change of variables y = π X (x) and using Lemma 5.5.The second relationship is proved by the same arguments and applying the chain rule.
We point out that, for triangulations verifying Hypothesis 2, the previous proposition generalizes to surfaces with boundary the Sobolev error estimates for Euclidean finite elements.

Discretization
This section aims to write down the discretized matching energy to compare the discrete problem defined on triangulations to the continuous one defined on the original surfaces.
For clarity, we detail the definition of all discrete norms and functionals in the framework of finite elements.Although it may seem redundant to a reader accustomed to this kind of formalization, we prefer to define the several operators to prevent misinterpretation of notations.

Notations
Let T be a triangulation in the sense of Definition 5.1.We denote by N v , N e , and N t the number of vertices, edges, and triangles of T , respectively.The family {v i } i=1,...,Nv denotes the vertices of the triangulation.For every k = 1, . . ., N t , we denote by {v k i } i=1,2,3 the vertices of the triangle the center of the edge linking v k i to v k j , and v k 0 = 1 3 3 i=1 v k i the center of mass of the triangle (see Figure 3).Analogously we denote by {f k i } 0≤i≤3 the values of the function f at location {v k i } 0≤i≤3 of T k .

P 0 and P 1 triangular finite elements
Let us start with the following definition, Definition 6.1.For a given triangulation T , we denote by P 0 (T ) (resp.P 1 (T )), the set of functions that are constant on the interior of each triangle and null on their edges (resp.the set of the continuous functions that are affine on each triangle).
The elements of P 0 (T ) (resp.P 1 (T )) are completely described by their values {f k 0 } k at the center of mass {v k 0 } k (resp.their values On the other hand, if f ∈ P 0 (T ), as f is null on the edges of each triangle, we can not compute f k 0 by the values at the vertices, and we can only set We denote by p 0 the L 2 projection of P 1 (T ) on P 0 (T ) where for f ∈ P 1 (T ), the function p 0 (f ) is the unique element of P 0 (T ) such that In other words, the operator p 0 replaces the affine approximation of a signal on each triangle with a constant approximation using the value at the center of mass.
A basis for P 1 (T ) (barycentric basis) is given by the family {φ} i=1,••• ,Nv with φ i ∈ P 1 (T ) and φ i (v j ) = δ ij (Kronecker's delta), for every i, j = 1, • • • , N v .Then, every f ∈ P 1 (T ) can be written as Remark that there exists a bijection between P 1 (T ) and R Nv , defined by the following operator

Discrete functional norms
Depending on how the continuous signal is discretized, various methods may be used to compute the norm of the discrete signal.
L p norm of P 0 finite elements.Let p ≥ 1 and f : T → R a function in L p (T , R).The pth power of the discrete L p norm of f is simply defined as 3) The formula (6.3) is exact for signals that are (almost everywhere) constant on each triangle so that we have L p 0 [f, T ] = ∥f ∥ p L p (T ) for any f ∈ P 0 (T ).
L p norm of P 1 finite elements.As in this work, we are concerned by L p norms with p = 1, 2, we need to define the discrete operators via an exact formula on piecewise quadratic polynomials.We use the following approximation by the evaluation at the midpoints (see [2] p. 178-179): ). (6.4) where for a piecewise linear signal f ∈ P 1 (T ) we have The formula (6.4) is a polynomial of degree 2 for any k = 1, . . ., N t .The formula (6.4) is also exact if p = 1 and f has constant sign on each triangle (i.e.we have . Suppose the signal f ∈ P 1 (T ) has a changing sign on a given triangle T k .In that case, the function to integrate is not a polynomial on the entire triangle, and the formula is not exact anymore.Then, the triangle is decomposed into several sub-triangles with constant sign signal and the operator in (6.4) is computed as the sum of the respective operators computed (exactly) on each subtriangle.This is equivalent to performing a local triangulation refinement to ensure the signal has a constant sign on each triangle.
We suppose now that f ∈ P 1 (T ) and we define the discrete operators corresponding to the H 1 and BV norms.For every T k ∈ T and for every f ∈ P 1 (T ), the gradient of f on T k can be computed by where The previous relationship can be stated by writing f in the barycentric coordinates system and recalling that the gradients of the basis elements are perpendicular to the triangle edges.
In this framework, the gradient ∇ T f is constant on each triangle and, by convention, is null on the edges.We can now define the discrete operators of the gradient norms that are exact for every f ∈ P 1 .
Total variation.The total variation of f ∈ P 1 on T is given by H 1 norm.For f ∈ P 1 (T ), the square of the L 2 norm of the gradient is given by and by formula (6.5) we have ∥∇ T f ∥ ∈ P 0 (T ).

Discrete fvarifold norm
The definition of discrete fvarifolds can be posed by considering a fvarifold associated to (T , f ).We point out (T , f ) does not verify the conditions (smoothness) defining fshapes; however, the fshape definition given in Section 3 remains consistent with the properties of T .Then, with an abuse of language, we call (T , f ) discrete fshape, and the related fvarifold supported on it is defined as follows : Moreover, we approximate such a fvarifold by a discrete operator µ (T ,f ) to simplify the computation.Such an approximation is set in the same way in both P 0 and P 1 framework : where The discrete fvarifold norm is defined by We point out that by the formula (6.1), we have We end this section with a technical lemma that will be useful in the following section.It proves an error estimate, with respect to the fvarifold norm, between a discrete signal f (belonging to P 0 or P 1 ) and its projection f l on the surface (see Definition 5.7).Lemma 6.2.Let X be a surface verifying Hypothesis 1 and T a h-admissible triangulation for X verifying Hypothesis 2. Then we have and Proof.Let f ∈ P 0 (T ).By the change of variables x = π X (y), because of Hypothesis 2 and formula (5.5), for every f ∈ P 0 (T ) and Moreover, because of Definition 5.3(iii) and Lemma 3.3, µ (T out ,f ) converges towards the null fvarifold as h → 0. Then where N t denotes the number of triangles contained in T and V k the tangent space to T k .So, by (3.1), we obtain and, by the triangle inequality, we get ∥µ (X,f ℓ ) − µ (T ,f ) ∥ ≤ O(h) which proves (6.11).
and, by the same arguments, we get which proves (6.12)

Discretization of continuous signals
Finally, we present the discretization process (on admissible triangulations) of a signal defined on a surface.Let X be a surface verifying Hypothesis 1 and T be a h-admissible triangulation for X (for some h > 0), and f ∈ BV (X) (resp.H 1 (X), L 2 (X)).
According to Definition 5.3(i), by using the map Ext defined in (2.3), we can extend the manifold X to a larger suitable manifold X +η such that T ⊂ N h (X +η ).Moreover, because of Theorem 2.4, the signal f can be extended to a signal f defined on X +η which is W 1,∞ (X +η ) and with a small norm on X +η \ X −η .
Let T be a triangle contained in T , we set the following discretization process: • P 0 elements: as every P 0 element is null on each edge of T , we define the piecewise (a.e.) constant f h corresponding to f on T as where v 1 , v 2 , v 3 denote the three vertices of T .
• P 1 elements: the piecewise linear f h corresponding to f on T is defined as where v 1 , v 2 , v 3 denote the three vertices of T and P 1 is defined by (6.2).
This process defines a discrete signal f h on the triangulation whose different norms can be computed by the operators introduced above.We note that these operators represent an approximation of their continuous counterparts, which motivates the study of the relationship between the discrete and continuous minimizers discussed in the next section.

Discrete problem and Γ-convergence results
This section aims to show the Γ-convergence of the discretized problems to the continuous one when the triangulation is fine enough.The main consequence is that the discrete optimal solution computed on triangulations represents a good approximation of the optimum computed on the surface (Theorem 7.4).The proof is provided under the geometric assumptions established in Hypothesis 2. This implies, in particular, that the area error (see Lemma 5.5) is O(h), ensuring the geometric consistency of the approximation process.
Let X, Y be two surfaces verifying Hypothesis 1, and g ∈ BV (Y ).We denote by {T h } h , {Y h } h two sequences of h-admissible triangulations of X and Y , respectively, verifying Hypothesis 2. We recall that the control on the angle between continuous and discrete corresponding normal vectors prevents pathological cases like the Schwartz lantern.We finally denote by {g h } h the discretization of g on the sequence of triangulations {Y h } h obtained by the discretization process described in Section 6.4.
We point out that h denotes the family parameter and the triangulation diameter (the maximum diameter of the triangles contained in T h ).Then, the Γ-convergence results hold for the typical case of increasingly fine triangulations.Finally, we note that the regularity condition (5.2) is fundamental in our proof to ensure global error estimate in the finite elements framework.
Then, for every h, the discrete energy is defined by where the operators L 1 1 , V , Var are define in (6.4), (6.6), and (6.10), respectively.For every h > 0, the optimal signal on T h is defined by the following problem inf The main goal of this section is to prove that the minimizers of E h converge to the minimizers of the continuous problem First of all, we prove the existence result for the discrete problem : Proposition 7.1.For every h > 0, there exists at least one solution to problem (7.2).
Proof.For every Then, every minimizing sequence is bounded in BV so that it weak- * converges in BV (up to a subsequence) to some function belonging to P 1 (T h ).The result follows from Lemma 3.4 and the lower semicontinuity of the total variation with respect to the weak- * convergence.
As the discrete and continuous energies are not defined in the same space, we introduce a suitable topology to compare the two spaces and generalize the classical definition of Γ-convergence (see [8]): Definition 7.2 (S-topology and Γ-convergence).Let X be a surface satisfying Hypothesis 1, f ∈ BV (X), and {T h } h be a sequence of admissible triangulations for X verifying Hypothesis 2. In this definition, {f h } h denotes a sequence of functions such that f h ∈ P 1 (T h ) for every h > 0.
• We say that {f h } h converges to f with respect to the S-topology where for any h > 0, f ℓ h is the projection of f h onto X (see Definition 5.7).• We say that (E h ) h Γ-converges to E if the following conditions hold: (i) Lower bound: for every f ∈ BV (X) and for every sequence (ii) Upper bound: for every f ∈ BV (X) there exists a sequence {f h } h such that f h S ⇀ f and We can now prove the main result of this work : Moreover, by Lemma 6.2, we have Now, as the BV norm is lower semicontinuous with respect to the L 1 convergence, we get Upper bound.Because of Theorem 2.4, we can assume that f ∈ C 1 (X) and get the general result by a diagonal argument.
The discretization process defines, for every h > 0, an extension X +η of X to extend f to f ∈ W 1,∞ (X +η ) (Theorem 2.4) and define the discrete function f h (see Section 6.4).Of course, the extension X +η can depend on h, but in the following, we prove some estimates for a fixed h > 0 and write X +η .In particular, as f is defined via a reflection symmetry with respect to the boundary of X (see Theorem 2.4), we have which implies that f h ∈ W 1,∞ (T h ) for every h, and sup The usual estimates for interpolation error, as h → 0, give This can be deduced from Theorem 3.1.6in [14] (applied with k = s = 0, q = m = 1, p = ∞) by considering the sum on all the triangles and by using the fact that the area of every triangle is bounded by π(h/2) 2 .We note, in particular, that these estimates hold for regular triangulations verifying (5.2).Now, according to Theorem 2.4, f coincides with f on X −η and its BV norm can be arbitrarily small on X \ X −η .Then, we get Moreover, as T h virifies Definition 5.3(ii) , X −η can be projected on a subset of T h , and by using (7.9) and (7.10), we get This proves in particular that f ℓ h → f strongly in L 1 (X), which means that f h S ⇀ f .Now, by Proposition 5.8 and Hypothesis 2, we have as h → 0. Similarly, the convergence of the total variation term follows from As in the case of the lower bound, the convergence of the fvarifold term follows by Lemmas 6.2 and 3.4, which finally proves The Γ-convergence result implies, in particular, the convergence of the minimizers : Proof.We consider the sequence {f h ∈ P 1 (T h )} h of minimizers of E h and, without loss of generality, we can also suppose that sup Similarly to (7.6), we have so that {f ℓ h } h is uniformly bounded in BV (X).Then, there exists a subsequence (not relabeled) such that {f ℓ h } h weak- * converges to f ∞ ∈ BV (X) and, by Lemmas 6.2 and 3.4, we get Then, by lower semicontinuity the BV norm with respect to the L 1 topology, we get The other inequality follows from the upper bound condition of Γ-convergence applied to a minimizer of E. This proves, in particular, that f ∞ minimizes E.

Remark 7.5 (H 1 model). The discrete problem for H 1 signals is defined by the discrete energy
Then, Theorems 7.3 and 7.4 still hold in this case, and their proofs can be adapted by considering as S-topology the L 2 convergence of the projection (i.e., f h S ⇀ f if f ℓ h → f strongly in L 2 (X)) and using the compact embedding in L 2 (X) and the compactness of the unit ball with respect to the weak topology of H 1 (X).Remark 7.6 (L 2 model).In this case, the discrete energy is defined on the set of P 0 finite elements by the following function to minimize on P 0 (T h ): We must suppose that γ f /γ W is large enough to apply Theorem 4.4.We note that in the case of W ′ -fvarifold norm (see Remarks 3.5 that gives the needed inequality.The other one follows from the upper bound condition of Γ-convergence.
In particular, we get min and, from Theorem 4.4, there exists

Appendix: Approximation theorem for BV functions on manifolds
Theorem.Let X be an orientable Riemannian compact manifold with boundary ∂X and let X 0 = X \ ∂X.For any f ∈ BV (X) there exists a sequence {f h } h ⊂ C 1 (X 0 ) such that Proof.Let {U i } i∈ 1,n be a finite atlas on X 0 and for any i ∈ 1, n let φ i : U i → V i be a local chart such that U i is compact and φ i is the restriction to U i of a C 1 diffeomorphisms from U i → V i (such an atlas exists since X is compact).Let {η j } j≥0 be a partition of unity such that supp(η j ) is compact for any j ≥ 0 and there exists a partition {J i } i∈ 1,n of N for which supp(η j ) ⊂ U i for any j ∈ J i which is locally finite on any U i (i.e., for any x ∈ U i , there exists an open set U i (x) ⊂ U i such that supp(η j ) ∩ U i (x) = ∅ for any j ∈ J i except on a finite number of j's).For any j ∈ J i we consider ϵ j such that d(φ i (supp(η j )), V c i ) > ϵ j and for ϵ = {ϵ j } j≥0 we consider the linear operator L ϵ : BV (X) → C 1 (X 0 ) defined by where ψ i : V i → U i is the inverse mapping of φ i .We recall the classical notation of differential geometry for pullbacks where for any function ℓ ∈ C c (U i ), ψ * i ℓ = ℓ • ψ i and for any v ∈ χ 1 c (U i ), ψ * i v = (dψ i ) −1 v • ψ i .We recall that χ 1 c (X 0 ) denotes the set of C 1 vector fields u : X → T X on X compactly supported in X 0 .Eventually, on every V i , we introduce α i dx the pullback of vol X ↾ U i by ψ i on V i such that for any ℓ ∈ C c (U i ), we have U i ℓvol X = V i (ψ * i ℓ)α i dx.Let δ > 0. We can assume that for any 1 ≤ i ≤ n and any j ∈ J i , we have ϵ j small enough so that Since {η j } j≥0 is a partition of unity, we have f = n i=1 j∈J i f η j and This first inequality is enough to prove an approximation result in a L 1 sense.We turn now to the control of the total variation part.Let u ∈ χ 1 c (X 0 ).We have the following decomposition using the integration by part formula (2.1) for equality (a) and the classical integration by part on R d for equality (b) However, for any δ 0, denoting | • | x the norm at x ∈ X induced by the metric, for x ∈ supp(η j ) and ϵ j small enough we have uniformly in u such that ∥u∥ ∞ ≤ 1 and j ∈ J i .Hence Moreover, for ϵ j small enough, for x ∈ supp(η j ) we can assume d ψ(φ i (x)−y) φ i u(ψ i (φ i (x) − y))ρ ϵ j (y)dy x ≤ δ so that Thus we have Let us consider now the B ij 's.We have For the B 2 ij terms, let us notice that sup V i |α i ψ * i u| < ∞ uniformly in u (since ∥u∥ ∞ ≤ 1) and (ψ * i f ∇(ψ * i η j )) * ρ ϵ j → (ψ * i f ∇(ψ * i η j )) in L 1 (R d , dx) so that for ϵ j sufficiently small, we can assume that |B 2 ij | ≤ δ2 −(j+1) .Summing along the indices, we get and, with (7.15), (7.17) and (7.18), we get eventually that, for sufficiently small values of the ϵ j 's, we have uniformly in u ∈ χ 2 c (X 0 ) satisfying ∥u∥ ∞ ≤ 1. Taking the supremum over such u, we get Since δ is arbitrary, we have shown that there exists a sequence (ϵ k ) k≥0 such that L ϵ k f ∈ C ∞ (X 0 ) and lim sup .21)

Figure 2 :
Figure 2: The surface X is a bent smooth star (solid grey), and two triangulations (black lines) are illustrated.

Figure 3 :
Figure 3: Labels of various points in the triangle T k .
Theorem 7.3.The sequence {E h } h (7.1) Γ-converges to E (7.3) with respect to the S-topology of Definition 7.2.Proof.Lower bound.Let {f h } h be a sequence of functions such that f h ∈ P 1 (T h ) for every h > 0 and f h Theorem 7.4 (Convergence of minimizers).Let {T h } h be a sequence of admissible triangulations for X verifying Hypothesis 2. Let {fh } h be a sequence of minimizers of E h (i.e., E h (f h ) = min f ∈P 1 (T h ) E h (f )).Then, {f ℓ h } h weakly- *converges in BV (X) (up to a subsequence) to a minimizer of E and lim and 4.3), Lemma 6.2 still holds and Lemma 3.4 can be proved under the hypothesis of a.e convergence of signals.For L 2 signals, it is sufficient to define the S-topology, f h where M X is defined in Definition 4.5.Because of (4.6), we get the lower semicontinuity for the lower bound of Γ-convergence.The upper bound follows similarly to the BV case, knowing that the discretization process guarantees the strong convergence in L 2 .Concerning the convergence of minimizers, Theorem 4.4 (in particular (4.4)) implies that the sequence of discrete minimizers {f h } h is uniformly bounded in L ∞ .Then, the sequence of measures {µ (X,f ℓ h ) } h ⊂ M X is tight and, because of the Prokhorov's theorem, it weak- * converges (up to a subsequence) to some µ ∞ ∈ M X .Now, the lower semicontinuity of Ẽ with respect to the weak- * convergence of measures and the fact that Ẽ is minimized by a fshape associated to a L 2 function (seeRemark 4.3), we get min