Four-million on, four-million off: even Fibonacci numbers

Project Euler #2 is a five-line problem with a forty-line worth of structure hiding behind it. The official answer fits in a u32, but the path from "loop the sequence" to "evaluate one Fibonacci number" passes through a parity argument, a linear recurrence in the even subsequence, and a closed-form summation identity — all four worth writing down.

Why Python

The problem's natural shape is infinite stream → take while bounded → sum the matches, and Python's generators map onto that idea more directly than anything compiled. itertools.takewhile, generator expressions, and sum together read like the spec; timeit is in the standard library, so I can show the benchmarks in the same file. Python's bignum integers also let me run the same code on a pathological 10^1000 limit without changing a line — useful when the post wants to compare asymptotic behaviour, not just the four-million case. None of that is a speed argument: for this problem the bound is so small that any of C, Rust, Go, or Python finishes inside a microsecond. The argument is legibility per line.

The one thing Python is bad at — recursion overhead — does bite the closed-form solver, and I'll be honest about that when the benchmark numbers come in.

The problem, restated

Project Euler problem #2 ("Even Fibonacci Numbers") asks for the sum of the even-valued terms of a Fibonacci-style sequence whose values do not exceed four million. Project Euler seeds the sequence with 1, 2 rather than the textbook 1, 1, so the first few terms are

1, 2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233, 377, 610, 987, ...

The seed shift is just an index relabel: 1, 2, 3, 5, 8, ... is the standard Fibonacci sequence with the first 1 chopped off, i.e. F_2, F_3, F_4, ... in the convention where F_1 = F_2 = 1. Both conventions produce the same set of values, so the answer is the same; only the indexing changes. I'll keep two parallel notations in this post: PE-indexing where f_1 = 1, f_2 = 2, and standard indexing where F_1 = F_2 = 1, F_3 = 2. Whenever I need a Fibonacci identity I'll use the standard convention because every reference uses it.

Find the sum

Project Euler's etiquette is that we don't post the numerical answer, so the post stops at the algorithm and the runtime; the script at the bottom tells you the number.

First instinct: just sweep

The most direct approach is exactly what the problem statement says aloud. Generate Fibonacci numbers, and for each one, check if it's even and within the limit; if both, add it to a running total.

from itertools import takewhile

def fibs_pe():
    a, b = 1, 2
    while True:
        yield a
        a, b = b, a + b

def solve_naive(limit=4_000_000):
    return sum(x for x in takewhile(lambda v: v <= limit, fibs_pe()) if x % 2 == 0)

fibs_pe() is an infinite generator. takewhile cuts it off at the limit. The if x % 2 == 0 filter selects even values. sum reduces.

This is fine — it's O(n) where n is the number of Fibonacci terms below the bound, and n is small. Below four million, n = 32. The whole loop runs in single-digit microseconds. There's no reason to make it faster for this bound. There are reasons to make it faster as an exercise, because they are the kind of structural observations that Project Euler exists to teach.

Two of them, in increasing order of cleverness.

A pattern in the parity

Print parity alongside the sequence and the answer falls out:

  i        F_i parity
  1          1    odd
  2          2   even
  3          3    odd
  4          5    odd
  5          8   even
  6         13    odd
  7         21    odd
  8         34   even
  9         55    odd
 10         89    odd
 11        144   even
 ...
 32    3524578   even
...next term 5702887 exceeds 4000000; stop.

The parity pattern is odd, even, odd, odd, even, odd, odd, even, ... — a fixed period of three. Why? Because Fibonacci's recurrence f_{n} = f_{n-1} + f_{n-2} lifts to a recurrence on parities:

Starting from (1, 2) ≡ (1, 0) and iterating that table reproduces O, E, O, O, E, O, O, E, ... forever. The state (prev, curr) mod 2 lives in a four-element space and one of the four states ((0, 0)) is unreachable from any nonzero seed, so the orbit lives in a three-cycle. Period three is forced by linear algebra over GF(2), not by the specific seed.

The first observation is therefore: in PE-indexing, the even terms are f_2, f_5, f_8, f_{11}, ... — every third term, starting at index 2. In standard indexing, they are F_3, F_6, F_9, ... — every third term, starting at index 3. So we can step three at a time and skip two-thirds of the work:

def even_fibs_pe(limit):
    a, b = 1, 2
    while a <= limit:
        if a % 2 == 0:
            yield a
        a, b = b, a + b

def solve_skip(limit=4_000_000):
    return sum(even_fibs_pe(limit))

This still walks every term — I'm only using the parity observation to drop the modulo from sum, not to skip iterations. To actually skip, we need a recurrence that operates on the evens alone.

A recurrence on the evens alone

Let E_k be the k-th even term, E_1 = 2, E_2 = 8, E_3 = 34, E_4 = 144, .... Claim:

Cheap sanity check: 4·8 + 2 = 34, 4·34 + 8 = 144, 4·144 + 34 = 610, 4·610 + 144 = 2584. Holds.

Now the proof, because the proof generalises. Use standard indexing for tidy algebra. We want a recurrence relating F_{n+6}, F_{n+3}, F_n for any n, because that's "every third term." Repeatedly expanding F_{m} = F_{m-1} + F_{m-2}:

F_{n+6} = F_{n+5} + F_{n+4}
        = (F_{n+4} + F_{n+3}) + F_{n+4}
        = 2 F_{n+4} + F_{n+3}
        = 2 (F_{n+3} + F_{n+2}) + F_{n+3}
        = 3 F_{n+3} + 2 F_{n+2}.

I want to eliminate F_{n+2} in favour of F_n and F_{n+3}. From F_{n+2} = F_{n+1} + F_n and F_{n+1} = F_{n+3} - F_{n+2} we get 2 F_{n+2} = F_{n+3} + F_n, so

Setting n = 3(k-1) gives F_{3(k+1)} = 4 F_{3k} + F_{3(k-1)}, i.e. E_{k+1} = 4 E_k + E_{k-1}. The constant 4 is what makes this special: most "every-dth-term" sequences over Fibonacci satisfy a_{k+1} = L_d · a_k + (-1)^{d+1} · a_{k-1} where L_d is the d-th Lucas number. For d = 3, L_3 = 4, and the sign is +1. The same machinery gives a_{k+1} = 7 a_k - a_{k-1} for every fourth term, a_{k+1} = 11 a_k + a_{k-1} for every fifth, and so on.

Translated into code:

def solve_recurrence(limit=4_000_000):
    e_prev, e_curr = 2, 8
    if limit < 2:
        return 0
    total = 2
    while e_curr <= limit:
        total += e_curr
        e_prev, e_curr = e_curr, 4 * e_curr + e_prev
    return total

No parity check, no modulo, no skipped terms — only addition and one multiplication-by-four (a left-shift) per step. The loop runs eleven times for limit = 4·10^6, because there are eleven even terms below four million and the loop trips once per term.

Trace:

  k        E_k      check 4*E_k + E_{k-1}
  1          2
  2          8
  3         34                    4*8 + 2
  4        144                   4*34 + 8
  5        610                 4*144 + 34
  6       2584                4*610 + 144
  7      10946               4*2584 + 610
  8      46368             4*10946 + 2584
  9     196418            4*46368 + 10946
 10     832040           4*196418 + 46368
 11    3524578          4*832040 + 196418
... 4*3524578 + 832040 = 14930352 exceeds 4000000; stop.

Eleven additions, eleven shifts, eleven comparisons. That's the whole computation.

Closed form: turning the sum into one Fibonacci value

There is a third observation, which collapses the loop entirely. Standard convention again. Among the textbook Fibonacci sum identities is this one:

Pretty clean. Quick proof from a more general identity, \sum_{i=1}^{m} F_i = F_{m+2} - 1. Split the standard sum into three residue classes modulo 3:

The first two sums (the odd-indexed contributions in our parity argument) can be related to the third using F_{3k-2} + F_{3k-1} = F_{3k}, which says "two odd terms sum to the next even." So \sum F_{3k-2} + \sum F_{3k-1} = \sum F_{3k} (over matching index ranges), and the total telescopes to 2 \sum F_{3k}, giving the identity above for m = 3n. I verified it numerically for n = 1..11:

  n     F_{3n}   sum so far   (F_{3n+2}-1)/2
  1          2            2                2 OK
  2          8           10               10 OK
  3         34           44               44 OK
  4        144          188              188 OK
  5        610          798              798 OK
  6       2584         3382             3382 OK
  ...

(I'm cutting the table off at n = 6 because the row for n = 11 is the answer to the puzzle.)

So if we know the largest n with F_{3n} ≤ 4·10^6, we can compute the sum by evaluating one Fibonacci number, F_{3n+2}, and dividing by two. Computing one Fibonacci number is O(log n) arithmetic operations using fast doubling:

def fib_pair(n):
    """Return (F_n, F_{n+1}) under standard convention F_1 = F_2 = 1.
    Fast doubling: F_{2k}   = F_k * (2 F_{k+1} - F_k)
                   F_{2k+1} = F_k^2 + F_{k+1}^2
    """
    if n == 0:
        return (0, 1)
    a, b = fib_pair(n >> 1)        # F_k, F_{k+1}
    c = a * (2 * b - a)            # F_{2k}
    d = a * a + b * b              # F_{2k+1}
    if n & 1:
        return (d, c + d)
    return (c, d)

def solve_closed_form(limit=4_000_000):
    n = 0
    while True:
        f3n_next = fib_pair(3 * (n + 1))[0]
        if f3n_next > limit:
            break
        n += 1
    f_3n_plus_2 = fib_pair(3 * n + 2)[0]
    return (f_3n_plus_2 - 1) // 2

Two notes about the implementation:

1. The "find the largest n" loop is the dumb linear version. For limits up to a few hundred digits this is fine; for huge limits you'd binary-search instead. I left it linear because I want to compare the closed form against the recurrence, and the linear search introduces an honest cost that the recurrence doesn't have.

2. fib_pair is recursive in Python, which carries a real per-call overhead. For n ≈ 33 the recursion depth is \log_2(33) ≈ 5, but each call still spawns two integer multiplications. We will pay for this in the benchmarks.

All four solvers in one file

Here is the complete euler2.py that runs every approach and asserts they agree:

"""
Project Euler #2 — sum of even Fibonacci terms whose value does not exceed
4,000,000.  Four solvers; assertion at the bottom checks they agree.
"""

from itertools import takewhile
from typing import Iterator

LIMIT = 4_000_000


# --- 1. The naive sweep -----------------------------------------------------
def fibs_pe() -> Iterator[int]:
    a, b = 1, 2
    while True:
        yield a
        a, b = b, a + b


def solve_naive(limit: int = LIMIT) -> int:
    return sum(x for x in takewhile(lambda v: v <= limit, fibs_pe()) if x % 2 == 0)


# --- 2. Skip via parity-aware filter ---------------------------------------
def even_fibs_pe(limit: int) -> Iterator[int]:
    a, b = 1, 2
    while a <= limit:
        if a % 2 == 0:
            yield a
        a, b = b, a + b


def solve_skip(limit: int = LIMIT) -> int:
    return sum(even_fibs_pe(limit))


# --- 3. The even-only recurrence -------------------------------------------
def solve_recurrence(limit: int = LIMIT) -> int:
    e_prev, e_curr = 2, 8
    if limit < 2:
        return 0
    total = 2
    while e_curr <= limit:
        total += e_curr
        e_prev, e_curr = e_curr, 4 * e_curr + e_prev
    return total


# --- 4. Closed-form via fast doubling --------------------------------------
def fib_pair(n: int) -> tuple[int, int]:
    if n == 0:
        return (0, 1)
    a, b = fib_pair(n >> 1)
    c = a * (2 * b - a)
    d = a * a + b * b
    if n & 1:
        return (d, c + d)
    return (c, d)


def solve_closed_form(limit: int = LIMIT) -> int:
    n = 0
    while True:
        f3n_next = fib_pair(3 * (n + 1))[0]
        if f3n_next > limit:
            break
        n += 1
    f_3n_plus_2 = fib_pair(3 * n + 2)[0]
    return (f_3n_plus_2 - 1) // 2


if __name__ == "__main__":
    a = solve_naive()
    b = solve_skip()
    c = solve_recurrence()
    d = solve_closed_form()
    assert a == b == c == d, (a, b, c, d)
    print("all four agree")
    print("answer =", a)

Running it locally:

$ python3 euler2.py
all four agree
answer = <redacted — readers run their own>

The assertion is the important line: four mutually independent algorithms produce the same number. That's much stronger evidence of correctness than a single solver returning what the OEIS says.

Benchmarks

bench.py runs each solver N times under timeit. Three limits: the actual PE bound (4·10^6), an over-the-horizon 10^18 to push the loop counts up, and a comically large 10^1000 to make the bignum arithmetic actually expensive enough to see asymptotic differences.

import timeit
import sys
sys.path.insert(0, "/labs-output")
from euler2 import solve_naive, solve_skip, solve_recurrence, solve_closed_form

def bench(fn, limit, n=200_000):
    t = timeit.timeit(lambda: fn(limit), number=n)
    return t / n

for limit in (4_000_000, 10**18, 10**1000):
    print(f"\n--- LIMIT = {limit if limit < 10**6 else f'10^{len(str(limit))-1}'} ---")
    n = 100_000 if limit < 10**100 else 200
    for fn in (solve_naive, solve_skip, solve_recurrence, solve_closed_form):
        try:
            t = bench(fn, limit, n=n)
            print(f"{fn.__name__:25s} {t*1e6:10.3f} us per call")
        except Exception as e:
            print(f"{fn.__name__:25s} {type(e).__name__}: {e}")

Output on the sandbox (Python 3.13.12, single core, no warm-up loop):

--- LIMIT = 10^6 ---
solve_naive                    8.547 us per call
solve_skip                     3.755 us per call
solve_recurrence               1.781 us per call
solve_closed_form             17.556 us per call

--- LIMIT = 10^18 ---
solve_naive                   19.570 us per call
solve_skip                    11.719 us per call
solve_recurrence               3.637 us per call
solve_closed_form             52.656 us per call

--- LIMIT = 10^1000 ---
solve_naive                 5045.577 us per call
solve_skip                  4652.642 us per call
solve_recurrence             596.370 us per call
solve_closed_form          16057.720 us per call

The headline result for the actual problem is solve_recurrence at 1.78 microseconds. That's 562,000 evaluations per second on this machine and well, well inside Project Euler's "one-minute rule." A single-shot wall-clock measurement (no timeit loop) puts it at about 2.6 microseconds for the very first call:

$ python3 single_shot.py
solve_naive       : 7.91 us per call (avg of 50000)
solve_recurrence  : 1.18 us per call (avg of 50000)

--- single shot ---
solve_recurrence  : 2650.0 ns single call
solve_naive       : 11946.0 ns single call

Honest finding: the closed-form solver is the slowest of the four at every limit we tested, including the 10^1000 stress test. That looks wrong, because asymptotically O(\log^2 n) Fibonacci should beat the linear recurrence. Two reasons it doesn't, in practice:

1. The "find the largest n" loop calls fib_pair once per increment, so the linear search dominates the logarithmic Fibonacci. Replacing it with a binary search on n would help.
2. Python recursion has a ~hundred-nanosecond floor per call. For n = 33, fast doubling is ~6 recursive calls, roughly a microsecond just in interpreter overhead, before any arithmetic.

If you needed the "even Fibonacci sum below 10^{10000}" version, the right move is a closed-form solver written in a language without per-call recursion overhead (Rust, C with GMP, Julia), with the index-search binary-searched. For the actual problem, solve_recurrence is the right tool.

A worked example, by hand

Suppose someone asks "what is the sum of even Fibonacci values up to 200?" Walk solve_recurrence step by step:

Result: 2 + 8 + 34 + 144 = 188. The check (F_{3·4+2} - 1)/2 = (F_{14} - 1)/2 = (377 - 1)/2 = 188 agrees. We did three additions, three shifts, four comparisons.

For the actual LIMIT = 4_000_000, the loop runs eleven of these rows instead of three, and the final e_curr at the time of exit is 14,930,352 — the next even Fibonacci past four million.

Edge cases I checked

• LIMIT = 0 and LIMIT = 1: the if limit < 2 short-circuit returns 0 in solve_recurrence. The other three return 0 by virtue of producing an empty stream. All agree.
• LIMIT = 2: just the first even term. All four return 2.
• LIMIT = 7: the next even term (8) hasn't appeared yet. All four return 2.
• LIMIT = 8: 8 is included. All four return 10.
• Massive limits: solve_recurrence and solve_naive have no precision issues because Python's ints are arbitrary precision; solve_closed_form is fine for the same reason.
• The fast-doubling recursion is bounded by log_2(3n+2), which for the limits we actually use stays under 20 — well below Python's default recursion limit.

A subtle one I want to call out: Project Euler's "1, 2" seed is not the same set of values as standard Fibonacci. It's the standard sequence with F_1 = 1 deleted. If you write f_1 = 1, f_2 = 1, f_3 = 2, ... and apply the same parity argument, the even terms appear at standard indices 3, 6, 9, 12, .... The values are the same; only the labels change. I picked the "standard" convention everywhere in the math sections because every textbook identity I cited is in that convention. If you read a different writeup that claims "every third term starting from index 2," that author is using PE-indexing — same observation.

When the closed form actually wins

The benchmark above said the closed form loses at every realistic limit. That's only true because of the linear search and the recursion overhead, both fixable. With those fixed (binary-search the index, switch to an iterative fast-doubling, or call a native bignum library), the closed form would beat the recurrence for limits where the recurrence's loop count starts to matter — say, LIMIT = 10^{10000}, where the recurrence does about 50,000 iterations and each iteration moves a 10,000-digit bignum.

For PE #2 specifically, the bound is so small that the constant factors win and the asymptotics never get a chance. Worth knowing: most Project Euler problems are stated at sizes where the clever algorithm and the asymptotically optimal algorithm are the same algorithm, but PE #2 is generous enough that all of them fit. The pedagogical value is in seeing the structure, not in shaving the last microsecond.

What I would do differently in a faster language

If I rewrote the recurrence in C, the only thing that changes is the constant: the eleven additions take a handful of nanoseconds, not microseconds. The structure of the algorithm is identical. There's no SIMD, no parallelism, no memory hierarchy in this problem to exploit; the bottleneck is the eleven serial dependencies in the recurrence. You can't add two evens in parallel without computing the recurrence twice.

Where languages would matter: if the bound were a billion-digit bignum, and you wanted the closed form, you'd want GMP. Python's bignums are fine but not the fastest available; for a 10-million-digit limit, a Rust solver linking rug (which wraps GMP) would be perhaps 5-10x faster, mostly in the multiplication kernel of fast doubling.

Reproducing this

Drop euler2.py and bench.py (both included in full above) into a directory and run:

$ python3 euler2.py
$ python3 bench.py

The first prints the answer. The second prints the benchmark table. No dependencies beyond a vanilla CPython 3.13 install.

References

• Project Euler problem #2: https://projecteuler.net/problem=2
• Fast doubling for Fibonacci: https://www.nayuki.io/page/fast-fibonacci-algorithms
• Lucas-number generalisation of "every-dth-term" recurrences: any standard reference on Fibonacci identities; see e.g. Knuth, The Art of Computer Programming Vol. 1 §1.2.8.
• The summation identity \sum_{i=1}^{m} F_i = F_{m+2} - 1 and its parity-class refinements: derivable from the matrix form [[1,1],[1,0]]^n, but easier to prove by induction on m.

Artefacts

The download bundle includes euler2.py, bench.py, tables.py, closed_form_check.py, and single_shot.py — every script referenced in the post. They run in any CPython 3.10+ install.