Problem 318: 2011 Nines

Project Euler Problem 318 Solution

EulerSolve provides an optimized solution for Project Euler Problem 318, 2011 Nines, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary For every pair of positive integers $$1\le p \lt q,\qquad p+q\le 2011,$$ consider $$\alpha=\sqrt p+\sqrt q.$$ We look at the even powers $\alpha^{2n}$, and $C(p,q,n)$ is the number of consecutive 9s at the beginning of the fractional part of $\alpha^{2n}$. Let $N(p,q)$ be the smallest $n$ such that $$C(p,q,n)\ge 2011.$$ The goal is to compute $$\sum_{p+q\le 2011} N(p,q),$$ but only for those pairs where the fractional part of $\alpha^{2n}$ approaches $1$. Mathematical Approach 1) Introduce the conjugate partner. Set $$\beta=\sqrt q-\sqrt p.$$ Then $$\alpha^{2n}+\beta^{2n}$$ is always an integer. The reason is that in the binomial expansions of $(\sqrt q+\sqrt p)^{2n}$ and $(\sqrt q-\sqrt p)^{2n}$, all odd-radical terms cancel and only integer terms remain. So if we define $$A_n=\alpha^{2n}+\beta^{2n}\in\mathbb Z,$$ then $$\alpha^{2n}=A_n-\beta^{2n}.$$ 2) Exactly when does the fractional part approach $1$? Because $p \lt q$, we have $\beta \gt 0$. The fractional part of $\alpha^{2n}$ approaches $1$ exactly when $\beta^{2n}\to0$, i.e. $$0\lt \beta \lt 1.$$ If $\beta\ge1$, then $\beta^{2n}$ does not decay to $0$, so the fractional part cannot converge to $1$. Thus the admissible pairs are precisely those with $$\sqrt q-\sqrt p \lt 1.$$ 3) Fractional part as $1-\beta^{2n}$....

Detailed mathematical approach

Problem Summary

For every pair of positive integers

$$1\le p \lt q,\qquad p+q\le 2011,$$

consider

$$\alpha=\sqrt p+\sqrt q.$$

We look at the even powers $\alpha^{2n}$, and $C(p,q,n)$ is the number of consecutive 9s at the beginning of the fractional part of $\alpha^{2n}$.

Let $N(p,q)$ be the smallest $n$ such that

$$C(p,q,n)\ge 2011.$$

The goal is to compute

$$\sum_{p+q\le 2011} N(p,q),$$

but only for those pairs where the fractional part of $\alpha^{2n}$ approaches $1$.

Mathematical Approach

1) Introduce the conjugate partner.

Set

$$\beta=\sqrt q-\sqrt p.$$

Then

$$\alpha^{2n}+\beta^{2n}$$

is always an integer. The reason is that in the binomial expansions of $(\sqrt q+\sqrt p)^{2n}$ and $(\sqrt q-\sqrt p)^{2n}$, all odd-radical terms cancel and only integer terms remain.

So if we define

$$A_n=\alpha^{2n}+\beta^{2n}\in\mathbb Z,$$

then

$$\alpha^{2n}=A_n-\beta^{2n}.$$

2) Exactly when does the fractional part approach $1$?

Because $p \lt q$, we have $\beta \gt 0$. The fractional part of $\alpha^{2n}$ approaches $1$ exactly when $\beta^{2n}\to0$, i.e.

$$0\lt \beta \lt 1.$$

If $\beta\ge1$, then $\beta^{2n}$ does not decay to $0$, so the fractional part cannot converge to $1$.

Thus the admissible pairs are precisely those with

$$\sqrt q-\sqrt p \lt 1.$$

3) Fractional part as $1-\beta^{2n}$.

For every admissible pair we have $0\lt \beta^{2n}\lt1$, hence

$$\alpha^{2n}=A_n-\beta^{2n}$$

lies just below the integer $A_n$. Therefore

$$\{\alpha^{2n}\}=1-\beta^{2n}.$$

This is the whole reason the problem becomes easy: the complicated-looking irrational power is controlled by the tiny positive quantity $\beta^{2n}$.

4) Translate "at least $K$ leading nines".

Let $x=\{\alpha^{2n}\}$. The fractional part begins with at least $K$ nines if and only if

$$x \gt 1-10^{-K}.$$

Since $x=1-\beta^{2n}$, this is equivalent to

$$\beta^{2n}\lt 10^{-K}.$$

Here $K=2011$ for the actual problem.

5) Take logarithms.

Because $0\lt \beta^2 \lt 1$, the quantity

$$\lambda=-\log_{10}(\beta^2)$$

is positive. The inequality above becomes

$$n\lambda \gt K.$$

So the minimal valid exponent is

$$N(p,q)=\left\lceil\frac{K}{\lambda}\right\rceil=\left\lceil\frac{K}{-\log_{10}(\beta^2)}\right\rceil.$$

This is exactly what the C++ function minimal_n computes, with a tiny tolerance to avoid floating-point boundary mistakes.

6) Worked example: $(p,q)=(2,3)$.

Here

$$\beta=\sqrt3-\sqrt2\approx0.3178372452,$$

$$\beta^2\approx0.1010205144,$$

and

$$\lambda=-\log_{10}(\beta^2)\approx0.9955904242.$$

For one leading 9 we need $K=1$, hence

$$N(2,3)=\left\lceil\frac{1}{0.9955904242}\right\rceil=2.$$

That matches the sequence shown in the statement:

$(\sqrt2+\sqrt3)^2=9.8989\ldots$ has no leading 9 in the fractional part,

$(\sqrt2+\sqrt3)^4=97.9897\ldots$ has one,

$(\sqrt2+\sqrt3)^6=969.9989\ldots$ has two,

$(\sqrt2+\sqrt3)^8=9601.9998\ldots$ has three.

So the checkpoints

$$N(2,3;K=1)=2,\qquad N(2,3;K=2)=3,\qquad N(2,3;K=3)=4$$

are exactly correct.

7) Small full example: $M=5,\ K=1$.

The admissible pairs are

$$ (1,2),\ (1,3),\ (2,3). $$

Their minimal values are

$$N(1,2)=2,\qquad N(1,3)=4,\qquad N(2,3)=2.$$

Therefore

$$S(5,1)=2+4+2=8,$$

which is the second checkpoint in the source code.

8) Final summation formula.

The complete answer is

$$\sum_{\substack{1\le p \lt q\\ p+q\le 2011\\ \sqrt q-\sqrt p\lt1}} \left\lceil\frac{2011}{-\log_{10}\!\left((\sqrt q-\sqrt p)^2\right)}\right\rceil.$$

No deeper number theory is needed after this reduction; the implementation simply iterates over all candidate pairs.

Algorithm

1) Loop over all pairs $(p,q)$ with $1\le p \lt q$ and $p+q\le M$.

2) Compute

$$\beta=\sqrt q-\sqrt p.$$

3) Skip the pair if $\beta\ge1$.

4) Compute $\lambda=-\log_{10}(\beta^2)$.

5) Add

$$\left\lceil\frac{K}{\lambda}\right\rceil$$

to the running total.

Complexity Analysis

The triangular region $p+q\le M$ contains

$$O(M^2)$$

pairs. Each pair needs a constant amount of work: two square roots, one logarithm, and a few arithmetic operations. So the total complexity is

$$O(M^2)$$

time and

$$O(1)$$

memory.

Checks And Final Result

The implementation checks

$$N(2,3;1)=2,\qquad N(2,3;2)=3,\qquad N(2,3;3)=4,$$

and

$$S(5,1)=8.$$

For the full problem $(M,K)=(2011,2011)$, the final answer is

$$709313889.$$

Problem 318 source code

C++

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

namespace {

using u64 = std::uint64_t;

struct Options {
    int max_sum = 2011;
    int target_nines = 2011;
    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, "--max-sum=", options.max_sum) ||
            parse_int_after_prefix(arg, "--target-nines=", options.target_nines)) {
            continue;
        }
        std::cerr << "Unknown argument: " << arg << '\n';
        return false;
    }
    return options.max_sum >= 3 && options.target_nines >= 1;
}

u64 minimal_n(const int p, const int q, const int target_nines) {
    const long double beta = std::sqrt(static_cast<long double>(q)) -
                             std::sqrt(static_cast<long double>(p));
    const long double beta2 = beta * beta;
    const long double lambda = -std::log10(beta2);
    const long double target = static_cast<long double>(target_nines);

    u64 n = static_cast<u64>(target / lambda);
    if (n == 0ULL) {
        n = 1ULL;
    }
    while (static_cast<long double>(n) * lambda < target - 1e-15L) {
        ++n;
    }
    while (n > 1ULL && static_cast<long double>(n - 1ULL) * lambda >= target - 1e-15L) {
        --n;
    }
    return n;
}

u64 solve(const int max_sum, const int target_nines) {
    u64 total = 0ULL;
    for (int p = 1; p < max_sum; ++p) {
        for (int q = p + 1; p + q <= max_sum; ++q) {
            const long double beta = std::sqrt(static_cast<long double>(q)) -
                                     std::sqrt(static_cast<long double>(p));
            if (beta >= 1.0L) {
                continue;
            }
            total += minimal_n(p, q, target_nines);
        }
    }
    return total;
}

bool run_checkpoints() {
    if (minimal_n(2, 3, 1) != 2ULL || minimal_n(2, 3, 2) != 3ULL || minimal_n(2, 3, 3) != 4ULL) {
        std::cerr << "Checkpoint failed for sqrt(2)+sqrt(3) nines progression" << '\n';
        return false;
    }
    if (solve(5, 1) != 8ULL) {
        std::cerr << "Checkpoint failed for small max-sum=5 target=1" << '\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.max_sum, options.target_nines) << '\n';
    return 0;
}

Python

import math

def minimal_n(p, q, target_nines):
    beta = math.sqrt(q) - math.sqrt(p)
    beta2 = beta * beta
    lam = -math.log10(beta2)
    target = float(target_nines)
    
    n = int(target / lam)
    if n == 0:
        n = 1
    while n * lam < target - 1e-15:
        n += 1
    while n > 1 and (n - 1) * lam >= target - 1e-15:
        n -= 1
    return n

def solve(max_sum=2011, target_nines=2011):
    total = 0
    for p in range(1, max_sum):
        for q in range(p + 1, max_sum - p + 1):
            beta = math.sqrt(q) - math.sqrt(p)
            if beta >= 1.0:
                continue
            total += minimal_n(p, q, target_nines)
    return str(total)

if __name__ == '__main__':
    print(solve())

Java

public class Euler318 {
    static long minimalN(int p, int q, int targetNines) {
        double beta = Math.sqrt(q) - Math.sqrt(p);
        double beta2 = beta * beta;
        double lambda = -Math.log10(beta2);
        double target = (double) targetNines;

        long n = (long) (target / lambda);
        if (n == 0) {
            n = 1;
        }
        while (n * lambda < target - 1e-15) {
            n++;
        }
        while (n > 1 && (n - 1) * lambda >= target - 1e-15) {
            n--;
        }
        return n;
    }

    public static String solve() {
        int maxSum = 2011;
        int targetNines = 2011;
        long total = 0;

        for (int p = 1; p < maxSum; ++p) {
            for (int q = p + 1; p + q <= maxSum; ++q) {
                double beta = Math.sqrt(q) - Math.sqrt(p);
                if (beta >= 1.0) {
                    continue;
                }
                total += minimalN(p, q, targetNines);
            }
        }
        return String.valueOf(total);
    }

    public static void main(String[] args) {
        System.out.println(solve());
    }
}