The largest 5th pivot may be the root of a 61st degree polynomial

Abstract: This paper introduces a number of new techniques in the study of the famous question from numerical linear algebra: what is the largest possible growth factor when performing Gaussian elimination with complete pivoting? This question is highly complex, due to a complicated set of polynomial inequalities that need to be simultaneously satisfied. This paper introduces the JuMP + Groebner basis + discriminant polynomial approach as well as the use of interval arithmetic computations. Thus, we are introducing a marriage of numerical and exact mathematical computations. In 1988, Day and Peterson performed numerical optimization on $n=5$ with NPSOL and obtained a largest seen value of $4.1325...$. This same best value was reproduced by Gould with LANCELOT in 1991. We ran extensive comparable experiments with the modern software tool JuMP and also saw the same value $4.1325...$. While the combinatorial explosion of possibilities prevents us from knowing there may not be a larger maximum, we succeed in obtaining the exact mathematical value: the number $4.1325...$ is exactly the root of a 61st degree polynomial provided in this work, and is a maximum given the equality constraints seen by JuMP. In light of the numerics, we pose the conjecture that this lower bound is indeed the maximum. We also apply this technique to $n = 6$, $7$, and $8$. Furthermore, in 1969, an upper bound of $4\frac{17}{18}\approx 4.94$ was produced for the maximum possible growth for $n = 5$. We slightly lower this upper bound to $4.84$.

• The paper identifies the numerically observed 5×5 growth factor, 4.13251707863247…, as the unique real root in [4,5] of an explicit 61st-degree integer polynomial.
• The authors combine JuMP optimization, Gröbner bases, discriminants, and exact case analysis to prove this value is optimal within the observed equality-constraint structure, while global optimality remains conjectural.
• The paper improves the rigorous n=5 upper bound from approximately 4.944 to 4.84 and derives new lower bounds for n=6, 7, and 8, including 5, 6.056953…, and 8.

This paper addresses a long-standing open problem in numerical linear algebra: determining the maximum growth factor for Gaussian elimination with complete pivoting on $5 \times 5$ real matrices. The authors, Chen, Edelman, and Urschel, establish that the numerically observed value $4.1325\ldots$ — reproduced independently by NPSOL in 1988, LANCELOT in 1991 and 2026, and JuMP in 2024 — is exactly the unique real root in $[4,5]$ of an explicit 61st-degree polynomial with integer coefficients. They prove this value is optimal within a specific equality-constraint structure suggested by numerical optimization, conjecture it is the global maximum, and improve the best known upper bound for $n=5$ from Tornheim's $4\tfrac{17}{18} \approx 4.9\overline{4}$ to $4.84$.

Background and problem formulation

For a matrix $A$, Gaussian elimination produces a factorization via iterated Schur complements; the growth factor

$$g(A) = \max\left\{ \|L\|_{\max},\ \frac{\|U\|_{\max}}{\|A\|_{\max}} \right\}$$

governs worst-case backward stability of the elimination. Complete pivoting selects at each step the largest-magnitude entry of the active submatrix as pivot. For completely pivoted matrices, maximizing growth reduces to maximizing the last pivot $|U_{nn}|$. The maximum growth factors are known exactly for $n = 1,\ldots,4$: the sequence $4.1325\ldots$0. The value for $4.1325\ldots$1 has resisted exact determination for decades despite repeated numerical attacks.

A structural observation underpins the algebraic approach: by the Tarski–Seidenberg theorem, the maximum growth factor over $4.1325\ldots$2 completely pivoted matrices is always an algebraic number. The authors prove it has degree at most $4.1325\ldots$3 — for $4.1325\ldots$4, a bound of $4.1325\ldots$5 that appears strikingly pessimistic against the degree-61 answer they obtain.

The JuMP + Gröbner basis + discriminant method

The central methodological contribution is a hybrid pipeline combining floating-point optimization with exact symbolic computation:

1. Numerical optimization with JuMP suggests both a candidate extremal matrix and, crucially, the equality structure of tight constraints at the optimum (which entries of the intermediate Schur complements attain their pivoting bounds).
2. Exact elimination: tight constraints are imposed as equalities; slack entries become variables solved symbolically from linear constraints.
3. Gröbner basis computation eliminates remaining nonlinearly dependent variables.
4. Discriminant polynomials perform the final exact optimization: if $4.1325\ldots$6 is defined implicitly by $4.1325\ldots$7, stationary points of $4.1325\ldots$8 as a function of $4.1325\ldots$9 satisfy $[4,5]$0, so candidate optima are roots of $[4,5]$1, a univariate polynomial in $[4,5]$2.

All steps after the first are exact — no rounding errors are involved. The authors emphasize that the combinatorial explosion of possible equality structures prevents exhaustive treatment; the numerical optimizer instead nominates which structure to analyze rigorously.

The 61st-degree polynomial and its optimality within the observed structure

For the elimination sequence observed consistently in JuMP runs (14 of 25 entries fixed at $[4,5]$3, ten tight pivoting constraints), the computation reduces to three polynomial equations $[4,5]$4 in variables $[4,5]$5. A Gröbner basis eliminating $[4,5]$6 yields a factorized polynomial; taking the discriminant of the operative factor with respect to $[4,5]$7 produces

$[4,5]$8

where $[4,5]$9 is the 61st-degree polynomial displayed in the paper. Descartes' rule of signs, applied after Möbius transformations mapping the intervals, confirms rigorously that $n=5$0 has exactly one real root in $n=5$1 and none in $n=5$2. That root is

$n=5$3

agreeing with Gould's quad-precision LANCELOT result to 30 digits. The extremal matrix, together with the auxiliary parameters $n=5$4 (each itself a root of a 61st-degree polynomial), is given explicitly.

Theorem 3.1 proves this value is maximal over all completely pivoted matrices satisfying the observed equality constraints. The proof systematically dismisses every other factor of the discriminant as extraneous: $n=5$5 cannot be maximal; degenerate cases ($n=5$6, forced values of $n=5$7 or $n=5$8) either blow up entries or cap growth at $n=5$9; and the remaining octic factor in $4\tfrac{17}{18} \approx 4.9\overline{4}$0 admits only one root in $4\tfrac{17}{18} \approx 4.9\overline{4}$1, which cannot yield sufficient growth. This case analysis is computer-assisted but rests on exact symbolic computations throughout.

The main limitation is candidly stated: because the super-exponentially many alternative constraint configurations are not ruled out, the global optimality remains a conjecture, not a theorem. The resolution "will likely require an impressive marriage of mathematics and computation, one that may not yet be possible."

Results for $4\tfrac{17}{18} \approx 4.9\overline{4}$2

Applying the same pipeline yields lower bounds for higher dimensions. For $4\tfrac{17}{18} \approx 4.9\overline{4}$3, JuMP's matrix has only six non-$4\tfrac{17}{18} \approx 4.9\overline{4}$4 entries, so no Gröbner basis is needed; the optimization collapses to maximizing $4\tfrac{17}{18} \approx 4.9\overline{4}$5 over $4\tfrac{17}{18} \approx 4.9\overline{4}$6, giving growth exactly $4\tfrac{17}{18} \approx 4.9\overline{4}$7. For $4\tfrac{17}{18} \approx 4.9\overline{4}$8, the JuMP matrix has at least $4\tfrac{17}{18} \approx 4.9\overline{4}$9 tight constraints, so after a Gröbner basis computation no optimization is required; the growth is the unique real root of the sextic

$4.84$0

namely $4.84$1 — notably larger than the value $4.84$2 reported by Day–Peterson and Gould, confirming the 2024 finding that the $4.84$3 maximum exceeds $4.84$4. For $4.84$5, numerics suggest growth $4.84$6, achieved by Sylvester's Hadamard matrix, though non-Hadamard maximizers also exist. If the lower bounds are true maxima, the sequence begins $4.84$7 — an irregular pattern whose form could scarcely have been predicted.

Improved upper bound of 4.84

The second contribution narrows the gap from above. Prior bounds were $4.84$8 (composing two $4.84$9 bounds), Cohen's $A$0, and Tornheim's $A$1 from 1969. The new bound combines analysis with interval arithmetic branch-and-bound:

• Reduction lemma: any completely pivoted normalized $A$2 matrix with third pivot $A$3 induces a completely pivoted normalized $A$4 matrix $A$5 satisfying $A$6 for bounded vectors $A$7, with $A$8 where $A$9 is $$g(A) = \max\left\{ \|L\|_{\max},\ \frac{\|U\|_{\max}}{\|A\|_{\max}} \right\}$$0's third pivot. Large-growth $$g(A) = \max\left\{ \|L\|_{\max},\ \frac{\|U\|_{\max}}{\|A\|_{\max}} \right\}$$1 blocks are thus incompatible with each other.
• Covers: interval-arithmetic pruning produces explicit unions of dyadic boxes — $$g(A) = \max\left\{ \|L\|_{\max},\ \frac{\|U\|_{\max}}{\|A\|_{\max}} \right\}$$2 with 1882 boxes and $$g(A) = \max\left\{ \|L\|_{\max},\ \frac{\|U\|_{\max}}{\|A\|_{\max}} \right\}$$3 with 19511 boxes — enclosing all high-growth normalized completely pivoted $$g(A) = \max\left\{ \|L\|_{\max},\ \frac{\|U\|_{\max}}{\|A\|_{\max}} \right\}$$4 matrices. These covers cluster around the three known $$g(A) = \max\left\{ \|L\|_{\max},\ \frac{\|U\|_{\max}}{\|A\|_{\max}} \right\}$$5 extremizers achieving growth $$g(A) = \max\left\{ \|L\|_{\max},\ \frac{\|U\|_{\max}}{\|A\|_{\max}} \right\}$$6.
• Infeasibility check: interval arithmetic verifies that inequality above cannot hold for any $$g(A) = \max\left\{ \|L\|_{\max},\ \frac{\|U\|_{\max}}{\|A\|_{\max}} \right\}$$7 in these covers when $$g(A) = \max\left\{ \|L\|_{\max},\ \frac{\|U\|_{\max}}{\|A\|_{\max}} \right\}$$8 would exceed $$g(A) = \max\left\{ \|L\|_{\max},\ \frac{\|U\|_{\max}}{\|A\|_{\max}} \right\}$$9.

Theorem 4.2 concludes $|U_{nn}|$0 for all completely pivoted $|U_{nn}|$1 matrices. The bound is not tight by design: the threshold constants can be lowered at the cost of runtime, but the technique "even when taken to its natural limit, certainly falls short of $|U_{nn}|$2" — closing the gap entirely requires further ideas. The computation used 64 threads for roughly 12 hours in Julia.

Limitations and open questions

Three limitations bear directly on the results. First, the optimality of $|U_{nn}|$3 is proven only conditional on the equality structure observed by JuMP; ruling out all other structures is computationally prohibitive, so Conjecture 3.2 remains open. Second, the upper-bound machinery depends on the reduction through $|U_{nn}|$4 blocks and on box covers at resolution $|U_{nn}|$5; finer resolutions are exponentially expensive in dimension, and the approach cannot reach the conjectured value even in principle. Third, the choice of thresholds $|U_{nn}|$6 and $|U_{nn}|$7 is admittedly arbitrary, tuned to runtime tolerance rather than to any structural feature of the problem. The specific question left open is whether the degree-61 lower bound equals the true $|U_{nn}|$8 maximum, and correspondingly whether the degree-6 lower bound is the true $|U_{nn}|$9 maximum.

Conclusion

The paper converts a four-decade-old numerical observation into an exact algebraic statement: the candidate $n = 1,\ldots,4$0 complete-pivoting growth factor is the root of an explicit 61st-degree integer polynomial, proven maximal within its constraint structure and supported by a rigorous root-count argument. The JuMP + Gröbner basis + discriminant methodology, together with interval-arithmetic upper bounding, provides a reproducible template applicable to larger dimensions, yielding new lower bounds for $n = 1,\ldots,4$1 and the first improvement since 1969 on the $n = 1,\ldots,4$2 upper bound. Whether the conjectured values are globally optimal remains the outstanding question.