Skip to main content

Exact Rational Arithmetic

Floating-point lies to you.

double x = 0.1 + 0.2;
std::cout << (x == 0.3);  // Prints 0 (false!)

The number 0.1 has no exact binary representation, for the same reason 1/3 has no exact decimal representation. Floating-point represents numbers as m x 2^e, and most decimal fractions don’t land on a power of two.

Rational arithmetic fixes this. 1/3 stays exactly 1/3.

The Representation

A rational number is a pair (numerator, denominator) kept in lowest terms:

template<std::integral T>
class rat {
    T num_;  // numerator (carries sign)
    T den_;  // denominator (always positive)

    void reduce() {
        T g = std::gcd(abs(num_), den_);
        num_ /= g;
        den_ /= g;
    }
};

Three invariants, always maintained:

  1. The denominator is positive (sign lives in the numerator)
  2. GCD(|num|, den) = 1 (always reduced)
  3. Zero is uniquely 0/1

Arithmetic

Addition needs a common denominator:

$$\frac{a}{b} + \frac{c}{d} = \frac{ad + bc}{bd}$$

Then reduce. In code:

rat operator+(rat const& rhs) const {
    return rat(num_ * rhs.den_ + rhs.num_ * den_,
               den_ * rhs.den_);
}

The constructor calls reduce() automatically.

Multiplication is simpler:

$$\frac{a}{b} \times \frac{c}{d} = \frac{ac}{bd}$$

Division multiplies by the reciprocal:

$$\frac{a/b}{c/d} = \frac{ad}{bc}$$

Exact Comparison

No floating-point fuzziness. Two reduced rationals are equal iff their numerators and denominators match:

bool operator==(rat const& rhs) const {
    return num_ == rhs.num_ && den_ == rhs.den_;
}

For ordering, cross-multiply (valid because denominators are positive):

$$\frac{a}{b} < \frac{c}{d} \iff ad < cb$$

The Mediant

The mediant of a/b and c/d is (a+c)/(b+d). It’s not the average. It has different, more interesting properties:

rat mediant(rat const& a, rat const& b) {
    return rat(a.numerator() + b.numerator(),
               a.denominator() + b.denominator());
}

If a/b < c/d, then a/b < mediant < c/d. The mediant is always in lowest terms when a/b and c/d are neighbors in the Stern-Brocot tree. And mediants generate all positive rationals exactly once.

The Stern-Brocot Tree

Start with 0/1 and 1/0 (representing infinity). Repeatedly take mediants:

Level 0:     0/1                     1/0
Level 1:     0/1       1/1           1/0
Level 2:     0/1   1/2   1/1   2/1   1/0
Level 3:  0/1 1/3 1/2 2/3 1/1 3/2 2/1 3/1 1/0

Every positive rational appears exactly once. The path from root to any node encodes its continued fraction. This connects to best rational approximations and Farey sequences.

GCD Ties Everything Together

Reducing fractions requires GCD. The algorithm is Euclid’s, from around 300 BCE:

T gcd(T a, T b) {
    while (b != 0) {
        a = a % b;
        std::swap(a, b);
    }
    return a;
}

The same algorithm works for any Euclidean domain. That’s not a coincidence. It’s a consequence of the algebraic structure.

Rational numbers form a field: every non-zero element has a multiplicative inverse (the reciprocal). The requirement that denominators be non-zero and fractions reduced comes from this algebraic structure, not from arbitrary convention.

Overflow

For large computations, numerators and denominators grow fast. Three options:

  1. Use int64_t and accept the limits
  2. Use arbitrary-precision integers
  3. Use interval arithmetic to track bounds

This implementation takes option one. Simplicity first.

Further Reading

  • Graham, Knuth, Patashnik, Concrete Mathematics, Chapter 4
  • Stern, “Uber eine zahlentheoretische Funktion” (1858)

Discussion