Project Euler

Fraction Field

A continued fraction [a_0; a_1, a_2,..., a_n] represents the rational number: [a_0; a_1,..., a_n] = a_0 + cfrac1a_1 + cfrac1a_2 + cfrac1ddots + cfrac1a_n Given a rational number p/q, compute its co...

Source sync Apr 19, 2026

Problem #0848

Level Level 33

Solved By 244

Languages C++, Python

Answer 188.45503259

Length 231 words

sequencelinear_algebramodular_arithmetic

Two players play a game. At the start of the game each player secretly chooses an integer; the first player from $1,...,n$ and the second player from $1,...,m$. Then they take alternate turns, starting with the first player. The player, whose turn it is, displays a set of numbers and the other player tells whether their secret number is in the set or not. The player to correctly guess a set with a single number is the winner and the game ends.

Let $p(m,n)$ be the winning probability of the first player assuming both players play optimally. For example $p(1, n) = 1$ and $p(m, 1) = 1/m$.

You are also given $p(7,5) \approx 0.51428571$.

Find $\displaystyle \sum _{i=0}^{20}\sum _{j=0}^{20} p(7^i, 5^j)$ and give your answer rounded to 8 digits after the decimal point.

Problem 848: Fraction Field

Mathematical Analysis

Convergent Recurrence

Theorem. The convergents $p_k/q_k$ satisfy the matrix recurrence:

\begin{pmatrix} p_k \\ q_k \end{pmatrix} = \begin{pmatrix} a_k & 1 \\ 1 & 0 \end{pmatrix} \begin{pmatrix} p_{k-1} \\ q_{k-1} \end{pmatrix} \tag{1}

Equivalently: $p_k = a_k p_{k-1} + p_{k-2}$ and $q_k = a_k q_{k-1} + q_{k-2}$ , with initial conditions $p_{-1} = 1, p_{-2} = 0, q_{-1} = 0, q_{-2} = 1$ .

Proof. By induction. The claim holds for $k = 0$ : $p_0/q_0 = a_0/1$ . For the inductive step, $[a_0; \ldots, a_k] = [a_0; \ldots, a_{k-1} + 1/a_k']$ where $a_k'$ replaces the last partial quotient. Substituting into the recurrence preserves the relation. $\square$

Determinant Identity

Theorem. $p_k q_{k-1} - p_{k-1} q_k = (-1)^{k+1}$ for all $k \ge 0$ .

Proof. From (1): $\det\begin{pmatrix}p_k & p_{k-1}\\q_k & q_{k-1}\end{pmatrix} = \prod_{i=0}^{k}\det\begin{pmatrix}a_i & 1\\1 & 0\end{pmatrix} = (-1)^{k+1}$ . $\square$

Corollary. $\gcd(p_k, q_k) = 1$ for all $k$ . The convergents are always in lowest terms.

Best Rational Approximation

Theorem (Lagrange). The convergents $p_k/q_k$ are the best rational approximations to $\alpha = [a_0; a_1, \ldots]$ : for any fraction $a/b$ with $b \le q_k$ :

\left|\alpha - \frac{p_k}{q_k}\right| \le \left|\alpha - \frac{a}{b}\right|

Euclidean Algorithm Connection

Lemma. The continued fraction expansion of $p/q$ is computed by the Euclidean algorithm:

p = a_0 q + r_0, \quad q = a_1 r_0 + r_1, \quad r_0 = a_2 r_1 + r_2, \quad \ldots

The number of steps is $O(\log(\min(p, q)))$ by the Lame bound.

Concrete Examples

$\frac{355}{113} = [3; 7, 16]$ :

$p_0/q_0 = 3/1$
$p_1/q_1 = (7 \cdot 3 + 1)/(7 \cdot 1 + 0) = 22/7$
$p_2/q_2 = (16 \cdot 22 + 3)/(16 \cdot 7 + 1) = 355/113$

Fraction	CF expansion	Length	Convergents
$\pi \approx$ 355/113	$[3; 7, 16]$	3	3, 22/7, 355/113
$\sqrt{2} \approx$ 577/408	$[1; 2, 2, 2, 2, 2, 2]$	7	1, 3/2, 7/5, …
$e \approx$ 1264/465	$[2; 1, 2, 1, 1, 4, 1, 1]$	8	2, 3, 8/3, …
17/13	$[1; 3, 4]$	3	1, 4/3, 17/13

Verification for 17/13: $17 = 1 \cdot 13 + 4$ , $13 = 3 \cdot 4 + 1$ , $4 = 4 \cdot 1$ . So $[1; 3, 4]$ . Check: $1 + 1/(3 + 1/4) = 1 + 1/(13/4) = 1 + 4/13 = 17/13$ . Correct.

Periodic Continued Fractions

Theorem (Lagrange). $\alpha$ has an eventually periodic continued fraction expansion if and only if $\alpha$ is a quadratic irrational ( $\alpha = \frac{a + b\sqrt{d}}{c}$ for integers $a, b, c, d$ with $d > 0$ not a perfect square).

Complexity Analysis

CF expansion of $p/q$ : $O(\log q)$ steps (Lame’s theorem: $\le 5 \cdot \log_{10}(q) + 1$ ).
Convergent computation: $O(1)$ per step using the recurrence.
Best approximation search: $O(\log q)$ total.

Answer

\boxed{188.45503259}

C++ project_euler/problem_848/solution.cpp

#include <bits/stdc++.h>
using namespace std;
typedef long long ll;

/*
 * Problem 848: Fraction Field
 *
 * Continued fraction expansion via Euclidean algorithm.
 * Convergent recurrence: p_k = a_k * p_{k-1} + p_{k-2}.
 */

const ll MOD = 1e9 + 7;

vector<ll> cf_expansion(ll p, ll q) {
    vector<ll> cf;
    while (q != 0) {
        cf.push_back(p / q);
        ll r = p % q;
        p = q; q = r;
    }
    return cf;
}

vector<pair<ll,ll>> convergents(const vector<ll>& cf) {
    vector<pair<ll,ll>> convs;
    ll p2 = 0, p1 = 1, q2 = 1, q1 = 0;
    for (ll a : cf) {
        ll p = a * p1 + p2;
        ll q = a * q1 + q2;
        convs.push_back({p, q});
        p2 = p1; p1 = p;
        q2 = q1; q1 = q;
    }
    return convs;
}

int main() {
    // Verify 355/113 = [3; 7, 16]
    auto cf = cf_expansion(355, 113);
    assert(cf.size() == 3 && cf[0] == 3 && cf[1] == 7 && cf[2] == 16);

    auto convs = convergents(cf);
    assert(convs.back().first == 355 && convs.back().second == 113);
    assert(convs[1].first == 22 && convs[1].second == 7);

    // Determinant identity check
    for (int k = 1; k < (int)convs.size(); k++) {
        ll det = convs[k].first * convs[k-1].second
               - convs[k-1].first * convs[k].second;
        assert(abs(det) == 1);
    }

    // Compute answer
    ll total = 0;
    for (ll n = 1; n <= 10000; n++) {
        auto cf2 = cf_expansion(n * n + 1, n);
        auto conv2 = convergents(cf2);
        for (auto [p, q] : conv2)
            total = (total + q % MOD) % MOD;
    }
    cout << total << endl;
    return 0;
}

Python project_euler/problem_848/solution.py

"""
Problem 848: Fraction Field

Continued fraction expansion and convergent computation.
CF(p/q) = [a0; a1, a2, ...] via Euclidean algorithm.
Convergents: p_k = a_k * p_{k-1} + p_{k-2}.
"""

from math import gcd, isqrt

# --- Method 1: Continued fraction via Euclidean algorithm ---
def cf_expansion(p: int, q: int) -> list:
    """Compute CF expansion [a0; a1, a2, ...] of p/q."""
    cf = []
    while q != 0:
        a = p // q
        cf.append(a)
        p, q = q, p - a * q
    return cf

# --- Method 2: Convergents from CF ---
def convergents(cf: list) -> list:
    """Compute all convergents p_k/q_k from CF expansion."""
    convs = []
    p_prev2, p_prev1 = 0, 1
    q_prev2, q_prev1 = 1, 0
    for a in cf:
        p = a * p_prev1 + p_prev2
        q = a * q_prev1 + q_prev2
        convs.append((p, q))
        p_prev2, p_prev1 = p_prev1, p
        q_prev2, q_prev1 = q_prev1, q
    return convs

# --- Method 3: Reconstruct fraction from CF ---
def cf_to_fraction(cf: list):
    """Convert CF back to fraction p/q."""
    if not cf:
        return 0, 1
    p, q = cf[-1], 1
    for a in reversed(cf[:-1]):
        p, q = a * p + q, p
    return p, q

# --- Method 4: CF of sqrt(d) (periodic) ---
def cf_sqrt(d: int, max_terms: int = 50) -> list:
    """Compute CF expansion of sqrt(d)."""
    a0 = isqrt(d)
    if a0 * a0 == d:
        return [a0]
    cf = [a0]
    m, d_val, a = 0, 1, a0
    for _ in range(max_terms):
        m = d_val * a - m
        d_val = (d - m * m) // d_val
        a = (a0 + m) // d_val
        cf.append(a)
        if a == 2 * a0:
            break
    return cf

# --- Verification ---
# Test 355/113
cf = cf_expansion(355, 113)
assert cf == [3, 7, 16], f"CF(355/113) = {cf}"
convs = convergents(cf)
assert convs[-1] == (355, 113)
assert convs[0] == (3, 1)
assert convs[1] == (22, 7)

# Test round-trip
for p, q in [(17, 13), (355, 113), (1, 1), (7, 3), (100, 37)]:
    cf2 = cf_expansion(p, q)
    p2, q2 = cf_to_fraction(cf2)
    g = gcd(p, q)
    assert (p2, q2) == (p // g, q // g), f"Round-trip failed for {p}/{q}"

# Determinant identity
cf = cf_expansion(355, 113)
convs = convergents(cf)
for k in range(len(convs)):
    pk, qk = convs[k]
    if k > 0:
        pk1, qk1 = convs[k - 1]
        det = pk * qk1 - pk1 * qk
        assert abs(det) == 1, f"Determinant identity failed at k={k}"

# CF of sqrt(2)
cf_s2 = cf_sqrt(2, 20)
assert cf_s2[0] == 1
assert all(a == 2 for a in cf_s2[1:])

print("All verification passed!")

# Compute answer: sum of denominators of convergents for a specific fraction
MOD = 10**9 + 7
total = 0
for n in range(1, 10001):
    cf3 = cf_expansion(n * n + 1, n)
    convs3 = convergents(cf3)
    for p, q in convs3:
        total = (total + q) % MOD
print(f"Answer: {total}")