## Efficient linear system solvers for mesh processing (2005)

Venue: | In IMA Conference on the Mathematics of Surfaces |

Citations: | 39 - 5 self |

### BibTeX

@INPROCEEDINGS{Botsch05efficientlinear,

author = {Mario Botsch and David Bommes and Leif Kobbelt},

title = {Efficient linear system solvers for mesh processing},

booktitle = {In IMA Conference on the Mathematics of Surfaces},

year = {2005},

pages = {62--83},

publisher = {Springer}

}

### Years of Citing Articles

### OpenURL

### Abstract

Abstract. The use of polygonal mesh representations for freeform geometry enables the formulation of many important geometry processing tasks as the solution of one or several linear systems. As a consequence, the key ingredient for efficient algorithms is a fast procedure to solve linear systems. A large class of standard problems can further be shown to lead more specifically to sparse, symmetric, and positive definite systems, that allow for a numerically robust and efficient solution. In this paper we discuss and evaluate the use of sparse direct solvers for such kind of systems in geometry processing applications, since in our experiments they turned out to be superior even to highly optimized multigrid methods, but at the same time were considerably easier to use and implement. Although the methods we present are well known in the field of high performance computing, we observed that they are in practice surprisingly rarely applied to geometry processing problems. 1

### Citations

2105 |
Numerical recipes in C: The art of scientific computing
- PRESS, FLANNERY, et al.
- 1988
(Show Context)
Citation Context ...w, on page 9). 2.2 Iterative Solvers In contrast to dense direct solvers, iterative methods are able to exploit the sparsity of the matrix A. Since they additionally allow for a simple implementations=-=[33]-=-, iterative solvers are the de-facto standard method for solving sparse linear systems in the context of geometric problems. A detailed overview of iterative methods with precious implementation hints... |

1488 |
Practical optimization
- Gill, Murray, et al.
- 1986
(Show Context)
Citation Context ...inear ones, like, e.g., the (semi-)implicit integration of non-linear geometric flows by solving a linear equation in each time step [13] or the Levenberg-Marquardt method for non-linear optimization =-=[17]-=-. Similarly, continuous energy functionals E (f) = � e (f, x) dx are approximated Ω by quadratic forms E (X) = XT QX, such that their minimizer surfaces cansefficiently be derived by solving the linea... |

847 | A fast and high quality multilevel scheme for partitioning irregular graphs
- Karypis, Kumar
- 1998
(Show Context)
Citation Context ...dividually. If the adjacency graph is connected, a small subset S of nodes, whose elimination would separate the graph into two components of roughly equal size, is found by one of several heuristics =-=[21]-=-. This graph-partitioning results in a matrix consisting of two large diagonal blocks (two connected components) and |S| rows representing their connection (separator S). Recursively repeating this pr... |

831 |
W.(1985): Multigrid Methods and Applications
- Hackbusch
(Show Context)
Citation Context ... This fact is exploited by multigrid methods, that build a fine-to-coarse hierarchy {M = M0, M1, . . . , Mk} of the computation domain M and solve the linear system hierarchically from coarse to fine =-=[19, 8]-=-. After a few (pre-)smoothing iterations on the finest level M0 the high frequencies of the error are removed and the solver becomes inefficient. However, the remaining low frequency error e0 = x ∗ −x... |

768 |
Methods of conjugate gradients for solving linear systems
- Hestenes, Stiefel
- 1952
(Show Context)
Citation Context ...ods only converge for a restricted set of matrices, e.g., for diagonally dominant ones. Non-stationary iterative solvers are more powerful, and for spd matrices the method of conjugate gradients (CG) =-=[20, 18]-=- is suited best, since it provides guaranteed convergence with monotonically decreasing error. For a spd matrix A the solution of Ax = b is equivalent to the minimization of the quadratic form φ (x) :... |

568 | Applied Numerical Linear Algebra
- Demmel
- 1997
(Show Context)
Citation Context ...completeness, the general case of a non-symmetric indefinite system is outlined in Sect. 2.5. More elaborate surveys on how to efficiently solve general large linear systems can be found in the books =-=[10, 29]-=-. 2.1 Dense Direct Solvers Direct linear system solvers are based on a factorization of the matrix A into matrices of simpler structure, e.g., triangular, diagonal, or orthogonal matrices. This struct... |

554 | A signal processing approach to fair surface design
- Taubin
- 1995
(Show Context)
Citation Context ...h an spd matrix Q. A very popular source of spd systems is the discrete Laplace-Beltrami operator ∆S [30], which is closely related to frequencies of scalar fields defined on a two-manifold surface S =-=[41]-=-. This operator has various applications in surface smoothing [13, 41], surface parameterization [32, 12], variational surface modeling [25, 7, 38, 44], mesh morphing [3, 4, 43], and shape analysis [3... |

525 |
Computer Solution of Large Sparse Positive Definite Systems
- George, Liu
- 1981
(Show Context)
Citation Context ...Cholesky factor is a dense lower triangular matrix (cf. Fig. 2, top row). However, an analysis of the factorization process reveals that a band-limitation of the matrix A will be preserved. Following =-=[15]-=-, we define the bandwidth β (A) in terms of the bandwidth of its ith row β (A) := max 1≤i≤n {βi (A)} with βi (A) := i − min 1≤j≤i {j | Aij �= 0} . If the matrix A has bandwidth β (A) then so has its f... |

521 |
der Vorst. Templates for the solution of linear systems: Building blocks for iterative methods
- Barrett, Berry, et al.
- 1995
(Show Context)
Citation Context ...rs are the de-facto standard method for solving sparse linear systems in the context of geometric problems. A detailed overview of iterative methods with precious implementation hints can be found in =-=[5, 36]-=-. Iterative methods compute a converging sequence x (0) , x (1) , . . . , x (i) of approximations to the solution x ∗ of the linear system, i.e., limi→∞ x (i) = x ∗ . In practice, however, one has to ... |

456 | Implicit fairing of irregular meshes using diffusion and curvature flow
- Desbrun, Meyer, et al.
- 1999
(Show Context)
Citation Context ...dle non-linear problems is their decomposition into a sequence of linear ones, like, e.g., the (semi-)implicit integration of non-linear geometric flows by solving a linear equation in each time step =-=[13]-=- or the Levenberg-Marquardt method for non-linear optimization [17]. Similarly, continuous energy functionals E (f) = � e (f, x) dx are approximated Ω by quadratic forms E (X) = XT QX, such that their... |

451 | Stable fluids
- Stam
- 1999
(Show Context)
Citation Context ... modeling [25, 7, 38, 44], mesh morphing [3, 4, 43], and shape analysis [31]. Besides from surfaces, the standard Laplace operator is for instance also used in image editing [40] and fluid simulation =-=[39]-=-. Finally, all linear problems Ax = b that cannot be solved exactly and hence are approximated in the least squares sense by using the normal equations A T Ax = A T b also result in spd linear systems... |

369 | Discrete differential-geometry operators for Triangulated 2-Manifold. Available on WWW at http://www.multires.caltech.edu/pubs/pubs.htm
- Desbrun, Meyer, et al.
(Show Context)
Citation Context ...ns. Such systems frequently occur when minimizing energy functionals of the form E (X) = X T QX with an spd matrix Q. A very popular source of spd systems is the discrete Laplace-Beltrami operator ∆S =-=[30]-=-, which is closely related to frequencies of scalar fields defined on a two-manifold surface S [41]. This operator has various applications in surface smoothing [13, 41], surface parameterization [32,... |

343 | An introduction to the conjugate gradient method without the agonizing pain. http://www.cs.cmu.edu/∼quake-papers/painless-conjugate-gradient.pdf
- Shewchuk
- 1994
(Show Context)
Citation Context ...In order to cancel out the effect of A’s eigenvalues on the search directions pi, those are chosen to be A-conjugate, i.e., orthogonal w.r.t. the scalar product induced by A: p T j Api = 0 for i �= j =-=[37]-=-. The computation of and minimization along these optimal search directions can be performed efficiently and with a constant memory consumption. The complexity of each CG iteration is mainly determine... |

285 | Interactive multiresolution modeling on arbitrary meshes
- Kobbelt, Campagna, et al.
- 1998
(Show Context)
Citation Context ...encies of scalar fields defined on a two-manifold surface S [41]. This operator has various applications in surface smoothing [13, 41], surface parameterization [32, 12], variational surface modeling =-=[25, 7, 38, 44]-=-, mesh morphing [3, 4, 43], and shape analysis [31]. Besides from surfaces, the standard Laplace operator is for instance also used in image editing [40] and fluid simulation [39]. Finally, all linear... |

281 | Computing discrete minimal surfaces and their conjugates. Experimental Mathematics
- Pinkall, Polthier
- 1993
(Show Context)
Citation Context ...[30], which is closely related to frequencies of scalar fields defined on a two-manifold surface S [41]. This operator has various applications in surface smoothing [13, 41], surface parameterization =-=[32, 12]-=-, variational surface modeling [25, 7, 38, 44], mesh morphing [3, 4, 43], and shape analysis [31]. Besides from surfaces, the standard Laplace operator is for instance also used in image editing [40] ... |

265 | Least squares conformal maps for automatic texture atlas generations
- Lévy, Petitjean, et al.
- 2002
(Show Context)
Citation Context ... Finally, all linear problems Ax = b that cannot be solved exactly and hence are approximated in the least squares sense by using the normal equations A T Ax = A T b also result in spd linear systems =-=[26]-=-. This large but still incomplete list of applications involving spd systems legitimates focusing on this special class of problems. Another important point to be considered is whether the linear syst... |

225 |
Reducing the bandwidth of sparse symmetric matrices
- Cuthill, McKee
- 1969
(Show Context)
Citation Context ...rresponding to the non-zeros of A, i.e., two nodes i, j ∈ {1, . . . , n} are connected by an edge if and only if Aij �= 0.sThe standard method for envelope minimization is the Cuthill-McKee algorithm =-=[9]-=-, that picks a start node and renumbers all its neighbors by traversing the adjacency graph in a breadth-first manner, using a greedy selection in the order of increasing valence. It has further been ... |

203 | A supernodal approach to sparse partial pivoting
- Demmel, Eisenstat, et al.
- 1995
(Show Context)
Citation Context ... conservative envelope is typically used, such that pivoting within this structure is still possible. A highly efficient implementation of a sparse LU factorization is provided by the SuperLU library =-=[11]-=-. 3 Laplace Systems Most of the example applications shown in Sect. 4 require the solution of linear Laplacian systems, therefore we analyze these matrices and compare different solvers for their solu... |

174 | Intrinsic parameterizations of surface meshes
- Desbrun, Meyer, et al.
- 2002
(Show Context)
Citation Context ...[30], which is closely related to frequencies of scalar fields defined on a two-manifold surface S [41]. This operator has various applications in surface smoothing [13, 41], surface parameterization =-=[32, 12]-=-, variational surface modeling [25, 7, 38, 44], mesh morphing [3, 4, 43], and shape analysis [31]. Besides from surfaces, the standard Laplace operator is for instance also used in image editing [40] ... |

174 | Laplacian surface editing
- SORKINE, COHEN-OR, et al.
- 2004
(Show Context)
Citation Context ...encies of scalar fields defined on a two-manifold surface S [41]. This operator has various applications in surface smoothing [13, 41], surface parameterization [32, 12], variational surface modeling =-=[25, 7, 38, 44]-=-, mesh morphing [3, 4, 43], and shape analysis [31]. Besides from surfaces, the standard Laplace operator is for instance also used in image editing [40] and fluid simulation [39]. Finally, all linear... |

143 |
Modification of the minimum degree algorithm by multiple elimination
- Liu
- 1985
(Show Context)
Citation Context ...he additional edges eij generated in this socalled elimination graph correspond to the fill-in elements Lij �= 0 = Aij. In order to minimize fill-in the minimum degree algorithm (MD) and its variants =-=[16, 27]-=- remove the nodes with smallest valence first from the elimination graph, since this causes the least number of additional pairwise connections. Many efficiency optimizations of this basic method exis... |

131 |
The evolution of the minimum degree ordering algorithm
- George, Liu
- 1989
(Show Context)
Citation Context ...he additional edges eij generated in this socalled elimination graph correspond to the fill-in elements Lij �= 0 = Aij. In order to minimize fill-in the minimum degree algorithm (MD) and its variants =-=[16, 27]-=- remove the nodes with smallest valence first from the elimination graph, since this causes the least number of additional pairwise connections. Many efficiency optimizations of this basic method exis... |

126 | Multiresolution modeling: Survey & future opportunities, Eurographics State of The Art Report - Garland - 1999 |

100 | Discrete fairing. In
- Kobbelt
- 1996
(Show Context)
Citation Context ..., x) dx are approximated Ω by quadratic forms E (X) = XT QX, such that their minimizer surfaces cansefficiently be derived by solving the linear systems QX = B, assuming proper boundary constraints B =-=[22]-=-. These examples motivate why a large class of geometric problems comes down to the solution of one or several linear systems. As a consequence, high performance linear system solvers are of major imp... |

97 | An intuitive framework for real-time freeform modeling
- Botsch, Kobbelt
(Show Context)
Citation Context ...encies of scalar fields defined on a two-manifold surface S [41]. This operator has various applications in surface smoothing [13, 41], surface parameterization [32, 12], variational surface modeling =-=[25, 7, 38, 44]-=-, mesh morphing [3, 4, 43], and shape analysis [31]. Besides from surfaces, the standard Laplace operator is for instance also used in image editing [40] and fluid simulation [39]. Finally, all linear... |

95 | Poisson matting
- Sun, Jia, et al.
- 2004
(Show Context)
Citation Context ..., 12], variational surface modeling [25, 7, 38, 44], mesh morphing [3, 4, 43], and shape analysis [31]. Besides from surfaces, the standard Laplace operator is for instance also used in image editing =-=[40]-=- and fluid simulation [39]. Finally, all linear problems Ax = b that cannot be solved exactly and hence are approximated in the least squares sense by using the normal equations A T Ax = A T b also re... |

86 | A general framework for mesh decimation
- Kobbelt, Campagna, et al.
- 1998
(Show Context)
Citation Context ...ur case the discrete computational domain M is an irregular triangle mesh instead of a regular 2D or 3D grid, the coarsening operator for building the hierarchy is based on mesh decimation techniques =-=[24, 14]-=-. The shape of the resulting triangles is important for numerical robustness, and the edge lengths on the different levels should mimic the case of regular grids. Therefore the decimation usually remo... |

75 |
Differential coordinates for local mesh morphing and deformation
- ALEXA
(Show Context)
Citation Context ... on a two-manifold surface S [41]. This operator has various applications in surface smoothing [13, 41], surface parameterization [32, 12], variational surface modeling [25, 7, 38, 44], mesh morphing =-=[3, 4, 43]-=-, and shape analysis [31]. Besides from surfaces, the standard Laplace operator is for instance also used in image editing [40] and fluid simulation [39]. Finally, all linear problems Ax = b that cann... |

64 |
Fair Morse Functions for Extracting the Topological Structure of a Surface Mesh
- Ni, Garland, et al.
- 2004
(Show Context)
Citation Context ...1]. This operator has various applications in surface smoothing [13, 41], surface parameterization [32, 12], variational surface modeling [25, 7, 38, 44], mesh morphing [3, 4, 43], and shape analysis =-=[31]-=-. Besides from surfaces, the standard Laplace operator is for instance also used in image editing [40] and fluid simulation [39]. Finally, all linear problems Ax = b that cannot be solved exactly and ... |

60 |
Computer solution of large linear systems
- Meurant
- 1999
(Show Context)
Citation Context ...completeness, the general case of a non-symmetric indefinite system is outlined in Sect. 2.5. More elaborate surveys on how to efficiently solve general large linear systems can be found in the books =-=[10, 29]-=-. 2.1 Dense Direct Solvers Direct linear system solvers are based on a factorization of the matrix A into matrices of simpler structure, e.g., triangular, diagonal, or orthogonal matrices. This struct... |

49 |
Comparative analysis of the Cuthill-Mckee and the reverse Cuthill-Mckee ordering algorithms for sparse matrices
- 25Liu, Sherman
- 1976
(Show Context)
Citation Context ...ks a start node and renumbers all its neighbors by traversing the adjacency graph in a breadth-first manner, using a greedy selection in the order of increasing valence. It has further been proven in =-=[28]-=- that reverting this permutation leads to better re-orderings, such that usually the reverse Cuthill-McKee method (RCMK) is employed. The result P T AP of this matrix re-ordering is depicted in the se... |

42 |
The cascadic multigrid method for elliptic problems
- Bornemann, Deuflhard
- 1996
(Show Context)
Citation Context ...and completely avoid V-cycles. In the latter case the number of iterations mi on level i must not be constant, but instead has to be chosen as mi = m γ i to decrease exponentially from coarse to fine =-=[6]-=-. Besides the easier implementation, the advantage of this cascading multigrid method is that once a level is computed, it is not involved in further computations and can be discarded. A comparison of... |

36 | Schröder,Multilevel Solvers for Unstructured Surface Meshes
- Aksoylu, Khodakovsky, et al.
(Show Context)
Citation Context ...o the next coarser one Mi+1 should additionally be restricted to remove a maximally independent set of vertices, i.e., no two removed vertices vj, vl ∈ Mi \ Mi+1 are connected by an edge ejl ∈ Mi. In =-=[2]-=- some more efficient alternatives to this standard DobkinKirkpatric hierarchy are discussed. In order to achieve higher performance, we do not change the simple way the hierarchy is constructed, but i... |

30 |
Local control for mesh morphing
- Alexa
- 2001
(Show Context)
Citation Context ... on a two-manifold surface S [41]. This operator has various applications in surface smoothing [13, 41], surface parameterization [32, 12], variational surface modeling [25, 7, 38, 44], mesh morphing =-=[3, 4, 43]-=-, and shape analysis [31]. Besides from surfaces, the standard Laplace operator is for instance also used in image editing [40] and fluid simulation [39]. Finally, all linear problems Ax = b that cann... |

28 |
der Vorst. Iterative solution of linear systems in the 20th century
- Saad, van
(Show Context)
Citation Context ...rs are the de-facto standard method for solving sparse linear systems in the context of geometric problems. A detailed overview of iterative methods with precious implementation hints can be found in =-=[5, 36]-=-. Iterative methods compute a converging sequence x (0) , x (1) , . . . , x (i) of approximations to the solution x ∗ of the linear system, i.e., limi→∞ x (i) = x ∗ . In practice, however, one has to ... |

27 |
Taucs: A library of sparse linear solvers. http://www.tau.ac.il/ stoledo/taucs
- Toledo, Chen, et al.
- 2003
(Show Context)
Citation Context ... chosen quite easily. For more details and implementation notes the reader is referred to the book of George and Liu [15]; a highly efficient implementation is publicly available in the TAUCS library =-=[42]-=-. 2.5 Non-Symmetric Indefinite Systems When the assumptions about the symmetry and positive definiteness of the matrix A are not satisfied, optimal methods like the Cholesky factorization or conjugate... |

26 |
Hierarchical Least Squares Conformal Map
- RAY, LEVY
- 2003
(Show Context)
Citation Context ...fic instance (values of A). Nevertheless, if iterative solvers are to be used, multigrid methods are the only way to achieve acceptable computing times for solving large systems, as has been shown in =-=[25, 34, 2]-=-. 2.4 Sparse Direct Solvers The use of direct solvers for large sparse linear systems is often neglected, since naïve direct methods have complexity O(n 3 ), as described above. The problem is that ev... |

24 |
Applied Numerical Linear Algebra (SIAM
- Demmel
- 1997
(Show Context)
Citation Context ...completeness, the general case of a non-symmetric indefinite system is outlined in Sect. 2.5. More elaborate surveys on how to efficiently solve general large linear systems can be found in the books =-=[10, 29]-=-. 2.1 Dense Direct Solvers Direct linear system solvers are based on a factorization of the matrix A into matrices of simpler structure, e.g., triangular, diagonal, or orthogonal matrices. This struct... |

22 |
Poisson shape interpolation
- XU, ZHANG, et al.
- 2005
(Show Context)
Citation Context ... on a two-manifold surface S [41]. This operator has various applications in surface smoothing [13, 41], surface parameterization [32, 12], variational surface modeling [25, 7, 38, 44], mesh morphing =-=[3, 4, 43]-=-, and shape analysis [31]. Besides from surfaces, the standard Laplace operator is for instance also used in image editing [40] and fluid simulation [39]. Finally, all linear problems Ax = b that cann... |

16 |
Hujun Bao, Baining Guo, Heung-Yeung Shum, Mesh editing with poisson-based gradient field manipulation
- Yu, Zhou, et al.
(Show Context)
Citation Context |

11 |
Freeform shape representations for efficient geometry processing, Shape Modelling International
- Kobbelt, Botsch
- 2003
(Show Context)
Citation Context ...to geometry processing problems. 1 Introduction In the field of geometry processing suitable data structures that enable the implementation of efficient algorithms are getting more and more important =-=[23]-=-, especially since the complexity of the geometric models to be processed is growing much faster than the steadily increasing computational power and available memory of today’s PC systems. Typical ex... |

4 |
A Multigrid Tutorial. SIAM, 2nd edition
- Briggs, Henson, et al.
- 2000
(Show Context)
Citation Context ... This fact is exploited by multigrid methods, that build a fine-to-coarse hierarchy {M = M0, M1, . . . , Mk} of the computation domain M and solve the linear system hierarchically from coarse to fine =-=[19, 8]-=-. After a few (pre-)smoothing iterations on the finest level M0 the high frequencies of the error are removed and the solver becomes inefficient. However, the remaining low frequency error e0 = x ∗ −x... |

2 |
Gmm++: a generic template matrix C++ library. http://www-gmm.insatoulouse.fr/getfem/gmm_intro
- RENARD, POMMIER
- 2005
(Show Context)
Citation Context ...olvers for Laplacian as well as for bi-Laplacian systems. All timings we report in this and the next section were taken on a 3.0GHz Pentium4 running Linux. The iterative solver from the gmm++ library =-=[35]-=- is based on the conjugate gradients method and uses an incomplete LDL T factorization as preconditioner. Our cascading multigrid solver performs preconditioned conjugate gradient iterations on each h... |