Table of Contents
Fetching ...

A Comparison of Sparse Solvers for Severely Ill-Conditioned Linear Systems in Geophysical Marker-In-Cell Simulations

Marcel Ferrari

Abstract

Solving sparse linear systems is a critical challenge in many scientific and engineering fields, particularly when these systems are severely ill-conditioned. This work aims to provide a comprehensive comparison of various solvers designed for such problems, offering valuable insights and guidance for domain scientists and researchers. We develop the tools required to accurately evaluate the performance and correctness of 16 solvers from 11 state-of-the-art numerical libraries, focusing on their effectiveness in handling ill-conditioned matrices. The solvers were tested on linear systems arising from a coupled hydro-mechanical marker-in-cell geophysical simulation. To address the challenge of computing accurate error bounds on the solution, we introduce the Projected Adam method, a novel algorithm to efficiently compute the condition number of a matrix without relying on eigenvalues or singular values. Our benchmark results demonstrate that Intel oneAPI MKL PARDISO, UMFPACK, and MUMPS are the most reliable solvers for the tested scenarios. This work serves as a resource for selecting appropriate solvers, understanding the impact of condition numbers, and improving the robustness of numerical solutions in practical applications.

A Comparison of Sparse Solvers for Severely Ill-Conditioned Linear Systems in Geophysical Marker-In-Cell Simulations

Abstract

Solving sparse linear systems is a critical challenge in many scientific and engineering fields, particularly when these systems are severely ill-conditioned. This work aims to provide a comprehensive comparison of various solvers designed for such problems, offering valuable insights and guidance for domain scientists and researchers. We develop the tools required to accurately evaluate the performance and correctness of 16 solvers from 11 state-of-the-art numerical libraries, focusing on their effectiveness in handling ill-conditioned matrices. The solvers were tested on linear systems arising from a coupled hydro-mechanical marker-in-cell geophysical simulation. To address the challenge of computing accurate error bounds on the solution, we introduce the Projected Adam method, a novel algorithm to efficiently compute the condition number of a matrix without relying on eigenvalues or singular values. Our benchmark results demonstrate that Intel oneAPI MKL PARDISO, UMFPACK, and MUMPS are the most reliable solvers for the tested scenarios. This work serves as a resource for selecting appropriate solvers, understanding the impact of condition numbers, and improving the robustness of numerical solutions in practical applications.
Paper Structure (50 sections, 5 theorems, 69 equations, 3 figures, 2 tables, 2 algorithms)

This paper contains 50 sections, 5 theorems, 69 equations, 3 figures, 2 tables, 2 algorithms.

Key Result

Theorem 1

Given a linear system and an approximate solution such that: Then the relative error $\|\mathbf{\hat{x}} - \mathbf{x}\| / \|\mathbf{x}\|$ of the approximate solution $\mathbf{\hat{x}}$ to the true solution $\mathbf{x}$ is bounded by: where $\kappa(\mathbf{A})$ is the condition number of the matrix $\mathbf{A}$ as defined in eq. eq:kappa.

Figures (3)

  • Figure 1: Visualisation of the sparsity structure of the matrix solved during the benchmark iterations. The image on the left shows a section of the matrix detailing the overall tri-banded structure, while the image on the right provides a detailed view of the main diagonal band.
  • Figure 4: Performance and accuracy results of evaluated sparse direct solvers. Each row presents the following data: optimal core count (* indicates no multithreading support or compiled without multithreading), total time to compute 30 benchmark iterations (including both symbolic and numerical factorization for each iteration), median iteration time, worst-case error for a benchmark linear system $\mathbf{Ax} = \mathbf{B}$, and tight error upper bound for the initial conditions system $\mathbf{A_0x_0} = \mathbf{B_0}$. Bootstrapped $95\%$ relative CIs are provided for total time and median iteration time when exceeding $1\%$. The values in bold represent the best recorded result.
  • Figure 5: Performance results with reused symbolic factorization. This table lists only the solvers that showed more than a $10\%$ performance improvement by reusing the symbolic factorization from the first benchmark iteration. Column definitions are the same as in table \ref{['tab:direct_full_benchmark']}.

Theorems & Definitions (5)

  • Theorem 1: Bounds on the relative error $\|\mathbf{\hat{x}} - \mathbf{x}\| / \|\mathbf{x}\|$
  • Theorem 2: Bounds on the relative error $\|\mathbf{\hat{x}} - \mathbf{x}\| / \| \mathbf{\hat{x}}\|$
  • Theorem 3: Tight upper bound of the error $\|\mathbf{\hat{x}} - \mathbf{x}\| / \| \mathbf{x}\|$
  • Theorem 4: Numerical singularity of a matrix
  • Theorem 5: Equivalent definitions of the $\mathcal{L}_2$ matrix norm