Wilkinson's Polynomial (1963)
A classic ill-conditioned root-finding case from James Wilkinson (1963): the roots are the integers 1 through 20, but in FP64 small coefficient errors misplace several of them. Exact arithmetic recovers all 20.
The Problem
w(x) = (x-1)(x-2)(x-3)...(x-20)
The roots are OBVIOUS: 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20
But expand to coefficient form and solve numerically in FP64, and several roots move substantially — some becoming complex.
What FP64 Does
- • Roots 16, 17 become COMPLEX numbers
- • Root 20 shifts by more than 10 units
- • Some roots disappear entirely
- • Residuals explode to 106 or higher
| True Root | FP64 Result | VLA Result | Status |
|---|---|---|---|
| 15 | 14.93 + 0.0i | 15 | |
| 16 | 16.73 + 2.81i | 16 | |
| 17 | 16.73 - 2.81i | 17 | |
| 18 | 19.50 + 0.0i | 18 | |
| 19 | 20.85 + 0.0i | 19 | |
| 20 | 8.91 + 0.0i | 20 |
4+
FP64 wrong roots
20/20
VLA exact roots
Why This Happens
Coefficient Explosion
When you expand (x-1)(x-2)...(x-20), the coefficients become MASSIVE:
- • x19 coefficient: -210
- • x10 coefficient: 101,623,020
- • x1 coefficient: 20!
- • Constant term: -2,432,902,008,176,640,000
The Sensitivity
Changing the x19 coefficient by just 2-23 (one bit in FP32) causes roots 16 and 17 to become complex conjugates.
FP64 has only ~15 decimal digits. The polynomial needs 20+ digits to represent its coefficients exactly, so in FP64 the recovered roots are well off the true values.
The VLA Solution
FP64 Approach
15-16 decimal digits of precision. Coefficients rounded. Errors compound through polynomial evaluation. Roots destroyed.
VLA Approach
77 decimal digits (512-bit integers). Coefficients stored EXACTLY. Exact polynomial evaluation. All 20 integer roots recovered perfectly.
Background
Wilkinson's polynomial is a long-standing textbook benchmark for numerical stability and conditioning.
The standard remedies are higher precision or symbolic computation. The question this demo explores is whether exact arithmetic can recover all 20 roots while staying GPU-accelerated — here it does.
Have Ill-Conditioned Problems?
If Wilkinson's polynomial breaks your solver, what else is silently failing? Let us validate your numerical pipelines.
Get in touch