A Modified Direct Eigenvalue Swapping Algorithm with Guaranteed Stability

1 Introduction

• The analysis restricts attention to swaps between two 2×2 diagonal blocks, the stated most challenging case for numerical stability.More general adjacent swaps involve 2×2, 3×3, or 4×4 matrices; treats 2×2 and 3×3 cases.Source: §1 · p. 1
• Earlier iterative swapping methods may fail on some ill-conditioned problems despite typically completing the operation in a few iterations.The introduction also situates direct swapping as an alternative that constructs the swapping orthogonal matrix from a QR factorization.Source: §1 · p. 2
• Direct swapping solves a Sylvester equation, then forms an orthogonal matrix from a basis for the solution-derived subspace.The method is used in LAPACK’s DLAEXC; its stability can be unsatisfactory on some ill-conditioned problems.Source: §1 · p. 2
• A complete rounding-error analysis remains absent, and Bai and Demmel’s analysis assumes exact QR even though its bound appears pessimistic in experiments.The paper identifies QR-factorization error as an unresolved part of the stability question.Source: §1 · p. 2
• The modified algorithm replaces Householder QR with a pivoted Givens-based strategy and is shown stable under mild assumptions using a mixed forward–backward error measure.The paper also gives an example of greater stability than the approaches in [1, 2].Source: §1 · pp. 2–3

2 Preliminaries

• The section gathers notation and lemmas for the subsequent analysis, including results on GECP and Givens rotations.It also includes a componentwise perturbation result for the SVD of a 2×2 matrix.Source: §2 · p. 3, §2.2 · pp. 4–5

2.1 Notation

• The analysis uses unit roundoff u and the standard floating-point model fl(x op y) = (x op y)(1 + δ), with |δ| ≤ u.The paper assumes u is small enough that γ_k ≪ 1 for constants k used in the analysis.Source: §2.1 · p. 3
• The notation defines vectorization and Givens rotations, along with matrix norms, condition numbers, and singular values used in the bounds.Applying the defined rotation to v zeros its second entry.Source: §2.1 · pp. 3–4

2.2 Auxiliary lemmas

• The auxiliary results include GECP componentwise backward stability and consequences for the computed solution and residual.For complete pivoting, the stated growth factor is ρ_4 = 4.Source: §2.2 · p. 4
• The lemmas also bound perturbations to computed Givens rotations and provide a componentwise backward-error result for them.The Givens-rotation perturbation result is attributed to Sai.Source: §2.2 · p. 4
• A separate lemma gives a componentwise perturbation result for the SVD of a 2×2 matrix under a perturbation-size condition.It assumes c + κ^-1(X) < 1 and ∥ΔX∥ ≤ c∥X∥.Source: §2.2 · p. 5

3 Pivoted Givens QR

• The factorization targets a matrix built from X and the scaled identity, with a nonsingular Z that need not be upper triangular.The QR strategy is introduced as a replacement for the Householder QR used in [1].Source: §3 · p. 5
• When the largest entry of |X| is in the first column, the method accumulates Q using Givens rotations instead of Householder QR.The passage gives a separate accumulation procedure when the largest entry is in the second column.Source: §3 · p. 5
• The modified direct swapping algorithm uses pivoted Givens QR, accumulating Q as the product of four Givens rotations.The matrix being factored is formed from X and the scaled identity; Z need only be nonsingular, not upper triangular.Source: §3 · p. 5

4 Main result

• The paper measures distance to the set of exact swapping matrices as the minimum norm distance from a matrix M to any Q in P(A).P(A) is the compact set of orthogonal matrices that swap A.Source: §4 · p. 5
• A computed matrix that is only approximately orthogonal cannot, in principle, be an exact solution to the reordering problem.This motivates an error measure that incorporates both perturbation and distance to exact solutions.Source: §4 · pp. 5–6
• Algorithm 2 solves a linear system by GECP, factors a matrix built from X and γI2 with pivoted Givens QR, then tests and applies the swap.An accepted swap is followed by reduction back to real Schur form.Source: §4 · p. 6
• The mixed forward–backward error is the smallest spectral-norm perturbation to A that places the approximate Q within distance θ of exact solutions for the perturbed problem.The error is monotonically decreasing in θ.Source: §4 · p. 6
• Algorithm 2 is stable in the norm of ˇA2,1 when A satisfies Theorem 4.1's assumption and ξ(A) is bounded away from u^-1.The corollary presents this as a stability guarantee under those conditions.Source: §4 · p. 7
• The proof bounds are aimed at O(u) error, but their constants may be too large for the practical algorithm's threshold and are not claimed to be sharp.In experiments, the actual norm of ˇA2,1 is reported as much smaller than the provided bound.Source: §4 · p. 7

5 Numerical experiments

• The comparison uses double-precision DLAEXC, MDS with pivoted Givens QR, and RDS, which refines when its stability threshold is exceeded.The algorithms use the same stability tests as DLAEXC, and the reference implementations are tested in their native languages.Source: §5 · pp. 7–8

5.1 Test matrices

• The three test matrices include a difficult case with ξ(A) ≈ 2.4 × 10^8 and |a2,1/a1,2| ≈ 3.2 × 10^12, a similarity transform of that case, and a case with separated eigenvalues but a small imaginary part.For the third matrix, ξ(A) ≈ 15.2 and the imaginary part of λ1 is on the order of 10^-7, much smaller than ∥A∥.Source: §5.1 · p. 8

5.2 Numerical results

• On test matrix 1, DLAEXC fails its stability test while MDS stays below the threshold; MDS also has smaller relative eigenvalue errors on matrices 1 and 2.The analysis attributes the contrasting threshold performance on matrix 1 to the QR factorization, with pivoted Givens QR more stable than Householder QR.Source: §5.2 · p. 9, Table 2 · p. 9
• On test matrix 3, MDS and DLAEXC perform comparably for both the swapped off-diagonal block norm and relative eigenvalue errors.For both algorithms, the relative error of λ1 is significantly larger than unit roundoff on all test matrices, even when the swapped off-diagonal block norm is near u∥A∥max.Source: §5.2 · p. 9
• RDS fails on its first attempt for test matrices 1 and 2, and one refinement step lowers the swapped off-diagonal block norm below threshold but appears less effective for eigenvalue stability.The reported final relative error indicates severe loss of precision for ˇλ2.Source: §5.2 · pp. 9–10, Table 3 · p. 9

• For the QR step, the analysis compares the computed factorization of the perturbed matrix with an exact factorization and accounts for the resulting factor's approximate orthogonality.The pivoted Givens QR analysis constructs a nearby orthogonal matrix and a perturbation to the computed X.Source: §6 · p. 10
• The analysis models the computed linear-system solution through residuals and perturbations, then bounds the smallest perturbation to the matrix data consistent with that residual.A Moore–Penrose pseudoinverse characterizes the minimum-norm perturbation.Source: §6 · pp. 10–11
• The bounds combine normwise and componentwise analyses of the residual, with the latter depending on the relative sizes of off-diagonal entries in the two diagonal blocks.The argument considers the case κ̂ > 4 after settling the other case, and relates singular vectors and values of the perturbed and computed solutions.Source: §6 · pp. 11–12
• Under Theorem 4.1's assumptions, the analysis concludes that the perturbation ∥∆A∥ is bounded.The supplied passage states the conclusion but does not include the bound's expression.Source: §6 · p. 12

7 Contribution to error from solving linear system

• Componentwise bounds for S yield the section’s main result, with the desired bound established directly when |â1,2| ≤ |â2,1| under the stated conditions.The analysis identifies two remaining cases when |â1,2| > |â2,1|.Source: §7 · pp. 12, 14
• The preliminary propositions characterize τ, bound the smallest entries of Û and V̂, and relate off-diagonal entries before the remaining cases are analyzed.The analysis uses permutation invariance to impose conditions on A without changing ϕ̂(A).Source: §7 · pp. 12–14

7.1 Case 1

• In Case 1, step-by-step GECP pivot bounds show that the bound in (7.12) holds under assumptions (7.10)–(7.11).The analysis tracks entry growth through the pivots and combines componentwise bounds to bound ϕ̂(A).Source: §7.1 · pp. 15, 18

7.2 Case 2

• Proposition 7.7 bounds entries in |T| for Case 2 under assumptions (7.29)–(7.30) and the section’s stated conditions.The detailed analysis considers |a4,3| > |a1,2|; its bounds include an extra η2 + 1 multiplier for |a4,3|, and η2 ≥ 7 in this setting.Source: §7.2 · pp. 18–19
• When |a1,2| ≥ |a4,3|, the bound (7.31) still holds.The passage states that this follows by applying the preceding analysis.Source: §7.2 · p. 21
• When |a1,2| = |a4,3|, choosing −a4,3 and a1,2 as the first and second pivots, or vice versa, bounds ˆϕ(A) by O((1 + ϵ_r)u||A||).In this case, s2,1, s3,1, and s4,1 are also bounded by O((1 + ϵ_r)u||A||).Source: §7.2 · p. 21

• The analysis compares finite-precision Givens QR of the computed matrix with exact Givens QR of a nearby matrix, constructing a nearby orthogonal factor.It assumes the largest-magnitude entry of X is in its first column and verifies the constructed perturbation and orthogonal factor satisfy the required condition.Source: §8 · pp. 21–23
• The QR-factorization contribution to the error bound is assembled from bounds on two terms, yielding the stated result for κ̂ > 4.The analysis bounds contributions associated with the perturbed linear-system solution and the difference between the orthogonal factors.Source: §8 · p. 23
• The first two Givens updates must use the specified ordering to ensure a perturbation ∆X exists that satisfies the required condition.Different orderings in another column still give an O(u) bound for the difference between the exact and computed orthogonal factors, but may change its constant.Source: §8 · pp. 21, 23

9 Conclusion

• The proposed pivoted Givens QR replaces Householder QR in DLAEXC for swapping adjacent 2×2 diagonal blocks.The analysis measures error using ∥ˇA2,1∥.Source: §9 · p. 23
• The authors report empirical evidence of an advantage over existing methods.The paper says the modified algorithm can be used confidently in practice with respect to error measured by ∥ˇA2,1∥.Source: §9 · p. 23
• Relative eigenvalue errors may remain large even when ∥ˇA2,1∥max is on the order of u∥A∥max.The authors suggest eigenvalue fragility may explain this and identify small relative eigenvalue errors as a challenging goal.Source: §9 · p. 23
• Stable reordering for the generalized eigenvalue problem remains an open challenge, greater than for the standard problem.The authors identify a provably stable generalized-problem algorithm as future work.Source: §9 · p. 23

Declaration of AI usage

• DeepSeek-V4-Flash was used exclusively for language refinement.The passage also identifies proofreading as a use.Source: Declaration of AI usage · p. 24
• DeepSeek-V4-Flash was used for proofreading.Source: Declaration of AI usage · p. 24
• The authors are responsible for the scientific contributions, including algorithm development and error analysis.Source: Declaration of AI usage · p. 24

A.1 Proof of Lemma 2.3

• The computed Givens parameters satisfy |ĉ − c| ≤ γ2k+4|c| and |ŝ − s| ≤ γ2k+4|s|.Source: §A.1 · p. 25

A.3 Proof of Lemma 2.5

• The invoked singular-value perturbation result assumes σ̃r ≥ α + δ and σr+1 ≤ α, with δ > 0.The excerpt defines ϵ using ∥(∆X)Ṽ1∥ and ∥(∆X)⊤Ũ1∥.Source: §A.3 · p. 27

A.4 Proof of Proposition 7.1

• The permutation argument yields Ŷ = P X̂ and R̂Y = P R̂X for the computed solutions and residuals.Source: §A.4 · p. 28

A.6 Proof of Proposition 7.3

• The proof establishes (ĉ²â1,2 + ŝ²â2,1)(ĉ²â2,1 + ŝ²â1,2) < 0 and the stated alternative bounds on |ŝ| and |ĉ|.The bounds are |ŝ| > √η1 |ĉ| or |ĉ| > √η1 |ŝ|.Source: §A.6 · p. 29

A.7 Proof of Proposition 7.5

• The proof uses equality of the traces of each original diagonal block and its corresponding swapped block.Source: §A.7 · p. 29

A.8 Proof of Proposition 7.6

• Under the stated assumptions, 3|a2,1| ≤ |a1,2| follows from Proposition 7.4 (2).Source: §A.8 · p. 30

A.9 Proof of Proposition 7.7

• The relative bound for |a4,3| in this proof uses assumption (7.30).Source: §A.9 · p. 31