Problem 64: Odd Period Square Roots

Project Euler Problem 64 Solution

EulerSolve provides an optimized solution for Project Euler Problem 64, Odd Period Square Roots, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary For each integer $N$ with $2 \le N \le 10{,}000$, we ignore perfect squares and study the simple continued fraction of $\sqrt{N}$. For every non-square $N$, that expansion has the form $\sqrt{N}=[a_0;\overline{a_1,a_2,\dots,a_r}]$, and the task is to count how many values of $N$ have odd period length $r$. The small example up to $N=13$ already shows what is being counted: the odd periods occur for $N=2,5,10,13$, so the count there is $4$. The full problem asks for the same test over the whole range up to $10{,}000$. Mathematical Approach The decisive observation is that the period is encoded by a very small integer state. Once the continued fraction of $\sqrt{N}$ is written through its complete quotients, the entire problem reduces to iterating an exact recurrence on integers and checking whether the resulting period length is odd. Complete Quotients and the Initial State Let $a_0=\lfloor \sqrt{N}\rfloor$. If $a_0^2=N$, then $\sqrt{N}$ is an integer and there is no periodic part at all, so such values are skipped immediately....

Detailed mathematical approach

Problem Summary

For each integer $N$ with $2 \le N \le 10{,}000$, we ignore perfect squares and study the simple continued fraction of $\sqrt{N}$. For every non-square $N$, that expansion has the form $\sqrt{N}=[a_0;\overline{a_1,a_2,\dots,a_r}]$, and the task is to count how many values of $N$ have odd period length $r$.

The small example up to $N=13$ already shows what is being counted: the odd periods occur for $N=2,5,10,13$, so the count there is $4$. The full problem asks for the same test over the whole range up to $10{,}000$.

Mathematical Approach

The decisive observation is that the period is encoded by a very small integer state. Once the continued fraction of $\sqrt{N}$ is written through its complete quotients, the entire problem reduces to iterating an exact recurrence on integers and checking whether the resulting period length is odd.

Complete Quotients and the Initial State

Let $a_0=\lfloor \sqrt{N}\rfloor$. If $a_0^2=N$, then $\sqrt{N}$ is an integer and there is no periodic part at all, so such values are skipped immediately. For a non-square $N$, we begin with

$$\alpha_0=\sqrt{N}=\frac{\sqrt{N}+0}{1}.$$

At every stage we write the current complete quotient in the standard quadratic-irrational form

$$\alpha_k=\frac{\sqrt{N}+m_k}{d_k},\qquad a_k=\lfloor \alpha_k \rfloor,$$

with initial values $m_0=0$, $d_0=1$, and $a_0=\lfloor \sqrt{N}\rfloor$. A continued-fraction step removes the integer part $a_k$ and inverts the remaining fractional part.

Deriving the Integer Recurrence

By definition of a simple continued fraction, the next complete quotient is

$$\alpha_{k+1}=\frac{1}{\alpha_k-a_k}.$$

Substituting $\alpha_k=(\sqrt{N}+m_k)/d_k$ and rationalizing the denominator gives

$$\alpha_{k+1}=\frac{d_k}{\sqrt{N}+m_k-a_kd_k} =\frac{d_k(\sqrt{N}+a_kd_k-m_k)}{N-(a_kd_k-m_k)^2}.$$

This suggests the standard definitions

$$m_{k+1}=d_ka_k-m_k,\qquad d_{k+1}=\frac{N-m_{k+1}^2}{d_k},$$

so that the next state again has the same form,

$$\alpha_{k+1}=\frac{\sqrt{N}+m_{k+1}}{d_{k+1}}.$$

Its integer part is then

$$a_{k+1}=\left\lfloor \frac{\sqrt{N}+m_{k+1}}{d_{k+1}} \right\rfloor =\left\lfloor \frac{a_0+m_{k+1}}{d_{k+1}} \right\rfloor,$$

because $a_0 < \sqrt{N} < a_0+1$. These are exactly the update formulas used by the implementations.

Why the Divisions Stay Exact

The recurrence never needs floating-point arithmetic after the initial floor. If $d_k\mid (N-m_k^2)$, then from $m_{k+1}=d_ka_k-m_k$ we have

$$m_{k+1}\equiv -m_k \pmod{d_k},\qquad m_{k+1}^2\equiv m_k^2 \pmod{d_k}.$$

Therefore $d_k\mid (N-m_{k+1}^2)$, so

$$d_{k+1}=\frac{N-m_{k+1}^2}{d_k}$$

is an integer at every step. Since $d_0=1$, exact divisibility holds inductively for the entire loop. This is why the C++, Python, and Java implementations can work purely with integers.

Why $a_k = 2a_0$ Marks the End of the Period

For square roots of non-squares, the continued fraction is periodic, and the closing state of one full period is the distinguished reduced surd

$$m_k=a_0,\qquad d_k=1.$$

If this state occurs, then

$$a_k=\left\lfloor \sqrt{N}+a_0 \right\rfloor=2a_0,$$

so reaching $(m_k,d_k)=(a_0,1)$ certainly produces $a_k=2a_0$.

The converse is just as important. In the reduced square-root process one has $0 \le m_k < \sqrt{N} < a_0+1$, so if $a_k=2a_0$, then

$$2a_0 \le \frac{a_0+m_k}{d_k} < \frac{2a_0+1}{d_k}.$$

The left inequality is impossible unless $d_k=1$. Once $d_k=1$, the formula becomes $2a_0=\lfloor a_0+m_k\rfloor$, and since $m_k$ is an integer with $m_k<a_0+1$, it follows that $m_k=a_0$. Thus the first appearance of $a_k=2a_0$ is exactly the first return to the closing state, which means one full period has been completed.

Worked Example: $\sqrt{23}$

Take $N=23$. Then $a_0=\lfloor \sqrt{23}\rfloor=4$. Iterating the recurrence produces

$$ (m_1,d_1,a_1)=(4,7,1),\quad (m_2,d_2,a_2)=(3,2,3),\quad (m_3,d_3,a_3)=(3,7,1),\quad (m_4,d_4,a_4)=(4,1,8). $$

Because $a_4=8=2a_0$, the period length is $4$, so

$$\sqrt{23}=[4;\overline{1,3,1,8}].$$

For contrast, $\sqrt{13}=[3;\overline{1,1,1,1,6}]$ reaches $2a_0=6$ after $5$ steps, so it contributes to the odd-period count.

From Period Parity to the Final Count

Once the period length $r(N)$ is known for one non-square $N$, the original problem becomes the count

$$\#\{\,N: 2 \le N \le 10{,}000,\ N \notin \{m^2 : m \in \mathbb{Z}\},\ r(N)\equiv 1 \pmod 2\,\}.$$

Each value of $N$ is completely independent. There is no interaction between different square roots, so the global answer is obtained by running the same period-finding routine once for every candidate $N$.

How the Code Works

Computing One Period

The implementation first computes $a_0=\lfloor \sqrt{N}\rfloor$ with an integer square-root operation. If $a_0^2=N$, it returns period length $0$ immediately. Otherwise it initializes the state with $m=0$, $d=1$, $a=a_0$, then repeatedly applies the three recurrence updates and increments a period counter.

The stopping rule is deliberately minimal: the loop ends at the first $a=2a_0$. Because that condition is mathematically equivalent to the closing state $(m,d)=(a_0,1)$, the implementation does not need a hash set or any other structure to remember earlier states.

Sweeping the Whole Interval

After period computation is available for a single $N$, the outer routine scans all integers from $2$ to the chosen limit, tests whether the returned period is odd, and increments the answer if it is. The default limit is $10{,}000$.

The C++, Python, and Java implementations also include a small built-in checkpoint from the statement: up to $13$, the correct count of odd periods is $4$. That single test validates the recurrence, the stopping condition, and the parity check together.

Complexity Analysis

For a fixed $N$, each iteration of the continued-fraction loop performs only $O(1)$ integer work, so the cost is proportional to the period length $r(N)$. Classical theory gives the worst-case upper bound $r(N)=O(\sqrt{N})$, which leads to

$$\sum_{N=2}^{L} O(\sqrt{N}) = O(L^{3/2})$$

for an overall limit $L$. With $L=10{,}000$, this is easily small enough for a fast computation.

Memory usage is $O(1)$. The algorithm stores only the current state $(m,d,a)$, the initial floor $a_0$, the current period length, and the global counter of odd periods.

Footnotes and References

Problem page: https://projecteuler.net/problem=64
Continued fractions: Wikipedia - Continued fraction
Periodic continued fractions: Wikipedia - Periodic continued fractions
Quadratic irrationals: Wikipedia - Quadratic irrational
Pell's equation: Wikipedia - Pell's equation

Problem 64 source code

C++

#include <cmath>
#include <iostream>
#include <string>

namespace {

struct Options {
    int limit = 10000;
    bool run_checkpoints = true;
};

bool parse_int_after_prefix(const std::string& arg, const std::string& prefix, int& value) {
    if (arg.rfind(prefix, 0U) != 0U) {
        return false;
    }
    const std::string tail = arg.substr(prefix.size());
    if (tail.empty()) {
        return false;
    }

    int parsed = 0;
    for (char c : tail) {
        if (c < '0' || c > '9') {
            return false;
        }
        parsed = parsed * 10 + static_cast<int>(c - '0');
    }
    value = parsed;
    return true;
}

bool parse_arguments(int argc, char** argv, Options& options) {
    for (int i = 1; i < argc; ++i) {
        const std::string arg(argv[i]);

        if (arg == "--skip-checkpoints") {
            options.run_checkpoints = false;
            continue;
        }
        if (parse_int_after_prefix(arg, "--limit=", options.limit)) {
            continue;
        }

        std::cerr << "Unknown argument: " << arg << '\n';
        return false;
    }
    return options.limit >= 2;
}

int period_length(const int n) {
    const int a0 = static_cast<int>(std::sqrt(n));
    if (a0 * a0 == n) {
        return 0;
    }

    int m = 0;
    int d = 1;
    int a = a0;
    int period = 0;

    do {
        m = d * a - m;
        d = (n - m * m) / d;
        a = (a0 + m) / d;
        ++period;
    } while (a != 2 * a0);

    return period;
}

int solve(const int limit) {
    int count = 0;
    for (int n = 2; n <= limit; ++n) {
        if ((period_length(n) & 1) == 1) {
            ++count;
        }
    }
    return count;
}

bool run_checkpoints() {
    if (solve(13) != 4) {
        std::cerr << "Checkpoint failed for limit=13" << '\n';
        return false;
    }
    return true;
}

}  // namespace

int main(int argc, char** argv) {
    Options options;
    if (!parse_arguments(argc, argv, options)) {
        return 1;
    }
    if (options.run_checkpoints && !run_checkpoints()) {
        return 2;
    }

    std::cout << solve(options.limit) << '\n';
    return 0;
}

Python

import math
import sys

def parse_int_after_prefix(arg, prefix):
    if not arg.startswith(prefix):
        return False, 0
    tail = arg[len(prefix):]
    if not tail:
        return False, 0
    
    try:
        parsed = int(tail)
        return True, parsed
    except ValueError:
        return False, 0

def period_length(n):
    a0 = int(math.isqrt(n))
    if a0 * a0 == n:
        return 0
    
    m = 0
    d = 1
    a = a0
    period = 0
    
    while True:
        m = d * a - m
        d = (n - m * m) // d
        a = (a0 + m) // d
        period += 1
        if a == 2 * a0:
            break
    
    return period

def solve(limit):
    count = 0
    for n in range(2, limit + 1):
        if period_length(n) % 2 == 1:
            count += 1
    return count

def main():
    limit = 10000
    run_checkpoints = True
    
    args = sys.argv[1:]
    for arg in args:
        if arg == "--skip-checkpoints":
            run_checkpoints = False
        else:
            success, value = parse_int_after_prefix(arg, "--limit=")
            if success:
                limit = value
            else:
                print(f"Unknown argument: {arg}", file=sys.stderr)
                sys.exit(1)
    
    if limit < 2:
        print("Limit must be at least 2", file=sys.stderr)
        sys.exit(1)
    
    if run_checkpoints:
        if solve(13) != 4:
            print("Checkpoint failed for limit=13", file=sys.stderr)
            sys.exit(2)
    
    print(solve(limit))

if __name__ == "__main__":
    main()

Java

class Euler64 {

    private static class Options {
        int limit = 10000;
        boolean runCheckpoints = true;
    }

    private static boolean parseIntAfterPrefix(String arg, String prefix, Options options) {
        if (!arg.startsWith(prefix)) {
            return false;
        }
        String tail = arg.substring(prefix.length());
        if (tail.isEmpty()) {
            return false;
        }

        int parsed = 0;
        for (int i = 0; i < tail.length(); i++) {
            char c = tail.charAt(i);
            if (c < '0' || c > '9') {
                return false;
            }
            parsed = parsed * 10 + (c - '0');
        }
        options.limit = parsed;
        return true;
    }

    private static boolean parseArguments(String[] args, Options options) {
        for (String arg : args) {
            if ("--skip-checkpoints".equals(arg)) {
                options.runCheckpoints = false;
                continue;
            }
            if (parseIntAfterPrefix(arg, "--limit=", options)) {
                continue;
            }

            System.err.println("Unknown argument: " + arg);
            return false;
        }
        return options.limit >= 2;
    }

    private static int periodLength(int n) {
        int a0 = (int) Math.sqrt(n);
        if (a0 * a0 == n) {
            return 0;
        }

        int m = 0;
        int d = 1;
        int a = a0;
        int period = 0;

        do {
            m = d * a - m;
            d = (n - m * m) / d;
            a = (a0 + m) / d;
            period++;
        } while (a != 2 * a0);

        return period;
    }

    private static int solve(int limit) {
        int count = 0;
        for (int n = 2; n <= limit; n++) {
            if ((periodLength(n) & 1) == 1) {
                count++;
            }
        }
        return count;
    }

    private static boolean runCheckpoints() {
        if (solve(13) != 4) {
            System.err.println("Checkpoint failed for limit=13");
            return false;
        }
        return true;
    }

    public static void main(String[] args) {
        Options options = new Options();
        if (!parseArguments(args, options)) {
            System.exit(1);
        }
        if (options.runCheckpoints && !runCheckpoints()) {
            System.exit(2);
        }

        System.out.println(solve(options.limit));
    }
}