Problem 293: Pseudo-Fortunate Numbers

Project Euler Problem 293 Solution

EulerSolve provides an optimized solution for Project Euler Problem 293, Pseudo-Fortunate Numbers, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary An even integer $n$ is called admissible when it is either a power of $2$ or, equivalently for even numbers, its distinct prime factors form a consecutive prefix of the prime list $2,3,5,7,\ldots$. For every admissible $n < 10^9$, define $M(n)$ as the smallest integer $m>1$ such that $n+m$ is prime. The task is to sum the distinct values of $M(n)$. Mathematical Approach 1. Exact shape of admissible numbers Because $n$ is even, the first prime $2$ must appear. If the distinct prime factors are consecutive, then every admissible number has the form $$n = 2^{a_1} 3^{a_2} 5^{a_3} \cdots p_k^{a_k}, \qquad a_i \ge 1.$$ The special case "power of $2$" is simply the case $k=1$. So admissible numbers are described completely by choosing a prime prefix and choosing a positive exponent for each prime in that prefix. This already gives a strong bound on the search space. The product of the first nine primes is $$2\cdot 3\cdot 5\cdot 7\cdot 11\cdot 13\cdot 17\cdot 19\cdot 23 = 223092870 < 10^9,$$ while including the next prime gives $$223092870 \cdot 29 = 6469693230 > 10^9.$$ Therefore no admissible number below $10^9$ can contain $29$ or any later prime as a distinct factor. The code keeps a slightly longer fixed prime list, but recursion automatically prunes everything beyond the limit. 2....

Detailed mathematical approach

Problem Summary

An even integer $n$ is called admissible when it is either a power of $2$ or, equivalently for even numbers, its distinct prime factors form a consecutive prefix of the prime list $2,3,5,7,\ldots$. For every admissible $n < 10^9$, define $M(n)$ as the smallest integer $m>1$ such that $n+m$ is prime. The task is to sum the distinct values of $M(n)$.

Mathematical Approach

1. Exact shape of admissible numbers

Because $n$ is even, the first prime $2$ must appear. If the distinct prime factors are consecutive, then every admissible number has the form

$$n = 2^{a_1} 3^{a_2} 5^{a_3} \cdots p_k^{a_k}, \qquad a_i \ge 1.$$

The special case "power of $2$" is simply the case $k=1$. So admissible numbers are described completely by choosing a prime prefix and choosing a positive exponent for each prime in that prefix.

This already gives a strong bound on the search space. The product of the first nine primes is

$$2\cdot 3\cdot 5\cdot 7\cdot 11\cdot 13\cdot 17\cdot 19\cdot 23 = 223092870 < 10^9,$$

while including the next prime gives

$$223092870 \cdot 29 = 6469693230 > 10^9.$$

Therefore no admissible number below $10^9$ can contain $29$ or any later prime as a distinct factor. The code keeps a slightly longer fixed prime list, but recursion automatically prunes everything beyond the limit.

2. Why the DFS generation is exact

The recursive generator works on the ordered prime list $(2,3,5,7,\ldots)$. Suppose we are currently at prime $p_i$ with a partial product built from the earlier prefix. The loop repeatedly multiplies by $p_i$:

$$\text{current},\ \text{current}\cdot p_i,\ \text{current}\cdot p_i^2,\ \text{current}\cdot p_i^3,\ldots$$

as long as the product stays below the limit. Every time at least one copy of $p_i$ has been chosen, the recursion may continue to the next prime $p_{i+1}$.

This does exactly what we need:

Positive exponents. Once a prime is included, it appears with exponent at least $1$.

Consecutive prefix only. The recursion may move from $p_i$ to $p_{i+1}$ only after choosing $p_i$, so later primes can never appear while skipping an earlier one.

No duplicates. Each exponent vector $(a_1,\ldots,a_k)$ corresponds to one unique recursion path.

So the DFS enumerates every admissible number once and only once, without scanning all even integers up to $10^9$.

3. Why the pseudo-Fortunate search starts at $3$

By definition, $M(n)$ is the smallest integer $m>1$ such that $n+m$ is prime. Since every admissible $n$ is even, any prime larger than $2$ has to be odd, so $m$ must be odd as well. Therefore the candidates are exactly

$$m = 3,5,7,9,\ldots$$

and the code tests them in increasing order until $n+m$ is prime.

For example, $16$ is admissible because it is a power of $2$. We cannot use $m=1$ because the definition requires $m>1$. Then

$$16+3 = 19,$$

which is prime, so

$$M(16)=3.$$

Likewise, $630 = 2\cdot 3^2\cdot 5\cdot 7$ is admissible because its distinct primes are $2,3,5,7$. Testing odd values gives

$$630+3=633,\ 630+5=635,\ 630+7=637,\ 630+9=639,$$

all composite, while

$$630+11=641$$

is prime, hence

$$M(630)=11.$$

4. Why we need a set of pseudo-Fortunate values

The problem asks for the sum of distinct values $M(n)$, not the sum over all admissible $n$. Different admissible numbers can share the same pseudo-Fortunate value, so the implementation inserts each result into a set first:

$$\mathcal{M} = \{M(n) : n \text{ admissible},\ n < 10^9\}.$$

The final answer is then

$$\sum_{m \in \mathcal{M}} m.$$

This is mathematically important: without deduplication we would answer a different question.

5. Primality testing is the main inner operation

After admissible numbers have been generated, the expensive repeated task is to test whether $n+m$ is prime. The code uses deterministic Miller-Rabin for 64-bit integers, so each test is exact in this range and still very fast. That makes the straightforward odd-step search practical.

6. Checkpoints and interpretation

The implementation verifies itself with three simple checks:

$$M(16)=3,\qquad M(630)=11,\qquad \text{solve}(10000)=\text{solve\_bruteforce}(10000).$$

The first two confirm the pseudo-Fortunate search itself. The third compares the efficient DFS enumeration against a brute-force admissibility test on a smaller limit, which is a strong check that the generator is complete and duplicate-free.

How the Code Works

prime_prefix() provides a fixed initial segment of the prime list. generate_admissible(...) performs the DFS over exponent choices and stores every admissible product in a set. For each stored $n$, pseudo_fortunate(n) checks odd $m=3,5,7,\ldots$ until $n+m$ is prime. These $M(n)$ values are inserted into another set, and the program sums that set at the end.

Complexity Analysis

Let $A$ be the number of admissible numbers below the limit, and let $T(n)$ be the number of odd candidates tested before finding $M(n)$. Then the total running time is roughly

$$O\!\left(A + \sum_{n \in \text{admissible}} T(n)\cdot C_{\text{prime}}\right),$$

where $C_{\text{prime}}$ is the cost of one primality test. Memory usage is $O(A + |\mathcal{M}|)$ for the two sets. The crucial gain comes from generating only admissible numbers rather than scanning all even numbers up to $10^9$.

Problem 293 source code

C++

#include <algorithm>
#include <cstdint>
#include <iostream>
#include <set>
#include <string>
#include <vector>

namespace {

using u64 = std::uint64_t;
using u128 = unsigned __int128;

struct Options {
    u64 limit = 1000000000ULL;
    bool run_checkpoints = true;
};

bool parse_u64_after_prefix(const std::string& arg, const std::string& prefix, u64& value) {
    if (arg.rfind(prefix, 0U) != 0U) {
        return false;
    }
    const std::string tail = arg.substr(prefix.size());
    if (tail.empty()) {
        return false;
    }
    u64 parsed = 0ULL;
    for (char c : tail) {
        if (c < '0' || c > '9') {
            return false;
        }
        parsed = parsed * 10ULL + static_cast<u64>(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_u64_after_prefix(arg, "--limit=", options.limit)) {
            continue;
        }
        std::cerr << "Unknown argument: " << arg << '\n';
        return false;
    }
    return options.limit >= 3ULL;
}

u64 mod_pow(u64 base, u64 exp, u64 mod) {
    u64 result = 1ULL % mod;
    u64 cur = base % mod;
    while (exp > 0ULL) {
        if ((exp & 1ULL) != 0ULL) {
            result = static_cast<u64>((static_cast<u128>(result) * cur) % mod);
        }
        cur = static_cast<u64>((static_cast<u128>(cur) * cur) % mod);
        exp >>= 1U;
    }
    return result;
}

bool is_prime(const u64 n) {
    if (n < 2ULL) {
        return false;
    }
    for (u64 p : {2ULL, 3ULL, 5ULL, 7ULL, 11ULL, 13ULL, 17ULL, 19ULL, 23ULL, 29ULL, 31ULL, 37ULL}) {
        if (n == p) {
            return true;
        }
        if (n % p == 0ULL) {
            return false;
        }
    }

    u64 d = n - 1ULL;
    int s = 0;
    while ((d & 1ULL) == 0ULL) {
        d >>= 1U;
        ++s;
    }

    for (u64 a : {2ULL, 325ULL, 9375ULL, 28178ULL, 450775ULL, 9780504ULL, 1795265022ULL}) {
        if (a % n == 0ULL) {
            continue;
        }
        u64 x = mod_pow(a, d, n);
        if (x == 1ULL || x == n - 1ULL) {
            continue;
        }
        bool witness = true;
        for (int r = 1; r < s; ++r) {
            x = static_cast<u64>((static_cast<u128>(x) * x) % n);
            if (x == n - 1ULL) {
                witness = false;
                break;
            }
        }
        if (witness) {
            return false;
        }
    }
    return true;
}

std::vector<int> prime_prefix() {
    return {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31};
}

void generate_admissible(const std::vector<int>& primes,
                         const int idx,
                         const u64 current,
                         const u64 limit,
                         std::set<u64>& out) {
    if (idx >= static_cast<int>(primes.size())) {
        return;
    }
    const u64 p = static_cast<u64>(primes[static_cast<std::size_t>(idx)]);
    u64 value = current;
    while (value <= (limit - 1ULL) / p) {
        value *= p;
        out.insert(value);
        generate_admissible(primes, idx + 1, value, limit, out);
    }
}

int pseudo_fortunate(const u64 n) {
    for (int m = 3;; m += 2) {
        if (is_prime(n + static_cast<u64>(m))) {
            return m;
        }
    }
}

u64 solve(const u64 limit) {
    std::set<u64> admissible;
    const auto primes = prime_prefix();
    generate_admissible(primes, 0, 1ULL, limit, admissible);

    std::set<int> pseudo_values;
    for (u64 n : admissible) {
        pseudo_values.insert(pseudo_fortunate(n));
    }

    u64 sum = 0ULL;
    for (int m : pseudo_values) {
        sum += static_cast<u64>(m);
    }
    return sum;
}

bool is_admissible_bruteforce(u64 n) {
    if ((n & 1ULL) == 1ULL) {
        return false;
    }
    std::vector<u64> factors;
    u64 x = n;
    for (u64 p = 2ULL; p * p <= x; ++p) {
        if (x % p != 0ULL) {
            continue;
        }
        factors.push_back(p);
        while (x % p == 0ULL) {
            x /= p;
        }
    }
    if (x > 1ULL) {
        factors.push_back(x);
    }
    if (factors.empty() || factors[0] != 2ULL) {
        return false;
    }
    u64 candidate = 2ULL;
    std::size_t idx = 0;
    while (idx < factors.size()) {
        if (factors[idx] != candidate) {
            return false;
        }
        ++idx;
        ++candidate;
        while (!is_prime(candidate)) {
            ++candidate;
        }
    }
    return true;
}

u64 solve_bruteforce(const u64 limit) {
    std::set<int> pseudo_values;
    for (u64 n = 2; n < limit; n += 2) {
        if (!is_admissible_bruteforce(n)) {
            continue;
        }
        pseudo_values.insert(pseudo_fortunate(n));
    }
    u64 sum = 0ULL;
    for (int m : pseudo_values) {
        sum += static_cast<u64>(m);
    }
    return sum;
}

bool run_checkpoints() {
    if (pseudo_fortunate(16ULL) != 3) {
        std::cerr << "Checkpoint failed for N=16 sample" << '\n';
        return false;
    }
    if (pseudo_fortunate(630ULL) != 11) {
        std::cerr << "Checkpoint failed for N=630 sample" << '\n';
        return false;
    }
    if (solve(10000ULL) != solve_bruteforce(10000ULL)) {
        std::cerr << "Checkpoint failed for brute cross-check at limit=10000" << '\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

def mod_pow(base, exp, mod):
    return pow(base, exp, mod)

def is_prime(n):
    if n < 2:
        return False
    for p in [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37]:
        if n == p:
            return True
        if n % p == 0:
            return False

    d = n - 1
    s = 0
    while (d & 1) == 0:
        d >>= 1
        s += 1

    for a in [2, 325, 9375, 28178, 450775, 9780504, 1795265022]:
        if a % n == 0:
            continue
        x = mod_pow(a, d, n)
        if x == 1 or x == n - 1:
            continue
        witness = True
        for _ in range(1, s):
            x = (x * x) % n
            if x == n - 1:
                witness = False
                break
        if witness:
            return False
    return True

def prime_prefix():
    return [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31]

def generate_admissible(primes, idx, current, limit, out):
    if idx >= len(primes):
        return
    p = primes[idx]
    value = current
    while value <= (limit - 1) // p:
        value *= p
        out.add(value)
        generate_admissible(primes, idx + 1, value, limit, out)

def pseudo_fortunate(n):
    m = 3
    while True:
        if is_prime(n + m):
            return m
        m += 2

def solve(limit=1000000000):
    admissible = set()
    primes = prime_prefix()
    generate_admissible(primes, 0, 1, limit, admissible)

    pseudo_values = set()
    for n in admissible:
        pseudo_values.add(pseudo_fortunate(n))

    return str(sum(pseudo_values))

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

Java

import java.util.*;

public class Euler293 {
    static long modPow(long base, long exp, long mod) {
        long result = 1L % mod;
        long cur = base % mod;
        while (exp > 0) {
            if ((exp & 1) != 0) {
                // To avoid overflow, use BigInteger for multiplication modulo
                result = java.math.BigInteger.valueOf(result)
                        .multiply(java.math.BigInteger.valueOf(cur))
                        .mod(java.math.BigInteger.valueOf(mod)).longValue();
            }
            cur = java.math.BigInteger.valueOf(cur)
                    .multiply(java.math.BigInteger.valueOf(cur))
                    .mod(java.math.BigInteger.valueOf(mod)).longValue();
            exp >>= 1;
        }
        return result;
    }

    static boolean isPrime(long n) {
        if (n < 2)
            return false;
        long[] smallPrimes = { 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37 };
        for (long p : smallPrimes) {
            if (n == p)
                return true;
            if (n % p == 0)
                return false;
        }

        long d = n - 1;
        int s = 0;
        while ((d & 1) == 0) {
            d >>= 1;
            ++s;
        }

        long[] bases = { 2, 325, 9375, 28178, 450775, 9780504, 1795265022 };
        for (long a : bases) {
            if (a % n == 0)
                continue;
            long x = modPow(a, d, n);
            if (x == 1 || x == n - 1)
                continue;
            boolean witness = true;
            for (int r = 1; r < s; ++r) {
                x = java.math.BigInteger.valueOf(x)
                        .multiply(java.math.BigInteger.valueOf(x))
                        .mod(java.math.BigInteger.valueOf(n)).longValue();
                if (x == n - 1) {
                    witness = false;
                    break;
                }
            }
            if (witness)
                return false;
        }
        return true;
    }

    static void generateAdmissible(int[] primes, int idx, long current, long limit, Set<Long> out) {
        if (idx >= primes.length)
            return;
        long p = primes[idx];
        long value = current;
        while (value <= (limit - 1) / p) {
            value *= p;
            out.add(value);
            generateAdmissible(primes, idx + 1, value, limit, out);
        }
    }

    static int pseudoFortunate(long n) {
        for (int m = 3;; m += 2) {
            if (isPrime(n + m)) {
                return m;
            }
        }
    }

    public static String solve() {
        long limit = 1000000000L;
        Set<Long> admissible = new HashSet<>();
        int[] primes = { 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31 };
        generateAdmissible(primes, 0, 1L, limit, admissible);

        Set<Integer> pseudoValues = new HashSet<>();
        for (long n : admissible) {
            pseudoValues.add(pseudoFortunate(n));
        }

        long sum = 0;
        for (int m : pseudoValues) {
            sum += m;
        }
        return String.valueOf(sum);
    }

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