Evaluation Algorithm of Root Canal Shape Based on Steklov Spectrum Analysis

In recent years, we have seen more and more interest in the field of medical images and shape comparison motivated by the latest advances in microcomputed tomography (μCT) acquisition, modelling, and visualization technologies. Usually, biologists need to evaluate the effect of different root canal preparation systems. Current root canal preparation evaluation methods are based on the volume difference, area difference, and transportation of two root canals before and after treatment. The purpose of root canal preparation is to minimize the volume difference and ensure the complete removal of the smear layer. Previous methods can reflect some general geometric differences, but they are not enough to evaluate the quality of root canal shape. To solve this problem, we proposed a novel root canal evaluation method based on spectrum and eigenfunctions of Steklov operators, which can be served as a better alternative to current methods in root canal preparation evaluation. Firstly, the ideal root canal model was simulated according to the root canal model before and after preparation. Secondly, the Steklov spectrum of the two models was calculated. Thirdly, based on the spectrum and the histogram of the Gaussian curvature on the surface, the weight of each eigenvalue was computed. Therefore, the Steklov spectrum distance (SSD), which measures shape difference between the root canals, was defined. Finally, the calculation method that quantifies the root canal preparation effect of root canals was obtained. Through experiments, our method manifested high robustness and accuracy compared with existing state-of-the-art approaches. It also demonstrates the significance of our algorithm's advantages on a variety of challenging root canals through result comparison with counterpart methods.


Introduction
In recent years, imaging technology represented by 3D digital geometry technology has widely been used in the medical field. e 3D model can display critical information such as the three-dimensional topological structure, geometric information, and anatomical structure of the lesion more intuitively. It lays a solid foundation for successive operations, such as simulating the operation, humancomputer interaction, quantitative morphometry, and telemedicine [1]. Researchers have paid more and more attention to the issues of segmentation, registration, difference evaluation, and classification retrieval of 3D medical models. With the popularity of measurement tools represented by μCT, the evaluation of anatomical structure before and after root canal preparation has become a hot topic in medical imaging [2][3][4][5][6][7][8][9]. e root canal is the part of the cavity in the middle of the tooth. e upper part of the cavity is broad, called the pulp cavity, and the lower part has a tubular structure from which the blood vessels of the dental nerve and the nutrient nerve are derived, as shown in Figure 1. In [10], Vertucci investigated 2,400 permanent teeth and categorized them into eight types. He studied the morphology and root canal of adult and established the evaluation criteria of treatment difficulty. From the point of view of dental anatomy, root canal cleansing and complete obturation with an inert filling material are very necessary when the root canal morphology of teeth is affected. Root canal therapy generally consists of three stages: root canal preparation, root canal disinfection, and root canal obturation. Root canal preparation is the critical step for root canal obturation therapy, including root canal cleaning and formation. Before root canal preparation, the location, length, shape, and number of root canal orifices need to be determined. Under the guidance of the dental operating microscope, ultrasonic root canal files and trumpet hand files are used to probe the root canal to the apical orifice. Root canals with small or blocked root canals should be cleaned and dredged with ethylenediaminetetraacetic acid (EDTA) or other solutions. Infected dentin should be removed and sterilized by root canal instruments to facilitate root canal obturation. e prepared root canals should contain the original root canal shape, retain the narrow part of the apex, and form the taper from the coronal to apical [3].
To improve the effect of root canal therapy, various root canal preparation systems have been developed to maintain the original morphology of root canal anatomy (including using different instruments and cooperating with different solutions to eliminate the smear layer of the root canal wall).
e results before and after the operation among groups were tested by the significance test, such as Mann-Whitney U test [3,6,7,9]. With the development of medical imaging measurement technology such as CT, more and higher resolution image data can be collected, which boost more accurate instrument evaluation. μCT is a nondestructive 3D imaging technology, which can clearly understand the internal microscopic structure of the sample. Hammad et al. [11] compared various instruments and considered that the root canal measurement results of the μCT are consistent with the actual results. Kierklo et al. [6] evaluated the root canal preparation operation based on 2D data; that is, the cross-sectional images collected by the μCT were compared before and after preparation. In [12], 3D models were reconstructed from images, and the volume and area of the model were used to evaluate the effect of preparation among different systems.
In the field of medical imaging, many studies of therapy evaluation are based on the idea of registration before comparison. As far as root canal models are concerned, there were many changes in the models before and after treatment, so the rigid registration method is unsuitable. e nonrigid registration method mainly includes two approaches based on geometric features or image statistics. Both of the two ideas need to extract features or statistical information from the model for registration. However, in most cases, it is difficult to obtain features accurately when the feature extraction is carried out between different models, or when the feature points are not distinct, and always accompanied by outliers. e outliners probably present immense barriers to the subsequent registration algorithm. Moreover, the geometric characteristics and the topology of the root canal model may change before and after preparation. It is also challenging to apply the idea of nonrigid registration to root canal preparation assessment.
Ørstavik [3] summarized the evaluation methods of root canal shaping ability of various instruments, including X-ray parallel projection, Bramante model, cross-sectional section method, CT, and μCT methods. e ideas of these methods are to obtain the volume, surface area, and transportation of root canal axis before and after preparation by integrating the slice images [13,14]. Specific methods are as follows: (1) Divide teeth with high h into N equal parts at the interval of δ � h/N, then get images I i N 1 (2) Measure the cavity area S i and the perimeter of the cavity of the root canal l i in I i (3) Calculate the total surface area with S � N− 1 i�1 δl i (4) Calculate volume with V � N i�1 δS i (5) e deviation of the root canal center and the average distance of each cross-sectional cavity are calculated, respectively, after registration, as shown in Figure 2 Existing methods of root canal preparation evaluation have obvious drawbacks. First, the differences in permeability and taper of different root canals were submerged behind a single geometric quantity, which cannot distinguish the root canals with similar volume and significant shape difference (shown in Figure 3). Second, the overall geometry and complex topological structure of root canals are not taken into account. ird, the cavity morphology before and after treatment is often not considered. A large extent of transformation (irregular to ellipsoid-like) occurs, which may cause intolerable errors in both rigid and nonrigid registration algorithms.
Every root canal has its own individual shape. Accurate knowledge of the canal anatomy is a prerequisite for effective endodontic treatment [15,16]. At present, the evaluation of the various instruments in root canal therapy is generally based on reexamining the volume percentage of the void (POV) of the root canal. is method typically requires a more extended period and cost. erefore, by measuring the model of root canal preparation, defining the "distance" between the shapes and calculating the distance between the current root canal and the idealized root canal according to the criteria conducive to the subsequent root canal obturation, it will be helpful to improve the efficiency of this work.
For the first time in [17], the method of studying geometric shape based on the spectrum was proposed. e paper answered the question of how much about the shape can be inferred from the knowledge of its geometry spectrum. By calculating the shape spectrum λ i and its eigenfunction ϕ i as shown in Figure 4, the shape surface can be uniquely determined and vice versa in ordinary cases.  In recent years, with the development of 4G/5G technology and CT technology, as well as the unprecedented expansion of the accuracy and breadth of medical 3D data acquisition, more and more researchers pay attention to the massive data of medical geometry. On the contrary, the shape analysis paradigm, theory, and practice based on the geometric spectrum have gradually become a research hotspot [18][19][20][21]. Based on the works, the author proposed a novel method based on the geometric spectrum for root canal preparation treatment. e spectrum of differential operator contains rich geometric and topological information of the shape, including curvature, area, and genus. e LBO spectrum and Steklov spectrum have the characteristics of independent of the coordinate system and scale. Although in rare highdimensional manifolds, there may be different shapes sharing the same spectrum. However, in low-dimensional manifold cases (such as root canals), the geometry of the shape can be uniquely determined by its spectrum. e traditional shape spectrum is based on the intrinsic Root canal ϕ 1 ϕ 2 ϕ 3 ϕ 10 ϕ 50 description of the Gauss curvature, which has an insufficient perception of concavity, convexity, and volume of shape, while the application of shape-based extrinsic geometric spectrum, such as Steklov spectrum, is increasing. It has gradually become a research hotspot [22][23][24][25][26].
In this work, we consider the problem of the root canal shape evaluation and investigate its relevance in a selection of literature from medical image and graphics. We propose a shape-aware evaluation method which is a precise evaluator for comparing different root canal preparation systems. Our key contributions can be summarized as follows: Firstly, despite being highly nonlinear and hard to compute, we indicate for the first time that the mapping between a 3D model of the root canal and its Steklov spectrum is addressable with modern mathematical tools Secondly, for the first time, an evaluation method of root canal shape based on an ideal root canal model is proposed, and an algorithm for generating an ideal root canal model is given At last, we evaluate the root canal preparation effect by defining the reasonable distance depending on the histogram of the Gaussian curvature and the Steklov spectrum of root canal models

Related Research
In recent years, more and more attention has been paid to the evaluation method of therapy based on 3D model reconstruction from data collected by μCT. e highresolution characteristics of μCT make it widely used in small-scale medical in vivo research. In [27], synchrotron radiation microcomputed tomography (SRμCT) was evaluated concerning the 3D morphology of the rat spinal cord microvasculature. In [28], Wang et al. used parallel programming techniques to accelerate existing image reconstruction algorithms for 3D optoacoustic tomography. In [29], Kobayashi-Velasco et al. compared conebeam computed tomography (CBCT) with periapical radiograph (PR) in the diagnosis of alveolar and root fractures.
Also, more and more applications of μCT in the evaluation of root canal therapy have been made. e purpose of root canal preparation is to clear the infected root canal anatomical space. Under ideal conditions, on the premise of keeping the centerline unchanged, the inner wall of the root canal should all be processed, and a better taper can be obtained to facilitate subsequent obturation. According to the different geometry of root canal preparation, the evaluation methods of root canal preparation can be divided into two categories. e first category focuses on volume and area changes before and after root canal preparation. In [10], the root canal was segmented into the coronal third, middle third, and apical third. Scanning electron microscopy (SEM) was used to observe the removal of the smear layer; then the images before and after root canal preparation were calculated and compared, which was used to evaluate the therapy. In [4], the geometry of the root canals was described by three different mathematical models: the elliptical model, the 1-ring model, and the 3-ring model, and then a method for analyzing the geometric characteristics of root canals was realized. In [30], Berutti et al. compared the canal curvature and axis transportation after instrumentation with three different systems. In [31], the critical section images of the coronal third, middle third, and apical third before and after root canal preparation were analyzed and compared by several nickeltitanium systems. Anbu et al. [32] observed the efficacy of various technique to fill root canals using spiral computed tomography (SCT). Orhan et al. [2] used both μCT and nano-CT to measure and analyze the volume ratio of filling voids at coronal third, middle third, and apical third and recommend that a voxel size of 11.2 μm can be used as a reliable cutoff value in μCT imaging. e second category of the evaluation method focuses on the transportation of spatial centerline and the balance of periodontal dentin thickness. In [9,16], root canal deviation and central positioning (or transportation) ability were considered to be essential parameters for evaluating the quality of canal preparation. In [33,34], the thickness of the mesial and distal root canal wall was scanned and measured, respectively, and then the displacement and axial centrality of root canals in each group were calculated by using Gambill's formula. In in vivo cases, the coordinate systems of canal measured before and after the operation are not consistent in general. How to register the two root canals remains a problem to be solved. Because of the particularity of teeth, it is not feasible to use registration points in industrial scenes for registration.
ere are two shortcomings in the above two methods: it significantly relies on knowledge and experience manipulators on understanding the geometry of the root canal. e accuracy and automation of the process are not high; the other is that the evaluation based on local sampling lacks a global shape assessment method suitable for subsequent obturation, so it is not conducive to accurately evaluate the quality of root canal preparation of various instruments.
3D shape comparison is a hotspot in the field of geometric analysis. In the field of manufacturing, it is often necessary to evaluate the accuracy of parts processing by shape comparison, such as blades and specialshaped bolts. In the medical field, the manufacture of porcelain fused to metal teeth and artificial bones also need to evaluate the difference of shape. In online massive shape retrieval scenarios, it is necessary to quickly match one or some kinds of shape according to certain retrieval conditions. In the interactive game industry, it is required to cluster and understand facial expressions. In [35], the mapping between two surfaces is expressed by the linear function map, and the difference between different shapes is evaluated by defining reference shapes. Fish et al. [36] introduced a meta-representation that represents the Computational and Mathematical Methods in Medicine essence of a family of shapes. Handels and Hacker [37] extended the volume-based representation of organs to enable the 3D visualization of organs' shape variability in the atlas. In [38], a functional map between two shapes based on the multiscale Laplace-Beltrami operator (LBO) shape descriptor was proposed, and the difference between shapes was evaluated by the basis function family of the spectrum to avoid the NP-hard problem of the registration algorithm. In [39], the comparative shapes were embedded in canonical domains and evaluated.
Shape comparison methods can be categorized into three types: first, the corresponding relationship between shapes is established, and then the geometric comparison is made; the other is to compress the geometric and topological information of three-dimensional image into a scalar, i.e., shape similarity, by using some global or local shape descriptors; and the third is to evaluate the shape difference by defining a specific norm from the global perspective, such as [40]. ese shape descriptors are all shape descriptors based on intrinsic geometry or geodesic equidistant transformation invariant. However, the intrinsic description of shape, such as the LBO spectrum, is unable to distinguish concave and convex shapes ( Figure 5). e conicity, concavity, and convexity of root canal shape and the overall permeability are the key factors to evaluate the effect of treatment. erefore, the extrinsic shape descriptor, which can capture volume difference, is suitable in the scenario. According to [41,42], the spectrum of the Steklov operator mapping the Dirichlet boundary to Neumann boundary is precisely extrinsic. e problem of calculating the spectrum of the Steklov operator of shape is generally transformed to solving the boundary value problem of the partial differential equation (PDE). Because of the complexity and discontinuity of the boundary, it is almost impossible to obtain the analytic solution.
erefore, finding the weak solution of the equation is the only feasible choice. Solving the PDE on the discrete point cloud model by the spectral method is an effective way to improve the accuracy and convergence speed. e basic idea is to map the original PDE to an integral equation. By using the orthogonal and recursive formulas of the basis function, the weak solution problem of PDE is simplified to an integral equation and then converted to the eigenvalue problem of linear equations [43]. erefore, the main obstacle is to select a fast convergent of the basis function family to calculate the shape spectrum for the root canal as a specific boundary condition.
e main contribution of this study is to overcome the shortcomings of the previous methods of evaluation of root canal preparation by surface area, volume, or eccentricity thickness of the central line. Based on strict mathematical theories, the proposed methods for comparing shape differences are given in terms of conicity, permeability, and smoothness of the whole root canal. is method does not require the establishment of calibration points for the shape of the root canal.
is thought avoids the precondition that two shapes must share the same topological structure in previous algorithms and is very suitable for the shape comparison of the root canal. It also breaks through the problems of different CT sampling rates and inconsistent coordinates of in vivo data acquisition and evaluates root canal preparation more accurately and comprehensively.

Theoretical Analysis
anks to the excellent characteristics of the Steklov spectrum in describing concave and convex shapes and volume changes and its robustness in the calculation, the authors propose an algorithm of root canal treatment evaluation based on the shape difference of the Steklov spectrum. e overall idea is as follows: (1) Firstly, the 3D models before and after root canal preparation M 0 and M 1 were acquired by μCT.
(2) To synthesize the changes of a cross section of M 0 and M 1 at the same height, the key cross section was selected, and the ideal root canal model M 1 was generated by lofting. (3) Calculate the Steklov shape spectrum of M 1 and M 1 , and then get the eigenvalues λ i and μ j . (4) Define the spectral distance Steklov spectrum discrete (SSD) to evaluate the difference. e distance is used to measure the shape difference between M 1 and M 1 so as to evaluate the quality of root canal preparation.

Mathematical Model of the Steklov Spectrum.
In the three-dimensional Euclidean space, the Steklov eigenvalue problem is defined as follows.
Let Ω ∈ R 3 be a compact Riemannian manifold, whose boundary Σ � zΩ satisfies the Lipschitz condition. e Steklov equation is used to find the harmonic function u(x) and real number λ ≠ 0. e Steklov problem on Ω is is the Laplace operator on functions on Ω and zu(x)/z n ⟶ is the outward derivative along Σ. In the second equation in equation (1), it can be interpreted as an eigenvalue problem that transforms the Dirichlet boundary condition into the Neumann boundary condition. It maps function : where D is the elliptic pseudodifferential operator of order one.
Hf is the harmonic extension of f on Ω and extends smoothly from the boundary to the interior of Ω [44,45]. e eigenfunctions ϕ 1 , ϕ 2 , . . . , corresponding to the eigenvalues constitutes the orthogonal basis in L 2 ( ) space. e trace of the eigen matrix formed by ϕ i is closely related to the intrinsic geometry and the extrinsic geometry of the shape itself.
Let Tr(e − tD ) be the heat kernel of the operator e − tD : f ⟶ u for the heat equation z t u + Pu � 0, when t ⟶ 0 + , and there is an asymptotic expansion of Tr(e − tD ) [46]: e coefficients a 0 , a 1 , . . . , a n− 1 , as well as b 0 , b 1 , . . . , b n− 1 , are the Steklov heat invariants which are determined by the Steklov spectrum. According to [26], when n � 3 and Ω ∈ R 3 , where K H (x) is the mean curvature at x.
where R M (x) is the Gaussian curvature of boundary Σ. From equations (3) to (6), it can be seen that the mean curvatures reflecting the intrinsic geometry and the Gaussian curvature of the extrinsic geometry are determined by the Steklov spectrum.
In summary, by calculating the eigenvalues and eigenfunctions of the Steklov operators, the intrinsic geometry and extrinsic geometry of shape can be transformed into the spectrum space for research, and then shape difference evaluation can be carried out. Our thought is based on PDE, which obtains the whole geometric structure from the boundary, which is a powerful mathematical tool for medical CT imaging and geological exploration imaging [47].
According to the literature [48], the Steklov spectrum and LBO spectrum of shape are invariant in isometric transformation. Moreover, it has strong robustness for model downsampling, and even nonuniform sampling has better immune effect for topological noise. What is more, the conjugate gradient descent method can be used to solve the above spectrum problem, and the running speed of the algorithm does not depend significantly on the number of edges and faces. e Steklov spectrum has good accuracy even for degenerated irregular meshes. erefore, it is very beneficial to realize the fast evaluation of root canal shape difference.

Discussion on the Numerical Solution of the Steklov Spectrum.
e numerical solution of PDE boundary value problems for high-dimensional manifolds, and two-manifold with higher genus often falls into local minima or diverges, leading to calculation failure. A series of thresholds or carefully adjusted penalty functions were often needed to intervene manually. e 3D topology and geometry of root canals are relatively fixed (tree-shaped cavity), gradually thinning from top to bottom. erefore, it is possible to avoid the shortcomings of the general algorithm and improve the automaticity and accuracy of the algorithm by thoroughly combining this priori knowledge.
us, the design of a particular shape difference evaluation algorithm based on relatively fixed root canal topology will be based on Computational and Mathematical Methods in Medicine the conicity, permeability, and smoothness of the whole root canal shape. e shortcomings of the previous local and single geometric evaluation method were broken through to provide a more scientific and comprehensive assessment method for root canal preparation effect. e calculation method is deduced as follows.
e numerical solution is based on potential theory. e potential theory is a robust basis for solving mathematical and physical equations by using the integral equation method. By using the linear superposition principle and fundamental solutions, various forms of integral expressions are obtained. en, by solving the boundary, the purpose of solving the definite solution problem of the differential equation corresponding to the physical problem is achieved.
In this paper, the boundary element method is used to solve the Steklov eigenvalue problem. is method is different from the finite element method, the variational method, and the radial basis function method. It has the characteristics of simplicity, fast convergence, and high accuracy. It is very suitable for the case of complex internal structure and gridding [49]. e basic idea of the boundary element method is that first, the differential equation is transformed into a boundary integral equation. After that, the discrete quadrature method is used for numerical calculation. In solving the discrete matrix elements of the integral equation, the solution needs to be obtained by a singular integral. In this paper, after derivation, the asymptotic expansion of the odd order can be obtained. When the accuracy can be achieved, the convergence speed of the eigenvalue problem can be significantly improved. At the same time, a posterior error estimate can be obtained based on residual estimation [50].
According to the potential theory of the differential equation, when y is a nonsingular point on Σ, the Steklov eigenvalue problem can be converted to the following integral equation with the Green formula [51]: where v(x) � zu(x)/zn x and U * (y, x) � 1/(4π|x − y|) is the fundamental solution of the following equation: us, the Steklov eigenvalue problem is transformed into the boundary element numerical calculation problem. e first term on the right-hand side of equation (7) is the double-layer potential, the second term is the single-layer potential, and the corresponding density functions are − u(y) and v(x).
According to the continuity of single-layer potential on the boundary, equation (7) still holds for x ∈ Σ. According to the discontinuity of the single-layer potential derivative and the nature of double-layer potential, the following relations can be obtained: where V : H 1/2 (Σ) ⟶ H 1/2 (Σ) is the single-layer potential operator.
where U * (y, x) � 1/(4π|x − y|) is the fundamental solution of Δu(x) � 0. e integral operator V maps the density function w(x) to the potential distribution on Σ.
K is the double-layer potential as follows: Since the right-hand side of equation (9) is not necessarily a positive definite matrix after being discretized into a matrix, the conjugate gradient method cannot be used to calculate its eigenvalues. According to [51], equation (9) can be transformed into another eigenvalue problem with the same solution through the Calderón projection: σ(x) ≔ lim ε⟶∞ 1 4πε 2 y∈Ω:|y− x|�ε ds y , for x ∈ Σ, for almost x ∈ Σ, σ(x) � 0.5 [51].
In equation (12), D is a hypersingular operator which is defined as K * : H − 1/2 (Σ) ⟶ H − 1/2 (Σ) is the adjoint double-layer potential, which is defined as the conormal derivative of V: e potential functions discussed above are integral operators, which are defined in both continuous space and infinite-dimensional function space. However, the root canal model is discretized, so these boundary operators should be further discretized into corresponding matrices which is discussed in Section 4.

Constructing an Ideal Root Canal Model as a Benchmark.
e ideal model serves as a benchmark for the evaluation of preparatory effect. For clearing the inner wall of the root canal, facilitating subsequent obturation and minimizing the transportation of centerline of the root canal, the prepared root canals should have the following criteria: (1) e central axis should be consistent with the original one (2) e new cross section should be enlarged evenly on the basis of the old one to ensure that all the inner walls are treated evenly (3) e volume of the cut inner cavity should be as small as possible on the premise of satisfying Step (2) (4) e transition between sections should be as smooth as possible Based on the four criteria, the second result is better than the first result in Figures 6(d) and 6(e).
Since most of the sections are nearly circular, we propose an ideal algorithm of root canal model generation based on the ideal inner diameter. e constructing steps are as follows: (1) Assuming that the height of the root canal is h, a root canal is divided into three segments: coronal, middle, and apical from top to bottom with heights h c , h m , and h a , respectively. For the coronal segment, three sections are intercepted at three heights of (h (2) To calculate the ideal inner diameter of the section, the following formulas are used: where A(·) represents the area of the given section. (3) For each section, the average distance d(·) between the actual boundary and the ideal boundary is calculated by using the following formula: where p i denotes the ith boundary point of the current section and c(·) is the barycentric coordinates of the boundary points of each section. N represents the number of boundary points of the current section. For each segment of the root canal, the cross section with the minimum d(·) is selected as the representative section of the segment. Obviously, when the section is circular, its d(·) equals zero exactly.

Spectrum Numerical Solution Algorithm.
e potential operator is defined in infinite-dimensional function space, and the actual root canal model is a triangular mesh composed of finite points. In this paper, the Galerkin method is used to convert the continuous equation to weak form. Specifically, by defining the basis function, some constraints are imposed on the original function space [42]: finding ϕ ∈ H d (Σ) in function space H, such as that for all u ∈ H d (Σ), a(ϕ, u) � f(u) ∈ R, where a(·, ·) represents a bilinear form, and it is also a bounded linear function. erefore, for the weak form of single-layer u n , it is expressed by a set of basis function e 1 , e 2 , . . . , e n : u n � n j�1 a e i , e j u j . (18) us, the original integral problem on the whole boundary is transformed into the issue of solving weak form linear equations.
Any piecewise linear hat function ϕ i : Σ ⟶ R on M equals one at their associated vertex and zero at all other vertices and is orthogonal at different subscripts. erefore, any square integrable function u on M can be expressed as Assume that w(x) ∈ Σ is the solution of en, the weak form of equation (20) is After Σ is discretized into the mesh M, w and f are also discretized into vectors , so the following linear system can be obtained: Similarly, the potential integral equations in equation (12) can be listed as a series of linear systems, namely, where d ij is defined as Computational and Mathematical Methods in Medicine and n x and n y denote the normal directions. Also in equation (24), k ij is defined as where k * ij is defined as From the above discussion, we then obtain that the Steklov eigenvalue problem on discrete meshes is e matrices D, K, K * , V, and G in equation (28) are symmetrical, and the summation and inverse matrices are symmetrical and positive definite matrices, which can be solved by using the conjugate gradient algorithm. Firstly, we rewrite equation (28) as follows: In this paper, the PCG method is used to solve equation (29). e main idea is to transform the linear system (29) into equivalent functional optimization problems. In view of the convergence of the steepest descent method and the insufficiency of the conjugate gradient method depending on the distribution of eigenvalues, the incomplete Cholesky decomposition conjugate gradient (ICCG) method not only keeps the matrix sparse but also makes the preoptimized matrix approximate identity matrix. erefore, the method accelerates the iteration process [52].
Assume that L T L is an approximate decomposition of A: where R is a residual matrix and can be used to control the sparsity of L. Since the elements of L should ensure that the preoptimal matrix is effective, L T L is an approximate decomposition of A. is paper adopts the following strategies.
Referring to [53], if the element a ij of A is not zero, then the element l ij of L is calculated with the complete Cholesky decomposition formula. Otherwise, if a ij � 0, then let l ij � 0. After finding the incomplete decomposition of A (let L be nonsingular), A � L T L can be taken as the preoptimal matrix. Since L is not singular, A is symmetrically positive definite. us, the process of solving the linear system AU � 0 can be converted to solving the following two triangular matrix linear systems: In this work, heuristic strategies were used to detecting the convergence rate at different resolutions to determine the optimal sparse degree of the direction. e specific algorithm is as follows.

Root Canal Shape Evaluation Based on the Spectral
Difference. In many kinds of literature, shape difference evaluation based on landmark points is the primary method, including geodesic distance based on point pairs and spectral distance based on point pairs. ese methods require that there are well-matched points between two models to be compared, and registration is needed before comparison.
Most of the shapes before and after root canal treatment have nonrigid deformation. It is almost impossible to find the corresponding points in the shape of the operation, so the above methods are not suitable for the problems studied in this paper.
In the related literature of shape differences, Ovsjanikov et al. [38] presented a novel representation of maps between pairs of shapes, and an energy measuring distortions induced by a map was defined on the function space to measure the difference of shapes. Huang et al. [54] proposed a pipeline for constructing an area-based shape difference operator on point clouds and showed numerically that the results are robust and informative.
Suppose the spectrum extracted from the two root canal models is λ � λ i and μ � μ i . A naive difference evaluation method is to calculate the Euclidean distance between λ and μ as However, the influence of different eigenvalues on shape is different since the first dozens of eigenvalues determined the shape.
Kac [17] proposed the spectral distance weighted by i − 1 but did not consider the influence factors of the regional characteristics of a specific shape.
Based on the spectral theory, we know that the lowfrequency eigenvalue reflects the global geometry (that is, the high-frequency signal is cancelled in a flat place) and the high-frequency eigenvalue reflects the local geometry (resonance occurs at the end of the shape). is principle was used in document [19,21,55,56] to realize the interchange of gestures and details of multiple animated characters. erefore, for each eigenvalue, its contribution to the difference between two shapes is proportional to the size of the variation range of the shape details. Intuitively, the curvature is the amount by which a geometric object such as a surface deviates from being a flat plane. e intensity of detail changes on the mesh can be quantified by the Gaussian curvature.
e Gaussian curvature of the triangular mesh surface can be calculated by using the following formula: where x ∈ V and A(x)is the area of the 1-ring region centered on x. θ i is the angle of one of the adjacent triangles. Examples of the Gaussian curvature are shown in Figure 7. Accordingly, we propose a spectral difference evaluation algorithm based on Gaussian curvature statistical distribution. Suppose the mesh model of a root canal is M � (V, E, F), N(v) represent the number vertex in V. We select the first l eigenvalues in the spectrum to evaluate the shape difference. e spectral difference evaluation result d(M, N) can be calculated by using Algorithm 1.
rough the above algorithm, we can quantify the weight of the difference between two spectrums in the overall shape difference comprehensively. erefore, we can get better experimental results.

Experiment and Results
Traditional evaluation indicators for the effectiveness of systems include volume and area difference, as well as the transportation of the root canal. In this paper, the conicity and permeability of root canals are evaluated from the spectral distance between shapes. To access the superiority of the algorithm, the accuracy and the robustness of the algorithm are evaluated.

Programming.
is paper mainly deals with the numerical solution of the integral equation. e most crucial problem is to find the eigenvalues of a matrix, i.e., to solve large-scale sparse linear system. erefore, we choose the popular eigen algorithm library to calculate. It provides a series of built-in solver libraries and encapsulates some external solver libraries, including LLT and LDLT used in this paper. e algorithm proposed in this paper was implemented through C++ programming with the hardware and software shown in Table 1.

Comparison with LBO-Based Method.
e Laplace-Beltrami operator is an intrinsic invariant operator whose spectrum also reflects the geometric and topological characteristics of some shapes, so it is often used for shape difference comparison, such as [57]. Most of the models processed by previous works are mesh models with low genus or high regularity. However, the mesh generated from the point cloud of the root canal obtained from μCT is often of poor quality. erefore, for some spectral methods with low robustness, there may be more problems in the aspect of accuracy.
In this paper, the Steklov spectrum method and LBO spectrum method are used to calculate three groups of root canal models, respectively, and the results are shown in Figure 8. Our algorithm has an apparent distinction between the apparent uneven root canal shape, the rough root canal, and the different topology, while the performance of the LBO method is not ideal. erefore, we infer that our algorithm is more suitable for the root canal model with complex geometry and topology in these cases.

Comparison with Volume/Area Evaluation Method.
Our method is to compare the prepared root canal with the idealized root canal and then calculate the overall geometric structure and topology. Shape difference based on the extrinsic spectrum accurately perceives the concavity and convexity of the shape. So our method is volume aware, area aware, and slice aware in real sense.
To evaluate the comparative effect of SSD, we infer that the best root canal shape is the smallest distance from the ideal root canal shape in theory, and the best therapy result is obtained. In [32], spiral computerized tomography (SCT) was reported to be a useful tool in various in vivo and laboratory studies. It was concluded that with SCT, three-dimensional volume measurements are possible without sectioning specimens, thus avoiding loss of material.
Specifically, the comparative method is as follows: A total of 60 extracted first molars were allocated into three groups randomly. Each group was treated with one instrument, respectively. e effect of root canal preparation was evaluated by the literature [2,14,58], which was based on the percentage of void (POV). e lower the POV, the better the obturation effect. e aim of this comparison is to assess the correctness of the volume and area evaluation method and our method. Statistical analysis was performed with ANOVA (Scheffé tests).
e results of experiments summarized in Table 2 show that there is no significant difference among the three groups of results obtained by the volume-based discriminant method, and the third group which is considered to be the best one is not consistent with the ground truth evaluation result. erefore, the approach is considered to be unsuccessful. e results obtained by the area-based discriminant method did not express a significant difference between the second and the third groups, and the best (second) group is also inconsistent with the ground truth evaluation result either.
e evaluation results of the SSD method have achieved good efficacy in two respects. First, only the SSD method is consistent with the ground truth evaluation result. Second, among the three methods, only SSD is capable of distinguishing the significant difference of System A from the other. erefore, it can be concluded that the new method can evaluate the effect of root canal preparation comprehensively Figure 7: An example of calculating the Gauss curvature of the root canal model.  and accurately. Hence, the experiment achieves the desired purpose in the dataset.

Robustness Evaluation.
In Figure 9, it shows an example of the robustness of our method to different sampling density. We compare the eigenvalues of a series of the root canal model with different resolution. It is observed that some geometric details in some regions are lost during the downsampling and the spectrum remains almost unaffected. According to the theory of the Steklov spectrum, the first dozens of eigenvalues and eigenfunctions are robust to downsampling for a mesh with 10k to 100k vertices. Moreover, the number of vertices is generally sufficient for shape difference evaluation.

Conclusions and Prospects
In this study, we addressed the problem of evaluating the operation result of a root canal by introducing a spectrumbased procedure called the Steklov spectrum distance, which measures the shape difference between the original shape and its ideal shape with its Steklov spectrum. We find it remarkable that the use of SSD simplifies the problem of evaluating the intrinsic and extrinsic difference between shapes, allowing to significantly improve the accuracy of standard pipelines comparing with the traditional volumearea literature. SSD appears to be a valuable tool to assess the efficacy of various root canal preparation systems.
Within the limitations of this study, one is that the error of SSD was seen in all groups, and it even diverges in the calculation of certain root canals. To solve these problems, we think that it can be optimized in two ways: one is to design a more expressive distance norm according to the specific shape of root canal; the other is to improve the numerical solution method of partial differential equation, using better presolver and matrix block parallel algorithm to improve the accuracy and efficiency of the algorithm. Another shortage to be improved is that, limited by the μCT itself, this study can only be performed in vitro, rather than applied to direct clinical surgical evaluation. In the future, we will continue to consider the related research based on the data sources of high-speed acquisition of CT images in vivo.

Data Availability
e reason why the authors choose not to disclose the data is that the research group is applying for patent at the same time, which may involve some intellectual property issues.

Conflicts of Interest
e authors declare that there are no conflicts of interest regarding the publication of this paper.  Serial number of eigenvalues i Figure 9: ree root canal meshes represented with 10K (blue), 50K (red), and 100K (yellow) vertices share similar Steklov spectrums.