1 Introduction

The importance of various types of symmetry is evident while solving problems such as shape segmentation, mesh repairing, shape matching, retrieving the normal forms of 3D models [12, 30, 60], inverse procedural modeling [31], shape recognition [14], shape understanding [34], shape completion [47, 49], and shape editing [10, 15, 20, 33, 54, 55]. The problem of detecting intrinsic symmetry is shown to be an NP-hard problem since it amounts to finding an intrinsically symmetric point for each point [37]. However, correspondences between the intrinsically symmetric points completely characterize the intrinsic symmetry of a shape since the intrinsic symmetry is a non-rigid transformation which can not be represented a matrix as opposed to the extrinsic symmetry which can be represented by a matrix [16]. We exploit the functional map approach that finds correspondences between functions instead of points [36]. Then, the point-to-point correspondences can be recovered in \(O(n\log (n))\), where n is the number of vertices. We extend this framework for the detection of intrinsic symmetry. The functional map framework has already been used for detecting intrinsic symmetry in previous works [24, 52]. The main task in these approaches was to determine the functional correspondence matrix which transforms a function to its intrinsic image. Various constraints have been enforced on this matrix. It is not known yet, how many constraints are sufficient. Further, they solved a non-linear optimization problem to estimate the functional correspondence matrix which makes the method slow. For the intrinsic symmetry detection problem, we show that the functional correspondences matrix is a diagonal matrix and a diagonal entry is \(+1\) (\(-1\)), if the corresponding eigenfunction is an even (odd) function. This closed-form solution makes our method faster.

We determine if a particular eigenfunction is even or odd based on our following result. An eigenfunction is an even (odd) function if its restriction to the shortest length geodesic between two intrinsically symmetric points is an even (odd) function. Therefore, we need to find a few accurate pairs of intrinsically symmetric points which we find using our following results. If we directly pair points based on the similarity between their heat kernel signatures [48], we may get false pairs. For example, a pair of points on the tip of the index finger and the tip of the ring finger of the same hand of a human model. The reason is that, if two neighboring points are subjected to the same strength heat sources, then their heat diffusion processes will also be similar because of the very small sizes of the fingers with respect to the body size. However, we observer that the sign of a low-frequency eigenfunction on the neighboring points is the same. Hence, we put high penalty for pairing two points if signs of first few low-frequency eigenfunctions are the same for both the points. The models of the benchmark datasets are obtained by applying an imperfect isometry, so the theory only holds approximately. Furthermore, some of the triangles may not be Delaunay triangles. Therefore, we may not get accurate correspondences using the original eigenfunctions. Hence, we transform the original eigenfunctions to make them near perfect even or odd functions. Following are our main contributions.

  1. 1.

    We propose a novel approach to find a few accurate pairs of intrinsically symmetric points based on the following property of eigenfunctions: the signs of low-frequency eigenfunction on neighboring points are the same.

  2. 2.

    We propose a novel and efficient approach for finding the functional correspondence matrix. We prove that the functional matrix for the intrinsic symmetry detection problem is a diagonal matrix and a diagonal entry is \(+1\;(-1)\) if the corresponding eigenfunction is an even (odd) function.

  3. 3.

    We propose a novel approach to determine the sign of an eigenfunction by showing that, if a manifold contains intrinsic symmetry and an eigenfunction is an even (odd) function, then its restriction to the shortest length geodesic between any two intrinsically symmetric points is an even (odd) function.

  4. 4.

    We transform the eigenfunctions to make them more invariant to self-isometry.

2 Related Works

Reflective symmetry detection methods are categorized by the types of the input data and the reflective symmetry present in the input data. The input data can be a digital image, point cloud, triangle mesh, etc. The main types of reflective symmetries are extrinsic and intrinsic. The well known methods for detecting the extrinsic symmetry are [21, 26, 28, 29, 40, 45,46,47, 51] and for detecting the intrinsic symmetry are [5, 6, 13, 17, 18, 22, 27, 37, 38, 41, 41, 42, 44, 45, 53, 56, 58, 60]. Furthermore, the intrinsic symmetry can be characterized as global and partial intrinsic symmetry. Our method finds global and partial intrinsic symmetry in triangle meshes. We discuss only the relevant works and suggest the readers to follow the excellent state-of-the-art report in [32] and the survey in [25].

Ovsjanikov et al. detected intrinsic symmetry using the global point signature (GPS) [37]. The main claim was that the GPS embedding transforms the problem from intrinsic to extrinsic symmetry detection. They showed that the GPS embedding was robust to the topological noises. They found pairs of symmetric points by comparing the GPS of one point to the signed-GPS of the other point. The time complexity of determining the sign (even or odd) of an eigenfunction is \(O(n\log (n))\) since they compared GPSs of all the possible pairs. Furthermore, the time complexity of the overall method is \(O(k^3 n\log (n))\) excluding the computation of eigenfunctions, where k is the number of eigenfunctions used. We propose a more efficient method which takes (\(O(kn\log (n))\)) for detecting intrinsic symmetry. Furthermore, we observe that the approach by [37] is sensitive to the sign flip and eigenfunction ordering. Whereas, our method is independent of sign flip and ordering since we determine the sign of each eigenfunction independently from the others. Mitra et al. used a voting based approach to detect intrinsic symmetry and then applied transformation in the voting space to deform the input model to have perfect extrinsic symmetry [30]. Xu et al. used a generalized voting scheme to find the partial intrinsic symmetry curve without explicitly finding the intrinsically symmetric point for each point [57]. Xu et al. efficiently found pairs of intrinsically symmetric points using a voting based approach [56]. They factored out symmetry based on the scale of symmetry. However, they needed to tune a parameter depending on how much the input shape is distorted [56]. The methods by Zheng et al. [13, 60] also used voting approach. These voting based methods do not utilize spatial coherency. Therefore, they may produce pairs of intrinsically symmetric points which may not be spatially continuous. Furthermore, they may have high complexity due to a large number of possible pairs for the voting.

In [41], the authors proposed a non-convex optimization framework to accurately detect full and partial symmetries of 3D models. However, the initialization severely affects the performance, and the complexity also is very high. Lipman et al. efficiently found the pairs of intrinsically symmetric points in point clouds and triangle meshes using the novel symmetry factored embedding technique [23]. However, their main bottleneck is the time complexity which is \(O(n^{2.5}\log (n))\). Kim et al. used anti-Möbius transformation to accurately find intrinsically symmetric pairs of points [16]. They first find a sparse set of pairs of intrinsically symmetric points. Then, they transform intrinsic symmetry into extrinsic symmetry using the Möbius transform. However, they required \(O(n^2)\) space for mid-edge flattening and \(O(|\mathcal {S}|^4)\) for finding the anti-Möbius transformation, where \(\mathcal {S}\) is the set of symmetry invariant points. Furthermore, false pairs in the first step may severely affect the overall performance.

3 Intrinsic Symmetry Detection

3.1 Background

Let \(\mathcal {M}\) be a compact and connected 2-manifold representing the input shape. Let \(L^2(\mathcal {M})=\{f:\mathcal {M}\rightarrow \mathbb {R}|\left\langle f,f\right\rangle _\mathcal {M}=\int _{\mathcal {M}}f^2(x)dx<\infty \}\) be the space of square integrable functions defined on \(\mathcal {M}\). The Laplace-Beltrami operator on a shape \(\mathcal {M}\) is defined as \(\Delta _\mathcal {M}f=-\text {div}_\mathcal {M}(\nabla _\mathcal {M}f)\) and admits an eigenvalue decomposition \( \Delta _\mathcal {M}\phi _i(x)=\lambda _i\phi _i(x),\forall x\in \mathcal {M}\). Here, \(0=\lambda _1\le \lambda _2\le \ldots \) are the eigenvalues and \(\phi _1,\phi _2,\ldots \) are the corresponding eigenfunctions. The eigenfunctions \(\phi _1,\phi _2,\ldots \) form a basis for the space \(L^2(\mathcal {M})\). Therefore, any function \(f\in L^2(\mathcal {M})\) can be represented as \( f(x)=\sum _{i=1}^{\infty }\left\langle f,\phi _i\right\rangle _\mathcal {M}\phi _i(x), \forall x\in \mathcal {M}. \) The functional map framework was first proposed in [36] for establishing point-to-point dense correspondence between two isometric shapes. The main idea was to establish correspondences between the functions, defined on the shapes, rather than the points. This idea reduced the time complexity to \(O(n\log n)\). Let \(\mathcal {M}\) and \(\mathcal {N}\) be two shapes. Let \(T_{\text {f}}:L^2(\mathcal {N})\rightarrow L^2(\mathcal {M})\) be a linear mapping between functions defined on these shapes. That is, if \(g:\mathcal {N}\rightarrow \mathbb {R}\) and \(f:\mathcal {M}\rightarrow \mathbb {R}\) are two corresponding functions then \(T_{\text {f}}(g)=f\). The mapping \(T_{\text {f}}\) is represented by a matrix \(\mathbf {C}\in \mathbb {R}^{k\times k}\) such that \(\mathbf {b}=\mathbf {Ca}\), where \(\mathbf {a}=\begin{bmatrix} a_1&a_2&\ldots&a_k\end{bmatrix}^\top \) and \(\mathbf {b}=\begin{bmatrix} b_1&b_2&\ldots&b_k\end{bmatrix}^\top \) are the representations of the functions g and f in the truncated bases \(\{\phi _i^\mathcal {N}\}_{i=1}^k\) and \(\{\phi _i^\mathcal {M}\}_{i=1}^k\), respectively. Therefore, the main goal is to find the matrix \(\mathbf {C}\) which completely characterizes the dense correspondence between the two shapes.

3.2 Functional Maps for Intrinsic Symmetry Detection

We extend the functional map framework for detecting the intrinsic symmetry which can be thought of as a shape correspondence problem where we have to find the correspondences between the points of the same shape rather than the points on the two different shapes. The functional map framework is applicable for two isometric shapes also. Therefore, we can use it to detect the intrinsic symmetry since a symmetric shape is a self-isometric shape [37]. The intrinsic symmetry \(T_{\text {p}}:\mathcal {M}\rightarrow \mathcal {M}\) of \(\mathcal {M}\) is defined as follows. If the points \(x\in \mathcal {M}\) and \(y\in \mathcal {M}\) are intrinsically symmetric, then \(T_{\text {p}}(x)=y\) and \(T_{\text {p}}(y)=x\). We first find the mapping between the functions and then use it to find the correspondences between the intrinsically symmetric points. Let us consider the space \(L^2(\mathcal {M})\) and let \(T:L^2(\mathcal {M})\rightarrow L^2(\mathcal {M})\) be a functional map which maps the functions defined on the same shape. Then, this functional map T completely characterizes the intrinsic symmetry \(T_{\text {p}}\) if \(T(g)=f\) and \(T(f)=g\), where \(f,g\in L^2(\mathcal {M})\) are intrinsically symmetric functions, i.e. \(f\circ T_\text {p}(x)=g(x),\;g\circ T_\text {p}(x)=f(x), \forall x\in \mathcal {M}\). Therefore, our goal is to find the matrix \(\mathbf {C}\) which characterizes the functional mapping T for the intrinsic symmetry detection problem. For the problem of finding correspondences between two shapes, various constraints have been imposed on the matrix \(\mathbf {C}\). Then, the matrix \(\mathbf {C}\) was the optimal solution of an optimization problem. However, we show that a closed form solution exists for the matrix \(\mathbf {C}\) for the problem of detecting the intrinsic symmetry, which we state as follows.

Theorem 1

Let \(T:L^2(\mathcal {M})\rightarrow L^2(\mathcal {M})\) be a mapping between the functions defined on a shape \(\mathcal {M}\) and T characterizes the intrinsic symmetry \(T_\text {p}\) of \(\mathcal {M}\), i.e. \(T(g)=f, T(f)=g,\;\forall f,g\in L^2(\mathcal {M})\) such that \(f\circ T_\text {p}(x)=g(x),\;g\circ T_\text {p}(x)=f(x), \forall x\in \mathcal {M}\). Then, the matrix \(\mathbf {C}\) representing T is a diagonal matrix. \(\mathbf {C}_{i,i}=+1\), if \(\langle T(\phi _i),\phi _i\rangle _\mathcal {M}=+1\), and \(\mathbf {C}_{i,i}=-1\), if \(\langle T(\phi _i),\phi _i\rangle _\mathcal {M}=-1\).

Proof

The functions fg belong to the space \(L^2(\mathcal {M})\). We can represent f and g as \(f(x)=\sum _{i=1}^{\infty }b_i\phi _i(x)\) and \(g(x)=\sum _{i=1}^{\infty }a_i\phi _i(x)\). Since T is a linear mapping, we have \( T(f)=T\left( \sum _{i=1}^{\infty }b_i\phi _i(x)\right) = \sum _{i=1}^{\infty }b_i T\left( \phi _i(x)\right) \). Since \(T\left( \phi _i(x)\right) \) is also a function in the space \(L^2(\mathcal {M} )\), it can be represented in the basis \(\{\phi _i\}_{i=1}^\infty \) as \(T\left( \phi _i(x)\right) =\sum _{j=1}^{\infty }c_{ij}\phi _j(x)\), where \(c_{ij}=\langle T(\phi _i),\phi _j\rangle _\mathcal {M}\). Therefore, we have that \(T(f) = \sum _{j=1}^{\infty } \sum _{i=1}^{\infty }c_{ij}b_i \phi _j(x)\). Since \(T(f)=g\), it follows that \(\sum _{j=1}^{\infty } \sum _{i=1}^{\infty }c_{ij}b_i \phi _j(x)=\sum _{j=1}^{\infty }a_j\phi _j(x)\). Therefore, \(a_j=\sum _{i=1}^{\infty }c_{ij}b_i\). Equivalently, we can write it as \(\mathbf {a}=\mathbf {C}\mathbf {b}\), where \(\mathbf {C}_{i,j}=c_{ij}=\langle T(\phi _i),\phi _j\rangle _\mathcal {M}\). According to [37], the eigenfunctions (corresponding to the non-repeating eigenvalues) are self-isometry invariant with sign ambiguity i.e. \(\phi _i\circ T_{\text {p}}(x)=\pm \phi _i(x)\), \(\forall x\in \mathcal {M}\). Furthermore, the functional map T completely characterizes the intrinsic symmetry. Therefore, \(T(\phi _i)=+\phi _i\), if \(\phi _i\circ T_{\text {p}}(x)=\phi _i(x)\), and \(T(\phi _i)=-\phi _i\), if \(\phi _i\circ T_{\text {p}}(x)=-\phi _i(x)\), \(\forall x\in \mathcal {M}\). Since the eigenfunctions \(\phi _1, \phi _2,\ldots \) form an orthogonal basis for the space \(L^2(\mathcal {M})\), we have that \(\langle \pm \phi _i,\phi _j\rangle _\mathcal {M}=0\) if \(i\ne j\). Therefore, \(\mathbf {C}_{i,j}=\langle T(\phi _i),\phi _j\rangle _\mathcal {M}=\langle \pm \phi _i,\phi _j\rangle _\mathcal {M}\). Therefore, \(\mathbf {C}_{i,j}=0\), if \(i\ne j\). Hence, the matrix \(\mathbf {C}\) is a diagonal matrix. Furthermore, \(\mathbf {C}_{i,i}=+1\), if \(\langle T(\phi _i),\phi _i\rangle _\mathcal {M}=+1\), and \(\mathbf {C}_{i,i}=-1\), if \(\langle T(\phi _i),\phi _i\rangle _\mathcal {M}=-1\).   \(\square \)

Therefore, the problem of determining whether \(\mathbf {C}_{i,i}=+1\) or \(\mathbf {C}_{i,i}=-1\) is equivalent to determining whether \(\phi _i\circ T_\text {p}(x)=+\phi _i(x)\) or \(\phi _i\circ T_\text {p}(x)=-\phi _i(x)\), \(\forall x\in \mathcal {M}\). It is observed that if \(\phi _i\circ T_\text {p}=+\phi _i\) then \(\phi _i\) is a symmetric or an even function and if \(\phi _i\circ T_\text {p}=-\phi _i\), then \(\phi _i\) is an anti-symmetric or an odd function in the intrinsic sense. We can not apply the definition of the vector space here, since the domain of the eigenfunctions is not a vector space. A function \(f:\mathbb {R}^2\rightarrow \mathbb {R}\) is an even function, if \(f(-\mathbf {x})=f(\mathbf {x}),\forall \mathbf {x}\in \mathbb {R}^2\). This definition is not valid for the functions defined on the manifolds, since if \(x\in \mathcal {M}\), then it may not always be true that \(-x\in \mathcal {M}\). However, we generalize the following property of vector spaces to the manifolds to determine the sign of eigenfunctions. Let \(f:\mathbb {R}^2\rightarrow \mathbb {R}\) be a function symmetric on \(\mathbb {R}^2\) and \(\ell =\{\mathbf {x}:\mathbf {x}=t\mathbf {x}_1+(1-t)\mathbf {x}_2, t\in [0,1]\}\) be the line segment joining the mirror symmetric points \(\mathbf {x}_1\) and \(\mathbf {x}_2\). Then, it is trivial to show that the restriction \(f_{\ell }:\ell \rightarrow \mathbb {R}\) of the function f on the set \(\ell \) is also symmetric. Here, we also observe that the set \(\mathbb {R}^2\) is symmetric about any of its coordinate axes and the set \(\ell \) is also symmetric. We formally generalize these results on manifolds as follows.

Theorem 2

Let \(\mathcal {M} \) be a compact and connected 2-manifold. Let there exist a self-isometry \(T_\text {p}:\mathcal {M}\rightarrow \mathcal {M}\) on \(\mathcal {M}\). Let \(x,y\in \mathcal {M}\) be two points which are intrinsically symmetric, i.e., \(T_\text {p}(x)=y\) and \(T_\text {p}(y)=x\). Let \(\gamma (t):[0,1]\rightarrow \mathcal {M}\) be the shortest length geodesic curve between the points x and y such that \(\gamma (0)=x\) and \(\gamma (1)=y\). Then, \(T_\text {p}(\gamma (t))=\gamma (1-t)\) and \(T_\text {p}(\gamma (1-t))=\gamma (t), \forall t\in [0,1]\).

Proof

Let \(\beta (t)=T_\text {p}(\gamma (t))\). We have to show that \(\beta (t)=\gamma (1-t)\). Since \(T_\text {p}\) is an isometry, according to (Proposition 16.3, [11], Chap. 3, p. 91 [35]) \(T_\text {p}\) maps a shortest length geodesic on \(\mathcal {M}\) to a shortest length geodesic on \(\mathcal {M}\). Therefore, \(\beta (t)\) is also a shortest length geodesic. Now, we have \(\beta (0)=T_\text {p}(\gamma (0))=T_\text {p}(x)=y\) and \(\beta (1)=T_\text {p}(\gamma (1))=T_\text {p}(y)=x\). Therefore, \(\beta (t)\) is the shortest length geodesic between the points x and y such that \(\beta (0)=y\) and \(\beta (1)=x\). Since there can only be a single shortest length geodesic curve between two points (except continuous symmetry, like sphere), both the geodesics \(\gamma (t)\) and \(\beta (t)\) trace the same path. However, their start and end points are flipped. Therefore, \(\beta (t)=\gamma (1-t)\Rightarrow T_\text {p}(\gamma (t))=\gamma (1-t)\). Since the self-isometry is an involution, i.e. \(T_\text {p}\circ T_\text {p}(x)=x,\forall x\in \mathcal {M}\), we have \(\gamma (t)=T_\text {p}\circ T_\text {p}(\gamma (t))=T_\text {p}(T_\text {p}(\gamma (t)))=T_\text {p}(\gamma (1-t))\Rightarrow T_\text {p}(\gamma (1-t))=\gamma (t)\).   \(\square \)

The intuitive is that if a shape is intrinsically symmetric, then the shortest length geodesic curve between any two intrinsically symmetric points is also intrinsically symmetric. This result helps us to determine the sign of the eigenfunctions of the Laplace-Beltrami operator. First, we show that the result \(\phi _i\circ T_\text {p}(x)=\pm \phi _i(x)\;\forall x\in \mathcal {M}\) holds true if we restrict the eigenfunctions on the shortest length geodesic curve between the intrinsically symmetric points. The restriction of \(\phi _i\) on a curve \(\gamma (t)\) is defined as \(\phi _i\circ \gamma (t):[0,1]\rightarrow \mathbb {R}\). Since, \(\phi _i\circ T_\text {p}=\pm \phi _i, \forall i\) such that i-th eigenvalue is non repeating, each eigenfunction is always either an even (sign\(=+1\)) or an odd (sign\(=-1\)) function. Hence, if restriction of the eigenfunction \(\phi _i\) on the shortest length geodesic between the intrinsically symmetric points has sign \(+1\) (\(-1\)), then the sign of \(\phi _i\) is also \(+1 \;(-1)\).

Proposition 1

Let \(x,y\in \mathcal {M}\) be two intrinsically symmetric points and \(\gamma (t):[0,1]\rightarrow \mathcal {M} \) be the shortest length geodesic curve between the points x and y. Then, \( \phi _i\circ \gamma (t)=\pm \phi _i\circ \gamma (1-t),\forall t\in [0,1].\)

Proof

Using Theorem 2, we proceed as \(\phi _i(\gamma (t))=\phi _i(T_\text {p}(\gamma (1-t)))=(\phi _i\circ T_\text {p})(\gamma (1-t))\). We know that \(\phi _i\circ T_\text {p}=\pm \phi _i\). Hence, \(\phi _i(\gamma (t))=\pm \phi _i(\gamma (1-t))\).    \(\square \)

We apply the above result to find whether an eigenfunction is even or odd as follows. We first determine a set of candidate pairs of intrinsically symmetric points. Then we find the shortest length geodesic curve between each pair. Then, for each eigenfunction \(\phi _i\), we determine if the restricted eigenfunction \(\phi _i\circ \gamma (t)\) is an even or an odd function for each pair.

3.3 Computation of Eigenfunction of Laplace-Beltrami Operator

Let \(\mathcal {T}=(\mathcal {V},\mathcal {F},\mathcal {E})\) be a triangle mesh, where \(\mathcal {V}\) is the set of n vertices, \(\mathcal {F}\) is the set of faces, and \(\mathcal {E}\) is the set of edges. We follow the method in [39] to find the eigenvalues and the corresponding eigenvectors of the Laplace-Beltrami operator in the discrete settings. The discrete Laplace-Beltrami operator is defined by the matrix \(\mathbf {L}=-\mathbf {A}^{-1}\mathbf {M}\). Both \(\mathbf {M}\) and \(\mathbf {A}\) are of size \(n\times n\) and are defined as follows.

$$ \mathbf {M}_{j,j^\prime }={\left\{ \begin{array}{ll} \frac{\cot (\alpha _{jj^\prime })+\cot (\beta _{jj^\prime })}{2} &{}\text { if } (j,j^\prime )\in \mathcal {E}\\ -\sum _{j^{\prime \prime }\ne j}\mathbf {M}_{j,j^{\prime \prime }} &{}\text { if } j=j^\prime \\ 0&{}\text { if } (j,j^\prime )\notin \mathcal {E}, \end{array}\right. } $$
figure a

\(\mathbf {A}=\text {diag}(a_1,a_2,\ldots ,a_n)\), and, \(a_j=\frac{1}{3}\sum _{j^\prime ,j^{\prime \prime }:(j,j^\prime ,j^{\prime \prime })\in \mathcal {F}}a_{jj^\prime j^{\prime \prime }}\) (the area of the shaded region in the inset figure). Here, \(a_{jj^\prime j^{\prime \prime }}\) is the area of the face \((j,j^\prime ,j^{\prime \prime })\). In the discrete settings, we denote eigenfunctions by \(\varvec{\phi }_i,i\in [k]\), and are the solutions of the generalized eigen-problem \(\mathbf {M}\varvec{\phi }_i=-\lambda _i \mathbf {A}\varvec{\phi }_i\). Here, \([k]=\{1,2,\ldots ,k\}\). We denote the value of \(\varvec{\phi }_i\) at the j-th point or vertex by \(\varvec{\phi }_i(x_j)\).

3.4 Detecting Pairs of Intrinsically Symmetric Points

In order to detect the intrinsic symmetry, according to Theorem 1, we need to find the matrix \(\mathbf {C}\) defined as \(\mathbf {C}_{i,i^\prime }=0\), if \(i\ne i^\prime \), \(\mathbf {C}_{i,i}=+1\), if \(\varvec{\phi }_i\) is an even function, and \(\mathbf {C}_{i,i}=-1\), if \(\varvec{\phi }_i\) is an odd function. According to Proposition 1, the eigenfunction \(\varvec{\phi }_i\) is an even function (odd) if its restriction on the shortest length geodesic between intrinsically symmetric points is also an even function (odd). Therefore, our first task is to find a few accurate candidate pairs of intrinsically symmetric points.

Fig. 1.
figure 1

Direct HKS matching vs. restricted HKS matching. (a)–(b) Pairs of intrinsically symmetric points using HKS similarity. It suffers from the fact that if two neighboring points are subjected to a same strength heat source, then their heat diffusions will be similar. (c) Pairs obtained by restricted HKS matching. We observe that the sign of low frequency eigenfunctions on two neighboring points is the same. Therefore, we assign a high cost for pairing two points having the same vectors of signs of eigenfunctions.

We find the heat kernel signature (HKS) feature points [48] on the given mesh \(\mathcal {T}\). HKS feature points are the local maxima of the function \(\sum _{i=1}^{k}e^{-\lambda _i t_h}\varvec{\phi }_i^2(x_j)\) defined on the shape \(\mathcal {T}\) for all vertices \(x_j\), \(j\in [n]\). We use \(k=13\) for our experiments. Let \(\{x_j\}_{j\in \mathcal {I}}\) be the set of HKS feature points, where, \(\mathcal {I}\subset [n]\), and \(|\mathcal {I}|=d\). We set \(t_h=\frac{4\ln 10}{\lambda _2}\) as defined in [48]. To find the pairs of intrinsically symmetric points, we do not directly match the HKS descriptors. The reason is that this strategy may fail very frequently. In Fig. 1(a), we show the detected HKS feature points and the pairs detected by direct matching of their HKS descriptors on the Kids model [9]. In Fig. 1(b), we show the zoomed pairs of symmetric points on the feet and the hands for better visualization. We observe that the tips of two fingers of the same hand got paired. The reason behind getting such matches is that if two neighboring points are subjected to a same strength heat source, then their heat diffusion processes will be similar. Therefore, their HKS descriptors will be similar. Hence, we have to assign a high cost for pairing two neighboring points. Furthermore, determining if two points are neighbors requires us to find the geodesic distance between each possible pair. This can be a costly process since there can be a large number of possible pairs. We propose a fast approach to determine if two points are neighbors based on the following observation.

Observation 1. Let \(x_{j}\) and \(x_{j^\prime }\) be any two neighboring points in \(\mathcal {T}\). Then, \(\mathbf {s}_{j}=\mathbf {s}_{j^\prime }\) for low frequency eigenfunctions, i.e., for small k.

Here, \(\mathbf {s}_{j}\) and \(\mathbf {s}_{j^\prime }\) are the j-th and \(j^\prime \)-th columns of the matrix \(\mathbf {S}=\begin{bmatrix} \mathbf {s}_1&\mathbf {s}_2&\ldots&\mathbf {s}_d\end{bmatrix}\in \{-1,+1\}^{k\times d}\). The i-th element of the vector \(\mathbf {s}_j\) is the sign of the i-th eigenfunction on the j-th point of the set \(\{x_j\}_{j\in \mathcal {I}}\). We give an intuitive understanding of this observation based on the nodal domains of the eigenfunctions. The nodal set of the eigenfunction \(\mathbf {\phi }_i\) is the set \(\mathcal {B}_i=\{x\in \mathcal {M}:\phi _i(x)=0\}\). A nodal domain is a component in the set \(\mathcal {M}\backslash \mathcal {B}_i\). The set \(\mathcal {M}\backslash \mathcal {B}_i\) is the collection of components or segments on the shape which are separated by the set \(\mathcal {B}_i\). The value of the eigenfunction \(\phi _i\) in any of its nodal domain is either positive or negative [43, 59]. Therefore, if two points lie in the same nodal domain, then the eigenfunction \(\phi _i\) will have the same sign on both the points. Now, the neighborliness of two points depends on the size or the area of the nodal domain. According to the Courant’s nodal domain theorem, the number of nodal domains of \(\phi _i\) is less than i [59]. Therefore, the size of the nodal domains remains significantly large for low frequency eigenfunctions. Hence, neighboring points remain in the same nodal domain for all the low frequency eigenfunctions. Therefore, the sign of all low frequency eigenfunctions remains the same on neighboring points. We choose eigenfunctions corresponding to the first 13 lowest eigenvalues and corresponding eigenvectors for our experiment.

Hence, we assign a high cost for pairing the points \(x_{j}\) and \(x_{j^\prime }\), if their sign vectors \(\mathbf {s}_{j}\) and \(\mathbf {s}_{j^\prime }\) are the same. Let \(\mathbf {H}=\begin{bmatrix} \mathbf {h}_1&\mathbf {h}_2&\ldots&\mathbf {h}_d\end{bmatrix}\in \mathbb {R}^{h\times d}\) be the heat kernel signatures matrix of the detected d HKS feature points, where we choose \(h=50\) steps. We define the affinity matrix \(\mathbf {W}\in \mathbb {R}^{d\times d}\) such that \(\mathbf {W}_{j,j^\prime }=\Vert \mathbf {h}_j-\mathbf {h}_{j^\prime }\Vert _2+q\psi (\Vert \mathbf {s}_j-\mathbf {s}_{j^\prime }\Vert _2),\forall (j,j^\prime )\in [d]\times [d],j\ne j^\prime \) and \(\mathbf {W}_{j,j}=q, \forall j\in [d]\). Here \(\psi (t)=0\), if \(t>0\) and \(\psi (t)=1\), if \(t=0\), and q is any large positive constant. We now pair these points such that if \(x_{j^\prime }\) is intrinsically symmetric point of the point \(x_j\), then \(x_{j}\) should be the intrinsically symmetric point of the point \(x_{j^\prime }\). We achieve this by representing the matching by a matrix \(\mathbf {\Pi }\in \{0,1\}^{d\times d}\), where \(\mathbf {\Pi }_{j,j^\prime }=1\) and \(\mathbf {\Pi }_{j^\prime ,j}=1\), if the points \(x_j\) and \(x_{j^\prime }\) form a pair and 0, otherwise. Now, we enforce the constraints \(\mathbf {\Pi }\mathbf {1}=\mathbf {1}\) and \(\mathbf {\Pi }^\top \mathbf {1}=\mathbf {1}\) to achieve one-to-one matching, where \(\mathbf {1}\) is a vector of size d with all elements equal to 1. We get many points which can not be paired. Therefore, we cap the number of pairs by c which we represent by the constraint \(\mathbf {1}^\top \mathbf {\Pi }\mathbf {1}=2c\). Further, to make it feasible, we modify the one-to-one matching constraints to \(\mathbf {\Pi }\mathbf {1}\le \mathbf {1}\) and \(\mathbf {\Pi }^\top \mathbf {1}\le \mathbf {1}\). Now, we frame the problem of pairing the points in the below optimization problem.

$$\begin{aligned} \underset{\mathbf {\Pi }\in \{0,1\}^{d\times d}}{\min }\sum _{j=1}^d\sum _{j^\prime =1}^d\mathbf {\Pi }_{j,j^\prime }\mathbf {W}_{j,j^\prime } \text {, subject to}&\mathbf {\Pi }\mathbf {1}\le \mathbf {1},\;\mathbf {\Pi }^\top \mathbf {1}\le \mathbf {1},\;\mathbf {1}^\top \mathbf {\Pi }\mathbf {1}=2c. \end{aligned}$$
(1)

We note that, problem (1) is equivalent to the below linear assignment problem.

$$\begin{aligned} \underset{\varvec{\pi }\in \{0,1\}^{d^2\times 1}}{\min }\mathtt {vec}(\mathbf {W})^\top \varvec{\pi } \text {, subject to}&\mathbf {C}_1\varvec{\pi }\le \mathbf {1},\;\mathbf {C}_2\varvec{\pi }\le \mathbf {1},\;\mathbf {c}_3^\top \varvec{\pi }=2c. \end{aligned}$$
(2)

Here, the vector \(\varvec{\pi }=\mathtt {vec}(\mathbf {\Pi })\), \(\mathbf {C}_2=\mathbf {1}^\top \otimes \mathbf {I}\), \(\mathbf {c}_3\) is the vector of size \(d^2\times 1\) with all elements equal to 1, and \(\mathbf {I}\) is the identity matrix of size \(d\times d\). The matrix \(\mathbf {C}_1\) is \(d\times d^2\) matrix and defined such that the j-th row is the circular shift, by jd elements on the right, of the row vector of size \(1\times d^2\) with the first d elements equal to 1 and the last \(d^2-d\) elements equal to 0. The time complexity of this problem is exponential in the number of variables. However, the size of our problem is very small. In all our experiments \(d\le 25\). We use the MATLAB function intlinprog to solve this problem which takes \(\approx 0.03\) s.

3.5 Determining the Sign of Eigenfunctions

Proposition 1 states that the eigenfunction \(\varvec{\phi }_i\) is an even (odd) function, if its restriction \(\varvec{\phi }_i\circ \gamma _{j}(t):[0,1]\rightarrow \mathbb {R}\) on the shortest length geodesic \(\gamma _j(t)\) between any two intrinsically symmetric points \(x_{j}\) and \(x_{j^\prime }\) is an even (odd) function. Let \(\{(x_j,x_{j^\prime })\}_{j=1}^{c}\), be the set of detected pairs of intrinsically symmetric points. We find the shortest length geodesic curve between two intrinsically symmetric points using [50] with approximate setting (Dijkstra’s algorithm), since the exact geodesic curve may not pass through the vertices of the mesh which may require us to perform interpolation for calculating the values of \(\varvec{\phi }_i\circ \gamma _{j}(t)\) for \(\gamma _{j}(t)\notin \mathcal {V}\). Let \(\mathbf {p}_{ij}\) be the restriction (vector of size equal to the number of vertices in the geodesic) of the eigenfunction \(\varvec{\phi }_i\) on the shortest length geodesic curve between the intrinsically symmetric points \(x_j\) and \(x_{j^\prime }\). Then, the sign \(s_i\) of eigenfunction \(\varvec{\phi }_i\) is equal to \(+1\), if \(\sum _{j=1}^c\mathbf {p}_{ij}^\top \mathtt {flip}(\mathbf {p}_{ij})>0\) and equal to \(-1\), if \(\sum _{j=1}^c\mathbf {p}_{ij}^\top \mathtt {flip}(\mathbf {p}_{ij})<0\). We do not consider the eigenfunction \(\varvec{\phi }_i\) if \(\sum _{j=1}^c\mathbf {p}_{ij}^\top \mathtt {flip}(\mathbf {p}_{ij})=0\). Equivalently, we define the diagonal entries of the functional correspondence matrix \(\mathbf {C}\) as \(\mathbf {C}_{i,i}=s_i\). In Fig. 2(a) and (c), we show the eigenfunctions \(\varvec{\phi }_2\) and \(\varvec{\phi }_3\), respectively, with the geodesic curves for two pairs of detected intrinsically symmetric points on a shape from Kids dataset [9]. In Fig. 2(b) and (d), we show the functions \(\varvec{\phi }_i\circ \gamma _{j}(t)\) for \(i=2,3\) and \(j=1,2\). We observe that \(\mathbf {C}_{2,2}=-1\) and \(\mathbf {C}_{3,3}=+1\).

One can directly use the values of an eigenfunction on the intrinsically symmetric points to find the sign instead of checking it on the geodesic between these points. However, this approach could be sensitive to the noise. If the value of an eigenfunction at the feature point has changed due to noise, then the point-based method will fail. Whereas, it is less likely that due to the noise the value of an eigenfunction will be changed at all the points on the geodesic. Our geodesic based method will detect the sign correctly due to averaging of signs.

Fig. 2.
figure 2

(a)–(b): The eigenfunction \(\varvec{\phi }_2\) on the Kids model which is an even function and its restrictions \(\varvec{\phi }_2\circ \gamma _j(t)\), \(j=1,2\), shown in red and blue colors (on hands and feet), which are also even functions. (c)–(d): The eigenfunction \(\varvec{\phi }_3\) on the Kids model which is an odd function and its restrictions \(\varvec{\phi }_3\circ \gamma _j(t)\), \(j=1,2\), shown in purple and orange colors (on hands and feet), which are also odd functions. (Color figure online)

3.6 Correcting the Eigenfunctions

The models of the benchmark datasets are obtained by applying an imperfect isometry, so the theory only holds approximately. Furthermore, some of the triangles may not be Delaunay triangles and the eigenfunctions are sensitive to the change in the triangulation of the mesh. Therefore, all the eigenfunctions may not be perfectly even or odd which may give the erroneous symmetry detected in Sect. 3.7. Consider Fig. 3(a), where the eigenfunction \(\varvec{\phi }_{10}\) is not perfectly even on the legs. We transform the eigenfunctions such that they preserve the pairs of intrinsically symmetric functions. We extend the framework in [19]. Let \(\mathbf {\Phi }=\begin{bmatrix}\varvec{\phi }_1&\varvec{\phi }_2&\ldots&\varvec{\phi }_k\end{bmatrix}\in \mathbb {R}^{n\times k}\), and \(\mathbf {D}=\text {diag}(\lambda _1,\lambda _2,\ldots ,\lambda _k)\in \mathbb {R}^{k\times k}\). Let \(\mathbf {\Phi }\mathbf {R}\) be the transformed basis obtained by applying the linear operator \(\mathbf {R}\) on the basis \(\mathbf {\Phi }\). Then, we impose the constraints \(\mathbf {R}^\top \mathbf {D}\mathbf {R}=\mathbf {D}\) and \(\text {off}(\mathbf {R}^\top \mathbf {D}\mathbf {R})=0\) so that the new eigenfunctions admit to the original eigenfunction decomposition problem as proposed in [19], where \(\text {off}(\mathbf {M})=\sum _{j}\sum _{j^\prime :j^\prime \ne j}\mathbf {M}_{j,j^\prime }^2\) for any matrix \(\mathbf {M}\).

Now, let \(f_j,g_j:\mathcal {M}\rightarrow \mathbb {R}\) be two functions such that \(f_j\) and \(g_j\) are intrinsic images of each other. That is, \(f_j\circ T_\text {p}(x)=g_j(x)\) and \(g_j\circ T_\text {p}(x)=f_j(x)\) are equivalent. Let \(\mathbf {f}_j\in \mathbb {R}^{n}\) and \(\mathbf {g}_j\in \mathbb {R}^{n}\) be the discrete versions of \(f_j\) and \(g_j\), respectively. Let \(\mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {f}_j\) and \(\mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {g}_j\) be the representations of the functions \(\mathbf {f}_j\) and \(\mathbf {g}_j\) in the transformed basis \(\mathbf {\Phi }\mathbf {R}\), respectively. We want the transformed basis \(\mathbf {\Phi }\mathbf {R}\) such that \(\mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {f}_j=\mathbf {C}\mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {g}_j\). Let \(\mathbf {F}=\begin{bmatrix}\mathbf {f}_1&\ldots&\mathbf {f}_c&\mid&\mathbf {g}_1&\ldots&\mathbf {g}_c \end{bmatrix}\in \mathbb {R}^{n\times 2c}\) and \(\mathbf {G}=\begin{bmatrix}\mathbf {g}_1&\ldots&\mathbf {g}_c&\mid&\mathbf {f}_1&\ldots&\mathbf {f}_c \end{bmatrix}\in \mathbb {R}^{n\times 2c}\) be the matrices representing 2c (bidirectional) pairs of intrinsically symmetric functions. We formulate the below optimization framework to find the transformation matrix \(\mathbf {R}\).

$$\begin{aligned} \nonumber&\underset{\mathbf {R}}{\min } \text { off}(\mathbf {R}^\top \mathbf {D}\mathbf {R})+\Vert \mathbf {R}^\top \mathbf {D}\mathbf {R}-\mathbf {D}\Vert _\text {F}^2\\ \text {subject to}&\mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {F}=\mathbf {C}\mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {G},\; \mathbf {R}^\top \mathbf {R}=\mathbf {I},\text {det}(\mathbf {R})=+1,\mathbf {R}\in \mathbb {R}^{k\times k}. \end{aligned}$$
(3)

Here, \(\mathbf {R}^\top \mathbf {R}=\mathbf {I}\) follows from the fact that the transformed basis \(\mathbf {\Phi } \mathbf {R}\) is an orthogonal basis. Here, the set \(\{\mathbf {R}\in \mathbb {R}^{k\times k}:\mathbf {R}^\top \mathbf {R}=\mathbf {I},\text {det}(\mathbf {R})=+1\}\) is the special orthogonal group \(\mathcal {SO}(k)\). Hence, we solve the below optimization problem.

$$\begin{aligned} \underset{\mathbf {R}\in \mathcal {SO}(k)}{\min }\text { off}(\mathbf {R}^\top \mathbf {D}\mathbf {R})+\Vert \mathbf {R}^\top \mathbf {D}\mathbf {R}-\mathbf {D}\Vert _\text {F}^2+\mu \Vert \mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {F}-\mathbf {C}\mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {G}\Vert _\text {F}^2. \end{aligned}$$
(4)
Fig. 3.
figure 3

Visualization of the eigenfunction correction.

We use the Riemannian-trust-region method, proposed in [1, 2], to solve this optimization problem. We use the manopt toolbox [7] for this purpose. We provide the Riemannian gradient and Hessian of this cost function in the supplementary material. We empirically found the optimal \(\mu \) to be equal to 1 in our experiments. We choose the functions \(\mathbf {f}_j\) and \(\mathbf {g}_j\) such that, \(\mathbf {f}_j=1\) at the point \(x_j\) and 0 everywhere else, and \(\mathbf {g}_j=1\) at the point \(x_{j^\prime }\) and 0 everywhere else. Here, \(x_j\) and \(x_{j^\prime }\) are intrinsically symmetric points. In Fig. 3(a) and (b), we show the effect of correction on (\(\varvec{\phi }_{10}\)). We observe that \(\varvec{\phi }_{10}\), which was not perfectly symmetric on the legs and the belly, becomes more symmetric. The large blue patch on the belly also got moved to center which was more towards left before correction. Here, the value of eigenfunction is color encoded, more blue implies more negative and more yellow implies more positive. In Fig. 3(c), we show the restriction of \(\varvec{\phi }_{10}\) on the geodesic between the two symmetric points which becomes more symmetric after the correction.

3.7 Dense Intrinsically Symmetric Correspondence

Let \(\mathbf {f}_{j}\) be the function such that \(\mathbf {f}_{j}=1\) at \(x_{j}\) and 0 elsewhere. Similarly, let \(\mathbf {g}_{j^\prime }\) be the function such that \(\mathbf {g}_{j^\prime }=1\) at \(x_{j^\prime }\) and 0 elsewhere. Let \(\mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {f}_{j}\) and \(\mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {g}_{j^\prime }\) be their basis representation. Then, if the point \(x_{j}\) and \(x_{j^\prime }\) are intrinsically symmetric then \(\mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {f}_{j}=\mathbf {C}\mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {g}_{j^\prime }\). Which is equivalent to \(\mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {F}=\mathbf {C}\mathbf {R}^\top \mathbf {\Phi }^\top \mathbf {G}\) if we consider all points. Now, if \(\mathbf {F}\) is equal to the identity matrix of size \(n\times n\), then \(\mathbf {G}_{j,j^\prime }=1\), if point \(x_j\) and \(x_{j^\prime }\) form a pair of intrinsically symmetric points, and 0 otherwise. Now following [36], the intrinsically symmetric point of \(x_j\) is the nearest neighbor of the j-th column of the matrix \(\mathbf {R}^\top \mathbf {\Phi }^\top \) among the columns of the matrix \(\mathbf {C}\mathbf {R}^\top \mathbf {\Phi }^\top \). The obtained correspondences are continuous as shown in [36]. Our method is invariant to the ordering of the eigenfunction since the sign of \(\varvec{\phi }_i\) and \(\mathbf {C}_{i,i}\) only depend on the eigenfunction \(\varvec{\phi }_i\). In Fig. 3(d), we show the detected symmetry on a geodesic on legs before and after correction.

Table 1. The total time for computing intrinsic symmetry for the methods MT [16], BIM[17], OFM [24], GRS [52], and the proposed approach on the SCAPE dataset [3].
Table 2. The correspondence rates and mesh rates for the methods MT [16], BIM [17], OFM [24], GRS [52], and the proposed approach on the SCAPE dataset [3].
Table 3. The correspondence rates and mesh rates for the methods MT [16], BIM [17], OFM [24], GRP [52], and the proposed approach on the TOSCA dataset [8].

4 Results and Evaluation

4.1 Time Complexity

Let n be the number of vertices and k be the number of eigenfunctions used. The feature points are the local maximums of \(\sum _{i=1}^{k}e^{-\lambda _i t_h}\varvec{\phi }_i^2(x_j)\). It requires us to find 2-ring neighborhoods of each vertex. We use the half-edge data structure which requires O(1) time. Hence, the overall time for finding the feature points is O(n). The optimization problem in Eq. (4) takes \(O(nk^2)\) when solved using Riemannian trust region method. We use the ANN library [4] to find the nearest neighbor for each column of the matrix \(\mathbf {R}^\top \mathbf {\Phi }^\top \in \mathbb {R}^{k\times n}\) among the columns of the matrix \(\mathbf {C}\mathbf {R}^\top \mathbf {\Phi }^\top \in \mathbb {R}^{k\times n}\) which takes time \(O(kn\log (n))\). Hence, the time complexity is \(O(kn\log (n))+O(n)+O(nk^2)\approx O(kn\log (n))\), since \(k<<n\). In our experiments \(k=13\) (empirical) and \(n\approx 15000\). The time complexity of computing the k smallest eigenvalues and corresponding eigenvectors of symmetric matrix is \(O(n^2k)\) which is common to all spectral decomposition based methods.

4.2 Comparison

Evaluation Metrics. We use the following evaluation metrics to compare the results of our method to that of the state-of-the-art methods as defined in [16]. Correspondence rate: Let \((x_j,x_{j^\prime }^\text {g})\) be the ground truth correspondence and \((x_j,x_{j^\prime }^\text {e})\) be the estimated correspondence, then the correspondence \((x_j,x_{j^\prime }^\text {e})\) is called true positive if the geodesic distance between the points \(x_{j^\prime }^\text {g}\) and \(x_{j^\prime }^\text {e}\) is less than \(\sqrt{area(\mathcal {T})/20\pi }\) as used in [16]. The correspondence rate is the fraction of true positive correspondences in the total estimated correspondences. Mesh rate: The mesh rate is the fraction of shapes for which the correspondence rate is more than \(75\%\) in the total shapes as used in [16]. Time Complexity: Total time required for computing symmetry for each shape in the given dataset. Datasets. We evaluate our approach on the SCAPE [3] and TOSCA [8] datasets. The SCAPE dataset contains 71 models. Each model in SCAPE dataset contains 12500 vertices and 24998 faces. The TOSCA dataset contains 80 models. On an average 20 ground truth intrinsically symmetric correspondences provided for each model in the datasets SCAPE and TOSCA. In the Fig. 4, we show a few results of the proposed approach on both the datasets. We have only shown the sparsely detected correspondences for better visualization.

Fig. 4.
figure 4

Results of our approach on the TOSCA [8] (first row) and the SCAPE [3] (second row) datasets. Detected correspondences (sparse) are shown in blue color. Correspondences in red color are the ones detected in Sect. 3.4.

Comparison Methods. We compare the results of our approach on the datasets SCAPE and TOSCA with the four methods Möbius transformation voting (MT) [16], Blended Intrinsic Maps (BIM) [17], Properly Constrained Orthogonal Functional Map (OFM) [24], and Group Representation of Symmetries (GRS) [52].

Discussions on the Comparison. In Table 1, we present the total time required for detecting the intrinsic symmetry in all the models of the TOSCA dataset for all the methods. We observe that our method is the fastest method on the TOSCA dataset. Our method takes around 6 s for each model whereas the method BIM takes around 270 s, the method OFM takes around 45 s, and the method GRS takes around 18 s. Our method takes 4.2 min to compute intrinsic symmetry in all the models of the SCAPE dataset. The possible reasons for our faster computation include finding the correspondence matrix using a closed form solution and determining the sign of eigenfunctions by computing the approximate shortest length geodesic curves between two intrinsically symmetric points. In Tables 2 and 3, we present the correspondence rate (Corr rate) and the mesh rate for all the methods for all the models of the SCAPE and the TOSCA datasets, respectively. The mesh rate for our method is equal to 100% and the correspondence rate is equal to 97.8% which is very close to the state-of-the-art correspondence rate 98% of the method [17]. However, the average computation time for each mesh is around 270 s for the method [17], whereas it is around 6 s for our method. We achieve the state-of-the-art performance on the SCAPE dataset.

Fig. 5.
figure 5

Partial intrinsic symmetry detection results on the dataset SHREC16 [9].

Effect of Holes. In Fig. 5, we show the detected intrinsic symmetry in the partial model from the SHREC16 [9] dataset. Here, the partial shape is obtained by making holes in the original shape such that it contains 90% area of the main shape. We observe that our method is invariant to significant holes.

5 Conclusions

We have presented a fast and an accurate algorithm for detecting intrinsic symmetry in triangle meshes. We showed that the functional correspondence matrix is diagonal and a diagonal entry is \(+1\; (-1)\) if the corresponding eigenfunction is even (odd). We showed that the restriction of an even (odd) eigenfunction on the shortest length geodesic between any two intrinsically symmetric points also is an even (odd) function. This result has helped us to derive a closed form solution to find diagonal entries of this matrix. We achieved state-of-the-art performance on the SCAPE dataset and second best on the TOSCA dataset. We achieved the best time complexity. Furthermore, our approach is invariant to the ordering of eigenfunctions and robust to the presence of holes in the input mesh. Our method is limited to the intrinsic reflective symmetry. It can not find the other types of symmetries such as rotational symmetry. We would like to extend our approach to more general symmetries. Our approach may fail to detect intrinsic symmetry in non-connected manifolds. As future work, we would like to extend the functional map to detect intrinsic symmetries in non-connected manifolds.