Exact Arithmetic: Big Integers, Rationals, and Floating-Point Traps
When \(10^{16}+1\) equals \(10^{16}\), and what to do about it: IEEE-754, catastrophic cancellation, exact rationals, and integer-only algorithms.
Leads to: 19.12 needs exact rational convergents; the Euler Lab’s Number Theory and Big-Number paths punish anyone who reaches for floats.
Learning Objectives
Click a status chip to cycle: Not started → In progress → Studied → Practiced → Needs review → Mastered.
- Explain the IEEE-754 double layout and derive why \(2^{53}\) is the last consecutively-representable integer.
- Diagnose catastrophic cancellation and rewrite an expression to avoid it.
- Decide correctly between int, Fraction, Decimal, and float for a given computation.
- Implement an exact integer square root and an exact rational comparison, avoiding all float contamination.
- Estimate the cost of big-integer arithmetic and explain when digit-count itself becomes the bottleneck.
Key Vocabulary
- IEEE-754 double
- 64 bits = 1 sign + 11 exponent + 52 stored mantissa bits; represents \(\pm(1.f)\times2^{e}\) with 53 bits of significand precision.
- Machine epsilon
- \(\varepsilon=2^{-52}\approx2.22\times10^{-16}\): the gap between 1 and the next representable double. Relative rounding error is \(\le\varepsilon/2\).
- Catastrophic cancellation
- Subtracting two nearly equal inexact numbers: the leading digits cancel, promoting tiny relative errors into huge ones.
- Arbitrary-precision integer
- Python’s int: exact, unbounded, cost grows with digit count (\(O(d^2)\) schoolbook, \(O(d^{1.585})\) Karatsuba for big \(d\)).
- Rational arithmetic
- fractions.Fraction: exact \(p/q\) kept in lowest terms via gcd - correct, but the denominators can explode.
- Integer square root
- math.isqrt(n): the exact \(\lfloor\sqrt n\rfloor\) with no floating point anywhere; the only safe way to test perfect squares for huge \(n\).
Intuition & Motivation
What a double actually is
A finite double is \((-1)^s\cdot(1.b_1b_2\ldots b_{52})_2\cdot 2^{e-1023}\). The significand has 53 bits of precision (52 stored plus the implicit leading 1). Hence:
Above \(2^{53}\) the representable integers thin out to every 2nd, then every 4th, and so on. This is why float(10**17)+1 == float(10**17) is True, and why a ‘sum of digits of \(2^{1000}\)’ computed via floats is not merely inaccurate but meaningless.
Catastrophic cancellation
Consider \(f(x)=\sqrt{x+1}-\sqrt{x}\) at \(x=10^{16}\). Both square roots are about \(10^8\) and agree to about 16 significant digits - so their difference in double precision keeps essentially zero correct digits. The fix is algebraic, not numerical:
which involves no subtraction of near-equal quantities and is accurate to full precision. Same trick as the stable quadratic formula: never subtract two nearly-equal inexact numbers if algebra can rewrite the expression as a sum or a quotient.
Choosing your number type
| Need | Type | Cost | Trap |
|---|---|---|---|
| Exact integers of any size | int | grows with digits | digit count itself can dominate |
| Exact \(p/q\) | fractions.Fraction | gcd per operation | denominators blow up; reduce or switch to modular |
| Exact decimal (money, digit puzzles) | decimal.Decimal | slow, configurable precision | still finite precision - set it deliberately |
| Continuous math, speed | float | 1 CPU instruction | 53 bits only; cancellation; \(0.1+0.2\ne0.3\) |
| Exact \(\lfloor\sqrt n\rfloor\) | math.isqrt | \(O(\log n)\) iterations | using sqrt instead - the single most common Euler bug |
Coding workbench
The Euler Lab’s Big-Number and Number Theory paths contain a whole family of problems (large powers, factorial digit sums, huge sums of roots) whose only real difficulty is refusing to use floats. Python hands you exact integers for free; the skill is noticing when you have silently left them.
- Testing math.sqrt(n)**2 == n or n**0.5 for large \(n\). Use math.isqrt.
- Using // on floats, or / on ints when you meant exact division: 10**20/10 is a float and already wrong.
- Accumulating a long sum of floats in arbitrary order and expecting the last digits to mean anything (see Kahan summation for the partial cure).
- Reaching for Fraction inside a hot loop: every operation does a gcd, and denominators grow multiplicatively. Prefer integer-only reformulations, or work modulo a prime.
- Believing Decimal is ‘exact’. It is exact for decimal fractions at a chosen precision; \(1/3\) is still truncated.
- The 100-digit-precision problems in any problem set are a hint, not a taunt: the intended tool is exact integer arithmetic or Decimal with a deliberately-set precision.
- When a formula requires \(\sqrt{2}\) to 100 places, do not compute it in floats and hope - use Decimal with getcontext().prec set well above the digits you need, or use integer square roots of \(2\cdot10^{2k}\).
- Big-integer multiplication is not free: multiplying two \(d\)-digit numbers is \(O(d^{1.585})\) (Karatsuba, which CPython uses above a threshold). A loop that repeatedly squares a growing integer is quadratic in disguise.
- Trefethen’s rule of thumb: an algorithm is backward stable if it gives the exact answer to a nearby problem. In discrete problems you can do better than stable - you can be exact.
Practice this in the Euler Lab
Computational problems that exercise exactly this technique. Each opens in the Euler Lab with a Python workbench, a progressive hint ladder, and answer checking. Tier A/B run at full scale in the browser.
Warm-up:
#16 Power Digit Sum (1%, tier A) #20 Factorial Digit Sum (1%, tier A) #34 Digit Factorials (2%, tier A) #33 Digit Cancelling Fractions (3%, tier A)
Applied:
#816 Shortest Distance Among Points (17%, tier B) #323 Bitwise-OR Operations on Random In (18%, tier B) #479 Roots on the Rise (19%, tier B)
Challenge:
#666 Polymorphic Bacteria (41%, tier C) #869 Prime Guessing (41%, tier C)
222 Project Euler problems in total are mapped to this lesson. Open the Euler Lab to filter them all.
Knowledge Check
Practical Exercise
(a) Show that the naive quadratic root \(x_1=\frac{-b+\sqrt{b^2-4ac}}{2a}\) suffers catastrophic cancellation when \(b\gt 0\) and \(4ac\ll b^2\), and derive a stable alternative. (b) You must sum \(\sum_{n=1}^{10^6} 1/n^2\) to as many correct digits as possible in double precision. In which order do you add, and why? (c) Give an example where Fraction is exactly right and one where it is a trap.
(a) With \(b\gt 0\) and \(4ac\) small, \(\sqrt{b^2-4ac}\approx b\), so the numerator \(-b+\sqrt{b^2-4ac}\) subtracts two nearly-equal quantities: the leading digits cancel and the relative error explodes. Stable route: compute the well-conditioned root first, \(q=-\tfrac12\big(b+\operatorname{sign}(b)\sqrt{b^2-4ac}\big)\) (an addition of same-sign quantities, no cancellation), then get the other from Vieta’s \(x_1x_2=c/a\):
(b) Add from the SMALLEST term upward (\(n=10^6\) down to 1). When you accumulate largest-first, the running sum is \(O(1)\) while the increments are \(O(10^{-12})\), so each tiny term is rounded away against a much larger accumulator, and you lose most of the tail. Summing small-first lets the small terms combine with each other at comparable magnitude before meeting the big ones. (The bullet-proof answer is math.fsum, which is exactly-rounded regardless of order; the exact answer here is close to \(\pi^2/6\approx1.6449340668\), and the two summation orders differ in the last few digits.)
(c) Right: comparing two rationals exactly, e.g. deciding whether \(p_1/q_1 \gt p_2/q_2\) for the convergents of 19.12 - or accumulating an exact probability from small denominators. (Even better: compare \(p_1q_2\) vs \(p_2q_1\) as integers - same answer, no gcd cost.) Trap: iterating a map like \(x\mapsto x/2 + 1/x\) as Fractions: the denominator roughly doubles in size each step, so after 50 steps you are doing arithmetic on numbers with \(10^{15}\) digits. Exactness is not free.
Lesson Summary
Retrieval Practice
Close the lesson and answer from memory before checking. This is deliberate, effortful recall - the single highest-yield study action.
A: \(2^{53}\), because the significand holds 53 bits; at \(2^{53}\) the spacing between representable doubles becomes 2, so \(2^{53}+1\) rounds back down.
A: Multiply by the conjugate: \(1/(\sqrt{x+1}+\sqrt x)\). For the quadratic, compute the same-sign root \(q=-\frac12(b+\mathrm{sign}(b)\sqrt{b^2-4ac})\) then use Vieta (\(x_2=c/q\)) for the other.
Completion Checklist
- I can explain the core ideas in my own words
- I worked the derivations/examples by hand
- I completed the interactive workbench(es)
- I passed the knowledge check
Source References
This lesson synthesizes and paraphrases concepts from the sources below. No copyrighted text, problem sets, or solutions are reproduced. Return to the originals for full depth.
- Numerical Linear Algebra (Trefethen & Bau, SIAM, 1997) foundational - Lectures 13–15 - Floating-point arithmetic, machine epsilon, conditioning and backward stability - the rigorous version of everything in this lesson.
- The Pragmatic Programmer (Thomas & Hunt, 2nd/20th Anniv. ed., 2019) current - Ch. 3 & Ch. 7 - Choose the right representation and make illegal states unrepresentable; a float where an integer belongs is a design bug, not a rounding bug.
- Calculus, Vol. I (Tom Apostol, 2nd ed., 1967) foundational - Ch. I.3 (the real number system) - Rigorous treatment of the real numbers and of limits - the mathematical objects that floating point can only approximate.