Numerically stable evaluation of closed-form expressions for eigenvalues of $3 \times 3$ matrices

TL;DR

This work tackles numerically unstable closed-form eigenvalue formulas for real $3\times3$ matrices by reexpressing the problem through four invariants: the trace $I_1$, the deviatoric invariants $J_2$ and $J_3$, and the discriminant $\Delta$. It develops stable, closed-form algorithms for evaluating these invariants (notably a backward-stable approach for $J_2$), analyzes forward error bounds via conditioning on the eigenbasis, and derives a backward-stable eigenvalue computation using a trigonometric formulation. Empirical benchmarks show that the $J_2$ algorithm achieves accuracy within the predicted stability bounds and can be about ten times faster than LAPACK on challenging cases, while $J_3$ and $\Delta$ are stable mainly for well-conditioned or orthogonal bases; ill-conditioned bases favor LAPACK. The work delivers an efficient, closed-form alternative for many practical matrices, with an open-source implementation (eig3x3) and clear directions for future theoretical guarantees and architecture-specific optimizations.

Abstract

Trigonometric formulas for eigenvalues of $3 \times 3$ matrices that build on Cardano's and Viète's work on algebraic solutions of the cubic are numerically unstable for matrices with repeated eigenvalues. This work presents numerically stable, closed-form evaluation of eigenvalues of real, diagonalizable $3 \times 3$ matrices via four invariants: the trace $I_1$, the deviatoric invariants $J_2$ and $J_3$, and the discriminant $Δ$. We analyze the conditioning of these invariants and derive tight forward error bounds. For $J_2$ we propose an algorithm and prove its accuracy. We benchmark all invariants and the resulting eigenvalue formulas, relating observed forward errors to the derived bounds. In particular, we show that, for the special case of matrices with a well-conditioned eigenbasis, the newly proposed algorithms have errors within the forward stability bounds. Performance benchmarks show that the proposed algorithm is approximately ten times faster than the highly optimized LAPACK library for a challenging test case, while maintaining comparable accuracy.

Figures (6)

• Figure 1: Typical approach for computing eigenvalues of $3 \times 3$ matrices via characteristic polynomial and its roots.
• Figure 2: Benchmark cases in this paper. The red squares represent the limiting path $\delta \to 0$.
• Figure 3: Numerical stability analysis for the invariant $J_2$.
• Figure 4: Numerical stability analysis for the invariant $J_3$.
• Figure 5: Numerical stability analysis for the discriminant $\Delta$.

Theorems & Definitions (24)

• Lemma 5.1
• proof
• Lemma 5.2
• proof
• Theorem 5.3
• proof
• Lemma 5.4
• proof
• Corollary 5.5
• proof