How large is the growth factor of random matrices under partial pivoting?

Growth factor measures entry blow-up during elimination

Gaussian elimination factors a matrix A into a lower triangular L and upper triangular U. In finite precision, the accumulated round-off error is controlled by the growth factor growth(A) = max(‖L‖_max, ‖U‖_max)/‖A‖_max, where ‖·‖_max is the maximum magnitude entry. If growth(A) is not too large, the LU factorization is stable; roughly log(n·growth(A)) bits of precision suffice.

Partial pivoting is the greedy row-swapping rule used by numpy.linalg.solve in Python and the backslash command in Matlab and Julia. At each elimination step the largest-magnitude entry in the current first column is moved into the first row. The worst-case growth factor under this rule can be exponentially large: a doubling argument gives an upper bound of 2^{n−1}, and Wilkinson constructed matrices in 1961 that achieve this in every dimension. The puzzle is that such blow-up is rarely observed in practice.

Trefethen's conjecture set the √n target

Trefethen and Schreiber ran empirical average-case analyses of Gaussian elimination, and Edelman's numerical experiments suggested the growth factor of a Gaussian matrix under partial pivoting grows like n^{1/2}. Edelman remarked that proving such a result would probably be very difficult. Trefethen stated the conjecture at the 1995 International Congress on Industrial and Applied Mathematics and in the 1997 text Numerical Linear Algebra, and later offered $1000 for a resolution.

Before this paper, mathematical progress came on two fronts. Trefethen claimed an unpublished proof for the last pivot |U_nn|, released 14 years later in parallel with this manuscript. Huang and Tikhomirov proved that a Gaussian matrix under partial pivoting has growth at most n^c with high probability for a large unspecified constant c. Sankar's smoothed analysis showed that any matrix plus a small Gaussian perturbation A+εG with ‖A‖_2≤1 has growth factor under partial pivoting bounded by (n/ε)^{O(log n)} with high probability.

Partial pivoting keeps Gaussian growth near √n

Theorem 1.4 states that for every α>0 there exist constants C_α and n_α such that P[growth(Π(G)) > √n·e^{C_α√(log n)}] < n^{−α} for all n≥n_α, where G is an n×n Gaussian matrix and Π(·) is the partial pivoting map. Because e^{√(log n)} = o(n^c) for every c>0, this implies Trefethen's Conjecture 1.1.

The 1/2 exponent is correct, not merely an upper bound. The paper notes that log‖G‖_max = o(log n), and by a cited result log‖U‖_max ≥ log(|det(G)|^{1/n}) ≥ (1/2)log n − O(1) with probability 1−o(1). Combining this lower bound with Theorem 1.4 yields Corollary 1.5: log growth(Π(G))/log n converges in probability to 1/2. So the growth factor is genuinely of order n^{1/2} up to subpolynomial factors, not smaller.

Orthogonal preconditioning produces larger growth than partial pivoting

The paper also proves an upper bound for a matrix left-multiplied by a Haar orthogonal matrix without pivoting. Theorem 1.2 states that for A∈GL_n(R) and Q∼Haar(O(n)), P[growth(QA) > 14n^{3/2}√(log(en))/δ] < δ for all δ∈(0,1). The 3/2 exponent is tight for any non-adaptive orthogonal preconditioner.

Theorem 1.3 characterizes the asymptotic distribution of the growth factor of a Gaussian matrix without pivoting: n^{−3/2}growth(G) converges in distribution to X∼InvExp(2/(3√π)), meaning Pr[X≤x] = e^{−2/(3√π x)} for x>0. Numerical experiments with one million samples at n=64 give a Kolmogorov–Smirnov distance of 2.8×10^{−3} between the empirical distribution and the inverse exponential, suggesting rapid practical convergence even for small n.

The comparison matters because orthogonal preconditioning is a reasonable alternative to partial pivoting. The paper states that this alternative produces a larger growth factor than one expects with partial pivoting on a Gaussian matrix, and its tail bounds decay inverse polynomially, which is slower than expected for partial pivoting. The scaled growth factor n^{−3/2}growth(G) is Θ(1) with probability n^{−O(1)}, whereas n^{−1/2}growth(Π(G)) is n^{o(1)} with overwhelming probability.

Smoothed analysis improves the perturbation bound to n^{O(log log n)}

Theorem 1.6 extends the argument to the smoothed setting. For every α>0 there exist constants C_α and n_α such that for all n≥n_α, every matrix A∈R^{n×n} with ‖A‖_2≤1 and every ε∈(0,1) satisfy P[growth(Π(A+εG)) > n^{C_α}ε^{−2/log 2·log log n}] < n^{−α}, where G is a Gaussian matrix and Π(·) is the partial pivoting map.

This shows that polynomially small perturbations produce n^{O(log log n)} growth. The paper states that the correct behavior is likely n^{O(1)}, so this is a significant improvement over the longstanding n^{O(log n)} bound from Sankar. The proof adapts the partial pivoting argument to a general mean by replacing A+εG with A/ε+G, which keeps unit variance but allows a possibly very large mean; the additional loss is tracked through a factor f_K(r) that accumulates the extra term.

The proof controls truncated inverses of Schur complements

The main technical tool is a bound on the inverse of a submatrix of a Haar orthogonal matrix. Lemma 2.1 states that for Q∼Haar(O(n)), fixed subsets S,T⊂[n] of size m, and fixed vectors x,y_1,…,y_k∈S^{m−1}(R) with k≤m, P[max_{j∈[k]}|y_j^T(Q_{S,T})^{−1}x| > t] ≤ (6/t)√((n−m+1)log(ek)/m) for all t>0. Repeated application of this lemma proves Theorem 1.2.

For the partial pivoting result, the paper proves six technical lemmas about singular values and Gaussian matrices, then controls the truncated inverse Σ_r(T_t) = Σ_{i=1}^{n−k} 1/σ_i(T_t)^2 of the leading t×t pivot submatrix. Proposition 5.5 shows that max_{k<n} Σ_r(Π(G)_{[k],[k]}) ≤ n·e^{C_α√(log n)} with probability at least 1−n^{−(α+1)} when r=⌈√(log n)⌉. Proposition 5.1 converts truncated inverse control into a growth bound. Combining the two with exponent α+1 gives Theorem 1.4.

The proof of Theorem 1.4 was developed through human–AI interaction. The AI usage statement says the proof of Conjecture 1.1 began as a human-only proof sketch for the bottom-right o(n) entries of U_ij and was improved to a proof for all entries through AI interaction. Conjecture 1.1 had been listed on the GitHub repository OpenProblemsInNLA as open, challenging, and interesting to the community, and remained unsolved prior to this manuscript.

Numerics hint at Gumbel fluctuations and possible narrowing

The paper reports numerical experiments that go beyond what the theorems prove. For n=4^k with k=3,…,7, it computes 7500 samples of n^{−1/2}growth(Π(G)). The empirical densities have rapidly decaying tails and appear to narrow somewhat as n increases. The paper says the evidence is not overwhelming, but empirically the fluctuations of n^{−1/2}growth(G), properly shifted and normalized, seem to resemble Gumbel distributions, which would not be too surprising since such distributions often appear when taking the maximum of many random variables.

The paper also notes that as n increases, the empirical density of n^{−1/2}growth(G) appears to narrow, leaving open the possibility that n^{−1/2}growth(G) converges to a constant. These are empirical observations, not proved results, and the paper explicitly flags the evidence as not overwhelming. They are worth separating from Theorem 1.4, which is a rigorous high-probability upper bound, and from Corollary 1.5, which is a rigorous convergence-in-probability statement about the exponent.

The guarantee is for Gaussian matrices and exact-arithmetic growth

Theorem 1.4 applies to n×n Gaussian matrices, meaning matrices whose entries are independent standard normal random variables, under the partial pivoting map Π(·). It is an asymptotic statement: the constants C_α and n_α depend on α, and the bound holds for all n≥n_α. The theorem does not claim a specific finite-n threshold, and the numerical figures are for n up to 4^7=16384.

A subtlety the paper flags is that the growth factor needed for error analysis is the one under floating point arithmetic, which can differ from the exact-arithmetic growth factor due to round-off accumulation or a change in pivot selection. The paper notes that the worst-case growth factor for partial pivoting is the same in exact and floating point arithmetic, and that the growth factor of a Gaussian matrix under partial pivoting is increased by at most a small multiplicative constant with high probability. The results in this paper are for exact arithmetic.

The smoothed analysis result in Theorem 1.6 applies to any fixed A∈R^{n×n} with ‖A‖_2≤1 plus an ε-sized Gaussian perturbation, for ε∈(0,1). It is also asymptotic in n. The paper states that the correct behavior is likely n^{O(1)}, so the n^{O(log log n)} bound is an improvement rather than a final answer. The orthogonal preconditioning results in Theorems 1.2 and 1.3 concern growth without pivoting, a different algorithmic setting from the partial pivoting results, and should not be read as bounds on growth(Π(G)).

Sources used in this answer