Programming for mathematics

Lens

MotivationH

A thousand cases in the blink of an eye

There is a famous story — told about Carl Friedrich Gauss as a schoolboy, though its details are uncertain — of a teacher who asked the class to add up all the numbers from $1$ to $100$ , hoping for a long quiet hour. Gauss answered almost at once: $5050$ .

Today a computer adds those hundred numbers in a millionth of a second. It could check the answer for every list from $1$ to $n$ , for every $n$ up to a million, before you finish this sentence. So why do mathematicians still bother with proofs?

Here is one reason. In 1752 Christian Goldbach suggested to Leonhard Euler that every odd number that is not prime can be written as a prime plus twice a square: $9 = 7 + 2 \cdot 1^2$ , $15 = 7 + 2 \cdot 2^2$ , $27 = 19 + 2 \cdot 2^2$ , and so on. It works for $9$ , $15$ , $21$ , and on through the hundreds and thousands. And then, at $5777$ , it fails.

And here is something stranger. Ask Python whether $0.1 + 0.2$ equals $0.3$ , and it answers False.

In this lesson you will learn to write short programs that do real mathematics — and learn exactly what a program's answer does, and does not, tell you.

IntuitionH

Programs are precise recipes

A program is a list of instructions precise enough for a machine to follow. You have met the main ingredients already in repetition and decisions; now we write them in a real programming language, Python.

A variable is a named box that holds a value. The value can change while the program runs.
A loop repeats some instructions, for example once for each number from $1$ to $100$ .
A conditional (an if) makes a decision: do this only when a condition is true.
A function packages a calculation under a name, so that you can use it again with different inputs.

DefinitionsH

Definitions

Definition— variable and assignment

A program variable is a name for a stored value that can change while the program runs. The assignment x = expression works out the value of the expression on the right, then stores it under the name on the left, replacing what was there before.

Definition— loop and conditional

The loop for i in range(a, b): runs its indented instructions once for each whole number i with $a \le i < b$ — notice that $b$ itself is not included. The conditional if condition: runs its indented instructions only when the condition is True.

Definition— program function

A program function, written with def, takes inputs (called arguments) and gives back an output with return. Unlike a mathematical function, it may also have other effects, such as printing.

Definition— binary representation

In base 2 (binary) the only digits are $0$ and $1$ , and the places count ones, twos, fours, eights, and so on, each place worth twice the one to its right:

1101_2 = 1 \cdot 8 + 1 \cdot 4 + 0 \cdot 2 + 1 \cdot 1 = 13.

A representation has no leading zeros if its first (leftmost) digit is $1$ ; the number $0$ is written $0$ .

Worked examplesH

Worked programs

Example— A loop with a running total

total = 0
for i in range(1, 101):
    total = total + i
print(total)        # prints 5050

Here is a trace — the value of each variable after each step — for the shorter loop range(1, 5):

step	`i`	`total` after the step
start	—	$0$
1	$1$	$1$
2	$2$	$3$
3	$3$	$6$
4	$4$	$10$

Example— A condition inside a loop

To add up only the multiples of $3$ below $20$ :

total = 0
for n in range(1, 20):
    if n % 3 == 0:          # % gives the remainder; == tests equality
        total = total + n
print(total)                # 3 + 6 + 9 + 12 + 15 + 18 = 63

Example— A function

def is_prime(n):
    if n < 2:
        return False
    for d in range(2, n):
        if n % d == 0:      # d divides n
            return False
    return True

print([p for p in range(30) if is_prime(p)])
# [2, 3, 5, 7, 11, 13, 17, 19, 23, 29]

The last line is a list comprehension: Python's version of set-builder notation $\{p : 0 \le p < 30,\ p \text{ is prime}\}$ (see set-builder notation).

Example— A sequence with two variables

The Fibonacci sequence starts $0, 1$ , and each further term is the sum of the two before it.

a, b = 0, 1
for step in range(10):
    print(a)
    a, b = b, a + b     # both new values are computed first, then stored

This prints $0, 1, 1, 2, 3, 5, 8, 13, 21, 34$ . The line a, b = b, a + b would be impossible as an equation (it would say $a = b$ and $b = a + b$ at once); as an assignment, it moves the pair one step along the sequence.

ExplanationH

Testing is not proving

Programs are wonderful for exploring. Add up the first few odd numbers:

1 = 1, \quad 1 + 3 = 4, \quad 1 + 3 + 5 = 9, \quad 1 + 3 + 5 + 7 = 16.

The answers are squares. A loop can check $1 + 3 + \dots + (2n - 1) = n^2$ for every $n$ up to a million in a moment. That is strong evidence, but it is not a proof: there are infinitely many $n$ , and a program can check only finitely many. Goldbach's conjecture above passed every test up to $5776$ , and then failed.

What a program can do decisively is find a counterexample. A single case where a “for all” statement fails refutes it for ever. And a program can suggest a pattern that we then prove. Here is a proof of the odd-number pattern, in two ways.

Theorem— Theorem 1 (sum of odd numbers)

For every whole number $n \ge 0$ : $\;1 + 3 + 5 + \dots + (2n - 1) = n^2$ .

Proof

By induction. For $n = 0$ the sum is empty and equals $0 = 0^2$ . Suppose the statement holds for some $n$ . The sum of the first $n + 1$ odd numbers is the sum of the first $n$ , plus the next odd number $2n + 1$ :

n^2 + (2n + 1) = (n + 1)^2.

So the statement holds for $n + 1$ as well. Starting from $0$ and stepping up one at a time, it holds for every $n$ .

By a picture. An $n \times n$ square of dots grows into an $(n + 1) \times (n + 1)$ square by adding an L-shaped border of $n + n + 1 = 2n + 1$ dots — the next odd number.

∎

The same kind of argument shows Gauss's formula $1 + 2 + \dots + n = \tfrac{n(n + 1)}{2}$ : write the sum forwards and backwards, and add — each of the $n$ columns adds up to $n + 1$ .

Note

A program that runs without an error message is not necessarily correct. If you write range(1, n) when you meant range(1, n + 1), the loop quietly stops one step early: it computes $1 + 2 + \dots + 99 = 4950$ instead of $5050$ , and nothing warns you. Checking a program against a case you can work out by hand, or against a proven formula, catches mistakes like this.

ProofH

Binary numbers

To write a number in binary, divide by $2$ again and again, writing down the remainders. Here is $45$ :

\begin{aligned} 45 &= 2 \cdot 22 + 1\\ 22 &= 2 \cdot 11 + 0\\ 11 &= 2 \cdot 5 + 1\\ 5 &= 2 \cdot 2 + 1\\ 2 &= 2 \cdot 1 + 0\\ 1 &= 2 \cdot 0 + 1 \end{aligned}

Reading the remainders from the bottom up gives $45 = 101101_2$ . Check: $32 + 8 + 4 + 1 = 45$ . In Python:

def to_binary(n):
    if n == 0:
        return "0"
    digits = ""
    while n > 0:
        digits = str(n % 2) + digits    # the new remainder goes in front
        n = n // 2                      # // divides and rounds down
    return digits

print(to_binary(45))    # 101101

Going the other way, read the digits from the left: start at $0$ , and for each digit, double what you have and add the digit. For $101101_2$ : $1, 2, 5, 11, 22, 45$ .

Theorem— Theorem 2 (binary representations exist and are unique)

Every whole number $n \ge 0$ has exactly one binary representation without leading zeros.

Proof

Existence. For $n = 0$ , the representation is $0$ . For $n \ge 1$ , divide by $2$ : $n = 2q + r$ with $r \in \{0, 1\}$ and $q < n$ . If $q = 0$ , then $n = r = 1$ , represented by $1$ . Otherwise $q$ , being smaller, already has a representation without leading zeros (this is the step that induction makes precise), and writing the digit $r$ after it gives a representation of $2q + r = n$ , still starting with $1$ .

Uniqueness. Suppose $n \ge 1$ had two representations without leading zeros. In each, every digit except the last one stands for an even amount (a multiple of $2$ ), so the last digit must be the remainder of $n$ on division by $2$ . So the last digits are equal. Removing them leaves two representations of the same smaller number $q = \tfrac{n - r}{2}$ — and by the same argument their last digits agree too. Repeating this, the two representations agree digit by digit, until both run out at the same moment, when the number reached is $0$ . (A representation without leading zeros of a number $\ge 1$ is never empty, and an empty list of digits stands for $0$ , so neither can run out before the other.)

∎

The same proof works in any base $b \ge 2$ , with the remainders $0, 1, \dots, b - 1$ as digits; base $10$ is the case you use every day.

ExplanationH

Why 0.1 + 0.2 is not 0.3

In base $10$ , the fraction $\tfrac{1}{3}$ never ends: $0.3333\ldots$ . In base $2$ , the same thing happens to $\tfrac{1}{10}$ :

\tfrac{1}{10} = 0.0001100110011\ldots_2,

with the block $0011$ repeating for ever. A computer's usual decimal numbers, called floating-point numbers, keep only $53$ significant binary digits. So Python stores not $\tfrac{1}{10}$ but the nearest number it can represent, which is exactly

\frac{3602879701896397}{36028797018963968} = 0.1000000000000000055511\ldots

The tiny errors in $0.1$ and $0.2$ do not cancel, and the sum is stored as $0.30000000000000004$ , which is not the stored version of $0.3$ . Nothing is broken: floating-point numbers are approximations, designed to be fast and almost always close enough. When you need exact answers, use exact arithmetic:

from fractions import Fraction
print(Fraction(1, 10) + Fraction(2, 10) == Fraction(3, 10))    # True
print(0.1 + 0.2 == 0.3)                                        # False

Whole numbers in Python are always exact, however large: 2 ** 100 is computed digit for digit.

Formal mathematicsF

Programs that come with proofs

Lean is a programming language too. Here is the summing loop, written almost exactly as in Python (:= is Lean's sign for assignment):

sumForLean 4Formally verified

/-- The sum `1 + 2 + ⋯ + n`, computed by a loop with a running total — the Lean version of the
Python program
```
total = 0
for i in range(1, n + 1):
    total = total + i
```
Here `total := total + i` is an **assignment**: replace the stored value of `total`. -/
def sumFor (n : Nat) : Nat := Id.run do
  let mut total := 0
  for i in [1:n+1] do
    total := total + i
  return total

Axioms used: propext, Quot.soundformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

sum_for_examplesLean 4Proved here from basic factsFormally verified

/-- Running the loop: `1 + ⋯ + 10 = 55` and `1 + ⋯ + 100 = 5050`. (`decide +kernel` asks Lean's
trusted kernel to run the program and compare the answers.) -/
theorem sum_for_examples : sumFor 10 = 55 ∧ sumFor 100 = 5050 := by
  decide +kernel

Axioms used: propext, Quot.soundformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

Assignment and equality really are different things:

assignment_exampleLean 4Proved here from basic factsFormally verified

/-- Assignment is not equality. The program "`x := 5`, then `x := x + 1`" ends with `x = 6`… -/
theorem assignment_example :
    (Id.run do
      let mut x := 5
      x := x + 1
      return x) = 6 := by
  rfl

Axioms used: No axiomsformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

no_solution_x_eq_x_add_oneLean 4Proved here from basic factsFormally verified

/-- … but as an *equation*, `x = x + 1` has no solution at all. -/
theorem no_solution_x_eq_x_add_one : ¬ ∃ x : Nat, x = x + 1 := by
  intro ⟨x, hx⟩
  omega

Axioms used: propext, Quot.soundformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

The same sum can be defined recursively: the sum up to $n + 1$ is the sum up to $n$ , plus $n + 1$ .

sumToLean 4Formally verified

/-- The same sum, defined **recursively**: the sum up to `0` is `0`, and the sum up to `n + 1` is
the sum up to `n`, plus `n + 1`. -/
def sumTo : Nat → Nat
  | 0 => 0
  | n + 1 => sumTo n + (n + 1)

Axioms used: No axiomsformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

loop_agrees_testedLean 4Proved here from basic factsFormally verified

/-- **Testing**: the loop and the recursive definition agree for each of the 40 inputs
`n = 0, 1, …, 39`. Lean checks this by running both programs on every one of these inputs. It
says nothing about `n = 40` or beyond. -/
theorem loop_agrees_tested : ∀ n, n < 40 → sumFor n = sumTo n := by
  decide +kernel

Axioms used: propext, Quot.soundformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

gauss_formulaLean 4Proved here from basic factsFormally verified

/-- **Proof for every `n`**: `2 · (1 + 2 + ⋯ + n) = n(n + 1)`, by induction. It holds for `n = 0`;
and if it holds for `n`, then `2 · sumTo (n + 1) = 2 · sumTo n + 2(n + 1) = n(n + 1) + 2(n + 1)
= (n + 1)(n + 2)`. -/
theorem gauss_formula (n : Nat) : 2 * sumTo n = n * (n + 1) := by
  induction n with
  | zero => rfl
  | succ k ih =>
    calc 2 * sumTo (k + 1) = 2 * sumTo k + 2 * (k + 1) := by
          rw [sumTo, Nat.mul_add]
      _ = k * (k + 1) + 2 * (k + 1) := by rw [ih]
      _ = (k + 1) * (k + 1 + 1) := by
          rw [← Nat.add_mul, Nat.mul_comm, Nat.add_comm k 2, Nat.add_comm 2 k]

Axioms used: propextformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

A loop that carries a running total can be written as a function with an extra argument, the accumulator, and then proved correct for every input:

sumAccLean 4Formally verified

/-- A loop written as a function with an **accumulator** `acc` (the running total): it adds
`n, n − 1, …, 1` to `acc`, one number per step. -/
def sumAcc : Nat → Nat → Nat
  | 0, acc => acc
  | n + 1, acc => sumAcc n (acc + (n + 1))

Axioms used: No axiomsformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

sum_acc_correctLean 4Proved here from basic factsFormally verified

/-- The accumulator loop is correct for **every** input: starting from any running total `acc`, it
ends with `acc + (1 + ⋯ + n)`. Proved by induction on `n`, for all `acc` at once. -/
theorem sum_acc_correct (n acc : Nat) : sumAcc n acc = acc + sumTo n := by
  induction n generalizing acc with
  | zero => rfl
  | succ k ih =>
    show sumAcc k (acc + (k + 1)) = acc + (sumTo k + (k + 1))
    rw [ih]
    omega

Axioms used: propext, Quot.soundformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

sumOddsLean 4Formally verified

/-- The sum of the first `n` odd numbers, `1 + 3 + 5 + ⋯ + (2n − 1)`. -/
def sumOdds : Nat → Nat
  | 0 => 0
  | n + 1 => sumOdds n + (2 * n + 1)

Axioms used: No axiomsformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

sum_odds_eq_squareLean 4Proved here from basic factsFormally verified

/-- The sum of the first `n` odd numbers is `n²` — for every `n`, by induction:
`n² + (2n + 1) = (n + 1)²`. -/
theorem sum_odds_eq_square (n : Nat) : sumOdds n = n * n := by
  induction n with
  | zero => rfl
  | succ k ih =>
    show sumOdds k + (2 * k + 1) = (k + 1) * (k + 1)
    rw [ih, Nat.add_mul, Nat.mul_add, Nat.mul_one, Nat.one_mul]
    omega

Axioms used: propext, Quot.soundformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

Binary numbers. Lean stores the digits as a list with the least significant digit first — the order in which repeated division produces them — so $13 = 1101_2$ is the list [1, 0, 1, 1].

fromBinaryLean 4Formally verified

/-- The value of a list of binary digits, **least significant digit first**: `[1, 0, 1, 1]` means
`1·1 + 0·2 + 1·4 + 1·8 = 13`. Each step: value of `d :: ds` is `d + 2 · (value of ds)`. -/
def fromBinary : List Nat → Nat
  | [] => 0
  | d :: ds => d + 2 * fromBinary ds

Axioms used: No axiomsformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

bitsAuxLean 4Formally verified

/-- Repeated division by `2`: the remainder is the next binary digit, and we continue with the
quotient. (The extra argument `fuel` limits the number of steps, which makes it easy for Lean to
see that the program stops; `n` steps are always enough.) -/
def bitsAux : Nat → Nat → List Nat
  | 0, _ => []
  | fuel + 1, n => if n = 0 then [] else (n % 2) :: bitsAux fuel (n / 2)

Axioms used: No axiomsformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

toBinaryLean 4Formally verified

/-- The binary digits of `n`, least significant first; `toBinary 0` is the empty list. -/
def toBinary (n : Nat) : List Nat := bitsAux n n

Axioms used: No axiomsformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

binary_examplesLean 4Proved here from basic factsFormally verified

/-- Examples, computed by running the programs: `13 = 1101₂` (digits `[1, 0, 1, 1]` from the
right), `45 = 101101₂`, and back again. -/
theorem binary_examples :
    toBinary 13 = [1, 0, 1, 1] ∧ toBinary 45 = [1, 0, 1, 1, 0, 1] ∧
      fromBinary [1, 0, 1, 1, 0, 1] = 45 := by
  decide

Axioms used: No axiomsformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

The conversion is correct for every number, not just the ones we tried:

from_bits_auxLean 4Proved here from basic factsFormally verified

/-- The heart of the correctness proof: if `n ≤ fuel`, the digits produced by repeated division
have value `n`. Induction on `fuel`: for `n > 0` the digits are `n % 2` followed by the digits of
`n / 2`, and `n % 2 + 2 · (n / 2) = n`. -/
theorem from_bits_aux : ∀ fuel n, n ≤ fuel → fromBinary (bitsAux fuel n) = n := by
  intro fuel
  induction fuel with
  | zero =>
    intro n hn
    show 0 = n
    omega
  | succ k ih =>
    intro n hn
    by_cases h0 : n = 0
    · subst h0
      rfl
    · show fromBinary (if n = 0 then [] else (n % 2) :: bitsAux k (n / 2)) = n
      split
      · contradiction
      · show n % 2 + 2 * fromBinary (bitsAux k (n / 2)) = n
        rw [ih (n / 2) (by omega)]
        omega

Axioms used: propext, Quot.soundformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

from_to_binaryLean 4Proved here from basic factsFormally verified

/-- **The conversion to binary is correct for every natural number**: the digits that `toBinary`
produces have value `n`. -/
theorem from_to_binary (n : Nat) : fromBinary (toBinary n) = n :=
  from_bits_aux n n (Nat.le_refl n)

Axioms used: propext, Quot.soundformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

And Theorem 2, existence and uniqueness, with “no leading zeros” stated precisely:

AllBitsLean 4Formally verified

/-- Every entry of the list is a binary digit, `0` or `1`. -/
def AllBits : List Nat → Prop
  | [] => True
  | d :: ds => d < 2 ∧ AllBits ds

Axioms used: No axiomsformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

NoLeadingZeroLean 4Formally verified

/-- No leading zeros: the most significant digit (the **last** entry, since we store the least
significant digit first) is not `0`. The empty list, which represents `0`, counts as fine. -/
def NoLeadingZero : List Nat → Prop
  | [] => True
  | [d] => d ≠ 0
  | _ :: d :: ds => NoLeadingZero (d :: ds)

Axioms used: No axiomsformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

value_ne_zeroLean 4Proved here from basic factsFormally verified

/-- A non-empty digit list without leading zeros has a non-zero value. -/
theorem value_ne_zero : ∀ (d : Nat) (ds : List Nat), NoLeadingZero (d :: ds) →
    fromBinary (d :: ds) ≠ 0
  | d, [], h => by
    show d + 2 * 0 ≠ 0
    have : d ≠ 0 := h
    omega
  | d, e :: es, h => by
    have ih := value_ne_zero e es h
    show d + 2 * fromBinary (e :: es) ≠ 0
    omega

Axioms used: propext, Quot.soundformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

binary_existsLean 4Proved here from basic factsFormally verified

/-- **Existence**: every natural number has a binary representation without leading zeros.
(Proved by induction, following the repeated-division algorithm: write `n = 2q + r` with
`r = n % 2`, represent `q`, and put `r` in front.) -/
theorem binary_exists : ∀ n, ∃ l, AllBits l ∧ NoLeadingZero l ∧ fromBinary l = n := by
  -- Induction on a bound `b` for `n`, so that we may use the statement for `n / 2`.
  have key : ∀ b n, n ≤ b → ∃ l, AllBits l ∧ NoLeadingZero l ∧ fromBinary l = n := by
    intro b
    induction b with
    | zero =>
      intro n hn
      exact ⟨[], trivial, trivial, by show 0 = n; omega⟩
    | succ k ih =>
      intro n hn
      by_cases h0 : n = 0
      · exact ⟨[], trivial, trivial, by show 0 = n; omega⟩
      · obtain ⟨l, hbits, hlead, hval⟩ := ih (n / 2) (by omega)
        refine ⟨(n % 2) :: l, ⟨Nat.mod_lt n (by decide), hbits⟩, ?_, ?_⟩
        · cases l with
          | nil =>
            -- Then `n / 2 = 0`, so `n = 1` and the single digit is `1`.
            have hq : n / 2 = 0 := by rw [← hval]; rfl
            show n % 2 ≠ 0
            omega
          | cons e es => exact hlead
        · show n % 2 + 2 * fromBinary l = n
          rw [hval]
          omega
  intro n
  exact key n n (Nat.le_refl n)

Axioms used: propext, Quot.soundformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

binary_uniqueLean 4Proved here from basic factsFormally verified

/-- **Uniqueness**: two binary representations without leading zeros that have the same value are
the same list. The last digit is forced (it is the remainder on division by `2`); remove it, and
the rest represent the same smaller number — so by induction they agree too. -/
theorem binary_unique : ∀ l₁ l₂ : List Nat, AllBits l₁ → AllBits l₂ →
    NoLeadingZero l₁ → NoLeadingZero l₂ → fromBinary l₁ = fromBinary l₂ → l₁ = l₂
  | [], [], _, _, _, _, _ => rfl
  | [], e :: es, _, _, _, h₂, hv => absurd hv.symm (value_ne_zero e es h₂)
  | d :: ds, [], _, _, h₁, _, hv => absurd hv (value_ne_zero d ds h₁)
  | d :: ds, e :: es, ⟨hd, hds⟩, ⟨he, hes⟩, h₁, h₂, hv => by
    have hv' : d + 2 * fromBinary ds = e + 2 * fromBinary es := hv
    have hde : d = e := by omega
    have hrest : fromBinary ds = fromBinary es := by omega
    have hlead₁ : NoLeadingZero ds := by
      cases ds with
      | nil => trivial
      | cons _ _ => exact h₁
    have hlead₂ : NoLeadingZero es := by
      cases es with
      | nil => trivial
      | cons _ _ => exact h₂
    rw [hde, binary_unique ds es hds hes hlead₁ hlead₂ hrest]

Axioms used: propext, Quot.soundformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

Finally, exact fractions behave exactly:

exact_tenthsLean 4Proved here from basic factsFormally verified

/-- With exact fractions, one tenth plus two tenths **is** three tenths. (Python's `Fraction` and
Lean's `Rat` both store numerator and denominator exactly; floating-point numbers do not.) -/
theorem exact_tenths : (1 / 10 + 2 / 10 : Rat) = 3 / 10 := by
  decide +kernel

Axioms used: propext, Classical.choice, Quot.soundformal/Lumenary/ProgrammingIntro.leanTry it in the Lean web editor

What is formalised: Theorem 1, Gauss's formula, the correctness of the accumulator loop and of the binary conversion for every input, and Theorem 2 (with $0$ represented by the empty list of digits), using core Lean only. The for loop sumFor is only tested against sumTo (for $n < 40$ ); proving it equal for every $n$ needs facts about how Lean runs loops that are beyond this lesson. Lean's Float type behaves like Python's floats — 0.1 + 0.2 == 0.3 evaluates to false there too — but its arithmetic is done by the computer's hardware, outside the logic, so Lean cannot prove facts about it; that is one reason careful mathematics uses exact numbers.

ComputationC

Experiments

Testing conjectures with loops, and a bug that runs finePythonExperimental evidenceCounterexampleExact computationComputationally tested

"""Testing conjectures with loops: confirmations, a famous failure, and a bug.

Loops make it easy to test a statement on thousands of cases. The experiment
checks two true formulas for n = 1..1000, then tests a conjecture Christian
Goldbach sent to Leonhard Euler in 1752 - every odd composite number is a
prime plus twice a square - and finds where it fails. Finally, a program
with an off-by-one bug runs without any error and prints a wrong answer.
Experimental evidence and counterexamples, not proofs.
"""

def sum_loop(n):
    """1 + 2 + ... + n with a loop and a running total."""
    total = 0
    for i in range(1, n + 1):
        total = total + i
    return total

def sum_loop_buggy(n):
    """Meant to compute 1 + ... + n, but range(1, n) stops at n - 1."""
    total = 0
    for i in range(1, n):
        total = total + i
    return total

def sum_of_odds(n):
    """1 + 3 + 5 + ... + (2n - 1)."""
    return sum(2 * k - 1 for k in range(1, n + 1))

def prime_sieve(limit):
    """is_prime[k] for 0 <= k <= limit (Sieve of Eratosthenes)."""
    is_prime = [False, False] + [True] * (limit - 1)
    for p in range(2, int(limit ** 0.5) + 1):
        if is_prime[p]:
            for multiple in range(p * p, limit + 1, p):
                is_prime[multiple] = False
    return is_prime

def prime_plus_twice_square(n, is_prime):
    """A pair (p, k) with n = p + 2*k*k and p prime, or None."""
    k = 0
    while 2 * k * k <= n:
        if is_prime[n - 2 * k * k]:
            return (n - 2 * k * k, k)
        k += 1
    return None

def goldbach_other_counterexamples(limit):
    """Odd composite numbers up to limit that are NOT a prime plus 2k^2."""
    is_prime = prime_sieve(limit)
    return [n for n in range(9, limit + 1, 2)
            if not is_prime[n] and prime_plus_twice_square(n, is_prime) is None]

def main():
    print("1. Testing two formulas for n = 1, 2, ..., 1000")
    ok1 = all(sum_loop(n) == n * (n + 1) // 2 for n in range(1, 1001))
    ok2 = all(sum_of_odds(n) == n * n for n in range(1, 1001))
    print("   1 + 2 + ... + n == n(n + 1)/2 every time:", ok1)
    print("   1 + 3 + ... + (2n - 1) == n^2 every time:", ok2)
    print()
    print("2. Goldbach (1752): every odd composite = prime + 2 * square?")
    is_prime = prime_sieve(10000)
    for n in [9, 15, 21, 25, 27, 33]:
        p, k = prime_plus_twice_square(n, is_prime)
        print(f"   {n} = {p} + 2 * {k}^2")
    bad = goldbach_other_counterexamples(10000)
    print(f"   odd composites up to 10000 that fail: {bad}")
    print(f"   5777 = 53 * 109 is composite: {53 * 109 == 5777}")
    print()
    print("3. A program can run without errors and still be wrong")
    print(f"   sum_loop(100) = {sum_loop(100)},"
          f" sum_loop_buggy(100) = {sum_loop_buggy(100)}")
    print("   The bug: range(1, n) stops at n - 1. No error message appears.")
    print()
    print("Epistemic status: part 1 is experimental evidence for 1000 cases;")
    print("the lesson proves both formulas for every n. Part 2 checks every odd")
    print("composite up to 10000: the counterexamples refute Goldbach's claim")
    print("for good, but the search says nothing about numbers beyond 10000.")

if __name__ == "__main__":
    main()

Runs entirely in your browser, in an isolated worker. Nothing is sent to a server.

The thousand successes in part 1 are evidence; Theorem 1 and the Lean proofs are what make the formulas certain. The two failures in part 2, $5777$ and $5993$ , settle Goldbach's question for good — they were found by hand in 1856 by Moritz Stern and his students, and a laptop now finds them instantly.

Binary by repeated division, and why 0.1 + 0.2 != 0.3PythonExact computationExhaustive finite checkComputationally tested

"""Binary numbers by repeated division, and why 0.1 + 0.2 != 0.3.

The experiment converts a number to binary by dividing by 2 again and again
and reading the remainders, compares the result with Python's bin() for
every number below 1024, and checks that the 1024 binary strings without
leading zeros of length at most 10 name each of 0..1023 exactly once. Then
it shows that a computer stores 0.1 as a nearby binary fraction, so
0.1 + 0.2 == 0.3 is False, while exact fractions get it right. Exact and
exhaustive-finite evidence.
"""

from fractions import Fraction

def division_steps(n):
    """The steps n = 2*q + r of repeated division, as (n, q, r) triples."""
    steps = []
    while n > 0:
        steps.append((n, n // 2, n % 2))
        n = n // 2
    return steps

def to_binary(n):
    """Binary digits of n as a string, most significant first ('0' for 0)."""
    if n == 0:
        return "0"
    digits = ""
    while n > 0:
        digits = str(n % 2) + digits
        n = n // 2
    return digits

def from_binary(digits):
    """The value of a string of binary digits."""
    value = 0
    for d in digits:
        value = 2 * value + int(d)
    return value

def strings_without_leading_zeros(max_length):
    """'0' and every binary string of length <= max_length starting with 1."""
    found = ["0"]
    for length in range(1, max_length + 1):
        for k in range(2 ** (length - 1), 2 ** length):
            found.append(format(k, "b"))
    return found

def binary_fraction_digits(x, count):
    """The first `count` binary digits after the point of 0 <= x < 1."""
    digits = ""
    for _ in range(count):
        x = 2 * x
        digit = int(x)
        digits += str(digit)
        x -= digit
    return digits

def main():
    print("1. 45 in binary, by repeated division by 2")
    for n, q, r in division_steps(45):
        print(f"   {n:>2} = 2 * {q:>2} + {r}")
    print(f"   remainders, read upwards: {to_binary(45)}"
          f"  (check: {from_binary(to_binary(45))})")
    print()
    print("2. Every number below 1024")
    same = all(to_binary(n) == format(n, "b") for n in range(1024))
    back = all(from_binary(to_binary(n)) == n for n in range(1024))
    names = strings_without_leading_zeros(10)
    values = sorted(from_binary(s) for s in names)
    print("   repeated division agrees with Python's bin():", same)
    print("   converting back gives the number again:     ", back)
    print(f"   {len(names)} strings without leading zeros name each of 0..1023"
          f" once: {values == list(range(1024))}")
    print()
    print("3. Decimals that are not exact in binary")
    print("   0.1 + 0.2 =", 0.1 + 0.2, "  so 0.1 + 0.2 == 0.3 is", 0.1 + 0.2 == 0.3)
    print("   1/10 in binary: 0." + binary_fraction_digits(Fraction(1, 10), 24)
          + "...")
    print("   the stored value of 0.1 is exactly", Fraction(0.1))
    tenths = 0.0
    for _ in range(10):
        tenths += 0.1
    print("   adding 0.1 ten times, one step at a time, gives", tenths)
    exact = Fraction(1, 10) + Fraction(2, 10) == Fraction(3, 10)
    print("   with exact fractions, 1/10 + 2/10 == 3/10 is", exact)
    print()
    print("Epistemic status: parts 1 and 2 are exact and exhaustive for the")
    print("numbers below 1024; the lesson proves existence and uniqueness of")
    print("binary representations for every whole number. Part 3 shows the")
    print("actual stored values: floating-point numbers are approximations.")

if __name__ == "__main__":
    main()

Runs entirely in your browser, in an isolated worker. Nothing is sent to a server.

Explore

Explore: an algorithm you can watch

The Goldbach experiment needs to know which numbers are prime, and it finds out with one of the oldest algorithms in the world, the Sieve of Eratosthenes: cross out $1$ , then the multiples of each prime in turn; the numbers that survive are the primes. Step through it here.

The Sieve of Eratosthenes

The same steps as the prime_sieve function in the experiment above.

prime
crossed out

Step 1 of 7

Here are the numbers from 1 to 120. Press “Next step” to start sieving.

Things to notice: why is it enough to sieve with the primes up to $10$ when the grid goes to $120$ ? (A composite number up to $120$ must have a factor at most $\sqrt{120} < 11$ .) In the experiment's prime_sieve, the line for multiple in range(p * p, limit + 1, p) starts crossing out at $p^2$ for the same reason you will see in the widget: smaller multiples of $p$ have already been crossed out.

Exercises

1. Programming for mathematics

Warm-upComputational

What is the value of x after these three lines run?

x = 3
x = x + 4
x = 2 * x

2. Programming for mathematics

Warm-upComputational

What does the Python line count = count + 1 do?

Choose one

It increases the value stored in count by

1

Nothing, because the equation

c = c + 1

has no solutionIt causes an error, because count appears on both sidesIt tests whether count equals count + 1

3. Programming for mathematics

Warm-upHuman

Write the binary number $110101_2$ in decimal.

4. Programming for mathematics

CoreHuman

Write $26$ in binary, without leading zeros. (Type just the digits, for example 101.)

5. Programming for mathematics

CoreComputational

What does this program print?

total = 0
for i in range(1, 6):
    if i % 2 == 1:
        total = total + i
print(total)

6. Programming for mathematics

CoreComputational

In Python, 0.1 + 0.2 == 0.3 evaluates to False. Why?

Choose one

0.1

and

0.2

cannot be stored exactly in binary, so the computer adds two slightly wrong numbersIt is a bug in Python that will be fixed in a later versionPython rounds every decimal to a whole numberPython adds from right to left, so the order is wrong

7. Programming for mathematics

StretchHuman

The uniqueness proof for binary representations starts: “the last binary digit of $n$ must be the remainder of $n$ on division by $2$ ”. Why is that true?

Choose one

Every other digit stands for an even amount (

2, 4, 8, \dots

), so only the last digit decides whether

n

is even or oddBecause we define the last digit to be the remainderBecause the last digit is the smallest digitBecause binary numbers always end in

1

8. Programming for mathematics

CoreComputational

Write to_binary(n) that returns the binary representation of a whole number n as a string without leading zeros, using repeated division by 2 (don't use Python's built-in bin). to_binary(0) should be "0".

Runs entirely in your browser, in an isolated worker. Nothing is sent to a server.

9. Programming for mathematics

StretchComputational

This function is meant to compute $1^2 + 2^2 + \dots + n^2$ . It runs without any error — but it is wrong. Find and fix the bug. Then check your fixed version against the formula $\tfrac{n(n + 1)(2n + 1)}{6}$ .

Runs entirely in your browser, in an isolated worker. Nothing is sent to a server.

10. Programming for mathematics

CoreFormal

In the Lean section, which of these machine-checked statements tells you that the sum $1 + 2 + \dots + 1000$ equals $\tfrac{1000 \cdot 1001}{2}$ ?

Choose one

gauss_formula : 2 * sumTo n = n * (n + 1) (proved by induction)loop_agrees_tested : ∀ n, n < 40 → sumFor n = sumTo n (proved by running both programs)sum_for_examples : sumFor 10 = 55 ∧ sumFor 100 = 5050sum_odds_eq_square : sumOdds n = n * n

Common misconceptionsH

Common misconceptions

Misconception

Reading = in a program as mathematical equality. x = x + 1 is not a false equation; it is an instruction that increases x by $1$ . In Python, equality is tested with ==, and x == x + 1 is always False.

Misconception

“The program ran without errors, so the answer is right.” A program does exactly what you wrote, not what you meant. range(1, n) stops at $n - 1$ ; a wrong formula is computed just as happily as a right one. Test programs on cases you can check by hand, and compare with what you can prove.

Misconception

“Computer arithmetic with decimals is exact.” Floating-point numbers are binary approximations: 0.1 + 0.2 gives 0.30000000000000004, and adding 0.1 ten times, one step at a time, gives 0.9999999999999999. Whole numbers and Fractions are exact; floats are not.

Misconception

“A conjecture that holds for a million cases is proved.” Goldbach's conjecture about odd numbers held for every case below $5777$ . A computer check is evidence, sometimes very strong evidence, but only a proof covers every case.

Misconception

“Binary is a different kind of number.” $101101_2$ and $45$ are the same number written in different bases, just as “forty-five” and “quarante-cinq” are the same number in different languages. Only the representation changes.

ViewpointsH

Other ways of seeing

Number theory— Counterexamples found by computer

In 1769 Euler conjectured that a fifth power cannot be the sum of only four fifth powers. In 1966, L. J. Lander and T. R. Parkin found by a computer search that

27^5 + 84^5 + 110^5 + 133^5 = 144^5,

refuting the conjecture in a two-sentence paper. Other famous questions have resisted every search. The Collatz problem — start with any whole number; halve it if it is even, otherwise triple it and add $1$ ; repeat — seems always to reach $1$ , has been checked by computer for an astronomical number of starting values, and is still unproved. Programs are the telescopes of number theory: they show what is out there, and proofs explain it (see experimental mathematics).

Logic— Programs and proofs

When can a computer settle a question completely? If there are only finitely many cases, checking all of them is a proof — that is how Lean's decide works, and how the Four Colour Theorem was proved with computer help by Kenneth Appel and Wolfgang Haken in 1976 (by reducing the problem to a long but finite list of cases). For infinitely many cases, a program can only test; then we need a proof, and systems like Lean check the proof itself. In the 1930s Alan Turing showed that some questions about programs cannot be answered by any program at all (see algorithms and efficiency).

Analysis— Living with approximation

Floating-point errors are tiny — about one part in $10^{16}$ for each operation — but they can grow. Subtracting two nearly equal numbers loses most of the correct digits, and a long calculation can pile up error upon error. Numerical analysis studies how big the error of a computation can be, and how to arrange calculations to keep it small (see numerical analysis). Its basic idea is the analyst's: an answer that is not exact can still be trusted, if we can prove a bound on how far it is from the truth.

History & global context

History and global context

The word algorithm comes from the name of Muḥammad ibn Mūsā al-Khwārizmī, who worked in Baghdad around 820 CE; his book on calculating with Hindu–Arabic numerals, translated into Latin, spread decimal place value to Europe.
Piṅgala's ancient Indian treatise on Sanskrit poetic metre (its date is uncertain, often put in the last centuries BCE) lists the patterns of short and long syllables in a systematic order that amounts to counting in binary.
Gottfried Wilhelm Leibniz published his account of binary arithmetic in 1703, and connected it with the hexagrams of the Chinese Yijing (I Ching), which he learned about from the missionary Joachim Bouvet.
In 1843 Ada Lovelace published notes on Charles Babbage's planned Analytical Engine, including a step-by-step method for computing Bernoulli numbers — often called the first published computer program.
Grace Hopper led the development of some of the first compilers in the 1950s, programs that translate human-readable instructions into machine code. Guido van Rossum released the first version of Python in 1991.

A thousand cases in the blink of an eye

Programs are precise recipes

Definitions

Worked programs

Testing is not proving

Binary numbers

Why 0.1 + 0.2 is not 0.3

Programs that come with proofs

Experiments

Explore: an algorithm you can watch

Exercises

Common misconceptions

Other ways of seeing

History and global context

Further reading