Source-linked AI summary
Solving Linear Systems in $\widetilde{O}(mn \log \frac{\kappa}{\epsilon})$ Bit Operations
Jonathan A. Kelner
TL;DR
Sparse linear systems pose a gap between exact-arithmetic iteration bounds and finite-precision bit complexity. The paper uses perturbation and lifting to construct a deterministic solver. For sparse, polynomially conditioned systems, it achieves near-quadratic time at inverse-polynomial accuracy.
Problem
Exact-arithmetic methods for sparse systems do not directly yield comparable bit-complexity guarantees, because intermediate values can grow and rounding can disrupt convergence.
Source: §1 · pp. 1–2Method
The algorithm perturbs the system so Dixon lifting uses trivial modular inversion and sparse iterations, then recovers a rational solution.
Source: §1 · p. 3, §2.1 · p. 5Results
eO(n^2) time solves sparse systems with polynomial condition number to inverse-polynomial accuracy.
Source: §1 · p. 2, §1.1 · p. 3Takeaways & Limitations
The result closes much of the gap between exact and finite-precision arithmetic for sparse linear systems.
Source: §1 · p. 2Ω(D) fixed-point output bits may be needed for residual below 1/2 even when the condition number is 1 and an input exponent encodes 2^-D.
Source: §4.4 · p. 15
Abstract
from arXiv · showhide
We give a deterministic algorithm that solves a nonsingular linear system $Ax=b$, where $A\in\mathbb{R}^{n\times n}$ has $m$ nonzero entries and condition number $\kappa$, to any relative residual tolerance $0<\epsilon\le1/2$ using $\widetilde{O(}mn\log(\kappa/\epsilon))$ bit operations for inputs with logarithmically many bits per entry. For sparse, polynomially conditioned systems with $m=\widetilde{O}(n)$, this gives an $\widetilde{O}(n^2)$ algorithm for any inverse-polynomial accuracy, improving on the algorithm of Peng and Vempala, as sharpened by Nie, whose running time in this regime is approximately $O(n^{2.2707})$ with the best current matrix multiplication exponent, and largely closing a gap between the idealized performance of the conjugate gradient method in exact arithmetic and the running time achievable with finite-precision computation that has persisted for over 70 years. The algorithm is surprisingly simple. For integer inputs, we apply Dixon's lifting algorithm to the perturbed system $(A+I/R)x=b$ for a suitable integer $R$. After scaling, the matrix of this system is $RA+I\equiv I\pmod R$, so its modular inverse is trivial, and each lifting step needs only a sparse matrix-vector product with $A$ on $O(\log R)$-bit numbers. Fast rational reconstruction then recovers the exact solution of the perturbed system, which is an $\epsilon$-accurate solution of the original one. Normalization and rounding extend the result to fixed-point and floating-point inputs, with floating-point outputs represented using short integer significands and a common encoded exponent. A computable certificate removes the need for prior knowledge of $\kappa$.
1 Introduction
Sparse linear systems expose a gap between exact-arithmetic iteration bounds and finite-precision bit complexity. The paper presents a deterministic solver with near-quadratic time for sparse, polynomially conditioned systems at inverse-polynomial accuracy.
1 Introduction
- Exact-arithmetic conjugate gradient can solve sparse systems in O(mn) operations, but finite-precision bit complexity can incur an extra factor of n.Intermediate exact rational values can grow linearly in bit length, while rounding can destroy the n-step termination guarantee.Source: §1 · p. 2
1.1 Main Result
- eO(n^2) time suffices for sparse systems with polynomial condition number and inverse-polynomial accuracy.More generally, polylogarithmic encodings, polynomial conditioning, and inverse-polynomial accuracy give eO(mn) time.Source: §1 · p. 2, §1.1 · p. 3
1 Introduction
- The method perturbs the system so modular inversion is trivial and lifting iterations use sparse matrix-vector products.It adapts Dixon’s integer-system algorithm and combines the iterations with faster rational recovery.Source: §1 · p. 3
1.1 Main Result
- The floating-point algorithm returns integer significands with a common exponent and requires no supplied condition-number bound.The returned significands and exponent have bit lengths bounded in terms of the theorem’s parameters.Source: §1.1 · p. 3
2 Background on Dixon’s Algorithm and Rational Reconstruction
The section reviews Dixon’s modular method for exact integer systems and rational reconstruction, then identifies costly steps in a straightforward implementation and the faster reconstruction used here.
2.1 Dixon’s Algorithm for Exactly Solving Systems with Integer Entries
- Dixon solves Mx=b by inverting M modulo R, successively lifting the solution modulo R^L, then reconstructing bounded rational coordinates from the final residue.The reconstruction is unique when the modulus is sufficiently large relative to the coordinate bounds.Source: §2.1 · pp. 4–5
2.2 Running Time of Dixon’s Algorithm
- Straightforward Dixon implementations can cost Õ(n^3): modular inversion costs Õ(n^ω), while a dense inverse can make each of Õ(n) lifting rounds quadratic.The entries used during lifting require only logarithmic precision under the stated polynomial-bounded input assumptions.Source: §2.2 · p. 6
- Replacing classical Euclidean reconstruction with a faster algorithm reduces Step 3 to Õ(n^2), while the algorithm also targets faster implementations of the first two steps.The faster reconstruction takes Õ(L) time per coordinate.Source: §2.2 · p. 6
2.3 Faster Reconstruction
- Scalar rational reconstruction recovers each exact solution coordinate from its modular residue when the theorem’s size bounds ensure uniqueness.Dixon’s method first computes a sufficiently precise p-adic approximation, then reconstructs the rational solution.Source: §2 · p. 4, §2.1 · p. 5, §2.3 · p. 7
3 Approximately Solving Systems with Integer Entries
For integer inputs, the paper gives a deterministic approximate solver whose bit complexity scales with matrix sparsity, dimension, input bit length, condition number, and requested accuracy. It applies Dixon lifting to a perturbed system, then bounds the residual and extends the approach to rounded inputs.
3 Approximately Solving Systems with Integer Entries
- For integer inputs, the deterministic solver returns rational x with relative residual at most ε using O~(mnτ) bit operations and O~(n²τ) space.Here τ = max(B, log κ, log(1/ε)); output coordinates have O~(nτ) bits.Source: §3 · p. 8
- Setting M = I + RA makes M ≡ I (mod R), so Dixon lifting avoids modular inversion and uses sparse products with A.The perturbation yields the perturbed system; lifting intermediates use O(log R)-bit numbers.Source: §3 · p. 8
3.1 Dixon’s Algorithm for 𝑴≡𝑰(mod 𝑹)
- For R = 2^s ≥ 2nT_A, SolvePerturbed returns (I + RA)^−1b in O~((mn + n²)s) bit operations, without requiring A to be nonsingular.The routine also uses O~((m + n²)s) bits of storage and reconstructs rational coordinates.Source: §3.1 · p. 10
3.2 Bounding the Approximation Error
- An exact residual test certifies the original system; doubling trial precision removes prior knowledge of κ while retaining the stated work and space bounds.The search begins at a precision based on input magnitudes and requested accuracy, and terminates once the tested precision is sufficient.Source: §3.2 · p. 11
3.3 Input and Output Rounding
- For nonsingular A0, the rounding certificate guarantees relative residual at most ε/2 when S ≥ 64κ²(A0)n³/ε.The construction rounds scaled inputs and the output; the sufficient scale bound follows from the certificate’s error analysis.Source: §3.3 · pp. 11–12
4 Floating-Point Inputs and Proof of the Main Theorem
The section converts floating-point inputs to a normalized integer system, solves it at a chosen precision, and maps the result back to the original scale. A certificate enables precision search without prior knowledge of the condition number, and the analysis establishes the residual guarantee and bit complexity.
4 Floating-Point Inputs and Proof of the Main Theorem
- The floating-point-input procedure normalizes scales, rounds entries to short integers, and converts the working solution back to the original scale.It first selects rounding precision from a supplied condition-number bound, then replaces that choice with certified precision search.Source: §4 · p. 12
4.1 Normalization and conversion to integers
- Normalization preserves the condition number, bounds the normalized input magnitudes, and transfers a relative-residual guarantee unchanged to the original system.The conversion rounds normalized entries to integers; subsequent lifting uses those exact integers.Source: §4.1 · p. 12
4.2 The algorithm with a known condition-number bound
- With a valid condition-number bound, one call to SolvePerturbed yields an original-system relative residual at most epsilon/2 after rounding and rescaling.The prescribed scales ensure the certificate holds, so no precision search or acceptance test is needed.Source: §4.2 · p. 13
4.3 Removing the condition-number bound
- Without a condition-number bound, the algorithm doubles the precision until its certificate passes, without computing or estimating the condition number.The test is sufficient rather than necessary, so rejection does not establish that the candidate's actual residual is too large.Source: §4.3 · p. 14
4.4 Correctness and complexity
- At precision h, each exact solve costs Õ(mnh) bit operations, and repeated input conversion adds only a logarithmic factor to input reading.The analysis sums the work over trials to obtain the stated total complexity.Source: §4.4 · p. 14
- For fixed-point inputs, with τ = max{B, log κ, E}, the algorithm uses Õ(mnτ) bit operations, Õ(n^2τ) storage, and Õ(τ) bits per output coordinate.The compact output representation avoids expanding large encoded powers of two, which can require Ω(D) bits even when κ = 1.Source: §4.4 · pp. 14–15
5.1 AI Disclosure
The disclosure describes the author's use of language models during exploration and writing. Their assistance included testing approaches, literature search, proofreading, and parts of the floating-point rounding analysis.
5.1 AI Disclosure
- The author used language models to explore and numerically test approaches, assist with literature search and proofreading, and help write parts of the rounding analysis.The author expected the floating-point rounding analysis to be messier than it proved to be.Source: §5.1 · p. 15