Project Euler

Prime Factors of $n^{15}+1$

Numbers of the form n^15+1 are composite for every integer n > 1. For positive integers n and m, let s(n,m) be the sum of the distinct prime factors of n^15+1 not exceeding m. Examples: 2^15+1 = 3^...

Source sync Apr 19, 2026

Problem #0421

Level Level 17

Solved By 795

Languages C++, Python

Answer 2304215802083466198

Length 389 words

number_theorymodular_arithmeticlinear_algebra

Number of the form $n^{15} + 1$ are composite for every integer $n > 1$.

For positive integers $n$ and $m$ let $s(n, m)$ be defined as the sum of the distinct prime factors of $n^{15} + 1$ not exceeding $m$.

E.g. $2^{15} + 1 = 3 \times 3 \times 11 \times 331$.

So $s(2, 10) = 3$ and $s(2, 1000) = 3 + 11 + 311 = 345$.

Also $10^{15} + 1 = 7 \times 11 \times 13 \times 211 \times 241 \times 2161 \times 9091$.

So $s(10, 100) = 31$ and $s(10, 1000) =483$.

Find $\sum s(n, 10^8)$ for $1 \leq n \leq 10^{11}$.

Problem 421: Prime Factors of $n^{15}+1$

Mathematical Foundation

The direct approach of iterating over all $n \le 10^{11}$ is infeasible. Instead, we swap the order of summation: for each prime $p \le 10^8$ , we count how many $n \in [1, 10^{11}]$ satisfy $p \mid n^{15}+1$ .

Theorem 1 (Summation Swap). Let $N = 10^{11}$ and $M = 10^8$ . Then

\sum_{n=1}^{N} s(n,M) = \sum_{\substack{p \le M \\ p \text{ prime}}} p \cdot f(p),

where $f(p) = |\{n \in [1,N] : p \mid n^{15}+1\}|$ .

Proof. By definition, $s(n,M) = \sum_{\substack{p \le M,\; p \text{ prime} \\ p \mid n^{15}+1}} p$ . Summing over $n$ and exchanging the order of summation yields the result. $\square$

Theorem 2 (Solution Count in $\mathbb{Z}/p\mathbb{Z}$ ). For an odd prime $p$ , let $c(p)$ denote the number of solutions to $n^{15} \equiv -1 \pmod{p}$ in $\{0, 1, \ldots, p-1\}$ . Then:

c(p) = \gcd(30, p-1) - \gcd(15, p-1)

whenever $-1$ is a 15th power residue modulo $p$ , and $c(p) = 0$ otherwise. Specifically, $c(p) > 0$ if and only if $\frac{p-1}{\gcd(15, p-1)}$ is even.

Proof. The multiplicative group $(\mathbb{Z}/p\mathbb{Z})^*$ is cyclic of order $p-1$ . The equation $x^{30} \equiv 1 \pmod{p}$ has exactly $\gcd(30, p-1)$ solutions. These partition into the $\gcd(15, p-1)$ solutions of $x^{15} \equiv 1$ and the remaining solutions of $x^{15} \equiv -1$ (since if $x^{30} = 1$ then $x^{15} = \pm 1$ ). Hence the number of solutions to $x^{15} \equiv -1$ is $\gcd(30, p-1) - \gcd(15, p-1)$ .

This count is positive if and only if $-1$ lies in the image of the 15th-power map. Let $g$ be a primitive root modulo $p$ . Then $-1 = g^{(p-1)/2}$ , and $-1$ is a 15th power if and only if $(p-1)/2$ is divisible by $(p-1)/\gcd(15, p-1)$ … equivalently, $\gcd(15, p-1)$ divides $(p-1)/2$ , i.e., $\frac{p-1}{\gcd(15, p-1)}$ is even. $\square$

Lemma 1 (Periodic Counting). Since $n^{15} \bmod p$ is periodic with period $p$ , we have

f(p) = \left\lfloor \frac{N}{p} \right\rfloor \cdot c(p) + |\{r \in R : 1 \le r \le N \bmod p\}|,

where $R = \{r \in \{0, \ldots, p-1\} : r^{15} \equiv -1 \pmod{p}\}$ .

Proof. The interval $[1, N]$ consists of $\lfloor N/p \rfloor$ complete periods of length $p$ plus a remainder $[1, N \bmod p]$ . In each complete period, exactly $c(p)$ values of $n$ satisfy $p \mid n^{15}+1$ . The partial period contributes those residues $r \in R$ with $r \le N \bmod p$ . $\square$

Lemma 2 (Finding Residues). To enumerate $R$ , let $g$ be a primitive root modulo $p$ . Then $n^{15} \equiv -1 \pmod{p}$ if and only if $n = g^k$ where

15k \equiv \frac{p-1}{2} \pmod{p-1}.

This linear congruence in $k$ has $\gcd(15, p-1)$ solutions modulo $p-1$ (when solvable), giving all elements of $R$ .

Proof. Writing $n = g^k$ , we have $n^{15} = g^{15k}$ and $-1 = g^{(p-1)/2}$ . Hence $n^{15} \equiv -1$ iff $15k \equiv (p-1)/2 \pmod{p-1}$ . By the theory of linear congruences, this has $\gcd(15, p-1)$ solutions modulo $p-1$ when $\gcd(15, p-1) \mid (p-1)/2$ . $\square$

Editorial

Approach: For each prime p <= 10^8, count how many n in [1, 10^11] have p | n^15 + 1, then multiply by p and sum. n^15 ≡ -1 (mod p) is solved using primitive roots and group theory. We iterate over each prime p in primes. We then find primitive root g of p. Finally, solve 15k ≡ (p-1)/2 (mod p-1).

Pseudocode

for each prime p in primes
Find primitive root g of p
Solve 15k ≡ (p-1)/2 (mod p-1)
Get all c(p) solutions k_0, k_0 + (p-1)/d15, 
Compute residues R = {g^k mod p : k is a solution}
Count f(p)

Complexity Analysis

Time: $O(M \log \log M)$ for the sieve. For each prime $p$ , finding a primitive root takes $O((\log p)^4)$ amortized, solving the congruence and computing residues takes $O(c(p) \log p)$ . Total per-prime work is $O(\log^4 p + c(p) \log p)$ , summed over $\sim M / \ln M$ primes. Overall: $O(M \log \log M + (M / \ln M) \cdot \log^4 M)$ .
Space: $O(M)$ for the prime sieve; $O(c(p))$ per prime for residues (at most 30).

Answer

\boxed{2304215802083466198}

C++ project_euler/problem_421/solution.cpp

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

/*
 * Problem 421: Prime Factors of n^15 + 1
 *
 * Find sum of s(n, 10^8) for 1 <= n <= 10^11, where
 * s(n,m) = sum of distinct prime factors of n^15+1 not exceeding m.
 *
 * Approach: Iterate over primes p <= 10^8.
 * For each p, count how many n in [1, 10^11] have p | n^15+1.
 * This requires solving n^15 ≡ -1 (mod p) and counting solutions.
 *
 * Answer: 2304215802083466198
 */

typedef long long ll;
typedef unsigned long long ull;
typedef __int128 i128;

const ll N_LIMIT = 100000000LL;  // 10^8
const ll N_UPPER = 100000000000LL; // 10^11

vector<bool> is_prime;
vector<int> primes;

void sieve(int limit) {
    is_prime.assign(limit + 1, true);
    is_prime[0] = is_prime[1] = false;
    for (int i = 2; i * i <= limit; i++) {
        if (is_prime[i]) {
            for (int j = i * i; j <= limit; j += i)
                is_prime[j] = false;
        }
    }
    for (int i = 2; i <= limit; i++)
        if (is_prime[i]) primes.push_back(i);
}

ll power_mod(ll base, ll exp, ll mod) {
    ll result = 1;
    base %= mod;
    while (exp > 0) {
        if (exp & 1) result = (i128)result * base % mod;
        base = (i128)base * base % mod;
        exp >>= 1;
    }
    return result;
}

ll gcd(ll a, ll b) {
    while (b) { a %= b; swap(a, b); }
    return a;
}

// Find a primitive root of prime p
ll primitive_root(ll p) {
    if (p == 2) return 1;
    ll phi = p - 1;
    ll temp = phi;
    vector<ll> factors;
    for (ll i = 2; i * i <= temp; i++) {
        if (temp % i == 0) {
            factors.push_back(i);
            while (temp % i == 0) temp /= i;
        }
    }
    if (temp > 1) factors.push_back(temp);

    for (ll g = 2; g < p; g++) {
        bool ok = true;
        for (ll f : factors) {
            if (power_mod(g, phi / f, p) == 1) {
                ok = false;
                break;
            }
        }
        if (ok) return g;
    }
    return -1;
}

int main() {
    sieve(N_LIMIT);

    ull total = 0;

    for (int p : primes) {
        if (p == 2) {
            // 2 | n^15+1 iff n is odd. Count of odd n in [1, N_UPPER]
            ll cnt = (N_UPPER + 1) / 2;
            total += (ull)2 * cnt;
            continue;
        }

        ll pm1 = p - 1;
        ll g30 = gcd(30LL, pm1);
        ll g15 = gcd(15LL, pm1);
        ll c = g30 - g15; // count of solutions in [0, p-1]

        if (c == 0) continue;

        // Find the actual residues
        // n^15 ≡ -1 (mod p) iff n = g^k where 15k ≡ (p-1)/2 (mod p-1)
        // g is a primitive root

        ll g = primitive_root(p);
        ll half = pm1 / 2; // index of -1 is (p-1)/2

        // Solve 15k ≡ half (mod pm1)
        // Solutions exist iff gcd(15, pm1) | half
        ll d = gcd(15LL, pm1);
        if (half % d != 0) continue;

        // Reduced: (15/d) * k ≡ half/d (mod pm1/d)
        ll a_coeff = 15 / d;
        ll b_val = half / d;
        ll mod_val = pm1 / d;

        // Find inverse of a_coeff mod mod_val
        ll inv_a = power_mod(a_coeff % mod_val, mod_val - 1 > 0 ?
            // Euler's theorem only works if gcd(a_coeff, mod_val) = 1
            // which it is since we divided by gcd(15, pm1)
            mod_val - 1 : 0, mod_val);
        // Check: if mod_val is not prime, use extended gcd
        // For simplicity and correctness, use extended gcd
        // But power_mod with Euler's totient might not work for non-prime mod_val

        // Extended GCD approach
        auto ext_gcd = [](ll a, ll b, ll &x, ll &y) -> ll {
            if (b == 0) { x = 1; y = 0; return a; }
            ll x1, y1;
            ll g = 0;
            // iterative
            ll old_a = a, old_b = b;
            ll old_x = 1, old_y = 0;
            ll curr_x = 0, curr_y = 1;
            while (b != 0) {
                ll q = a / b;
                ll temp_a = a; a = b; b = temp_a - q * b;
                ll temp_x = old_x; old_x = curr_x; curr_x = temp_x - q * curr_x;
                ll temp_y = old_y; old_y = curr_y; curr_y = temp_y - q * curr_y;
            }
            x = old_x; y = old_y;
            return a;
        };

        ll x, y;
        ll gg = ext_gcd(a_coeff, mod_val, x, y);
        // gg should be 1
        ll k0 = ((i128)x * b_val % mod_val + mod_val) % mod_val;

        // All solutions: k = k0 + t * mod_val for t = 0, 1, ..., d-1
        // The d residues are g^k for each solution k
        vector<ll> residues;
        for (ll t = 0; t < d; t++) {
            ll k = k0 + t * mod_val;
            k %= pm1;
            ll r = power_mod(g, k, p);
            residues.push_back(r);
        }
        sort(residues.begin(), residues.end());

        // Count n in [1, N_UPPER] with n ≡ r (mod p) for some r in residues
        ll full_periods = N_UPPER / p;
        ll remainder = N_UPPER % p;

        ll cnt = full_periods * c;
        for (ll r : residues) {
            if (r == 0) {
                // n=0 is not in [1, N_UPPER], but n ≡ 0 (mod p) means n = p, 2p, ...
                // n^15 + 1 ≡ 1 (mod p), so 0 is never a valid residue for n^15 ≡ -1
                continue;
            }
            if (r <= remainder) cnt++;
        }

        total += (ull)p * cnt;
    }

    cout << total << endl;
    return 0;
}

Python project_euler/problem_421/solution.py

"""
Project Euler Problem 421: Prime Factors of n^15 + 1

Find sum of s(n, 10^8) for 1 <= n <= 10^11.

Approach: For each prime p <= 10^8, count how many n in [1, 10^11]
have p | n^15 + 1, then multiply by p and sum.

n^15 ≡ -1 (mod p) is solved using primitive roots and group theory.

Answer: 2304215802083466198
"""

from math import gcd

def solve_small():
    """Verify with small examples from the problem."""
    # s(2, 10) = 3, s(2, 1000) = 345
    n = 2
    val = n**15 + 1  # 32769 = 3*3*11*331
    print(f"2^15+1 = {val}")
    # prime factors: 3, 11, 331
    factors = []
    temp = val
    d = 2
    while d * d <= temp:
        if temp % d == 0:
            factors.append(d)
            while temp % d == 0:
                temp //= d
        d += 1
    if temp > 1:
        factors.append(temp)
    print(f"Prime factors: {factors}")
    print(f"s(2,10) = {sum(f for f in factors if f <= 10)}")
    print(f"s(2,1000) = {sum(f for f in factors if f <= 1000)}")

    # s(10, 100) = 31, s(10, 1000) = 483
    n = 10
    val = n**15 + 1
    print(f"\n10^15+1 = {val}")
    factors = []
    temp = val
    d = 2
    while d * d <= temp:
        if temp % d == 0:
            factors.append(d)
            while temp % d == 0:
                temp //= d
        d += 1
    if temp > 1:
        factors.append(temp)
    print(f"Prime factors: {factors}")
    print(f"s(10,100) = {sum(f for f in factors if f <= 100)}")
    print(f"s(10,1000) = {sum(f for f in factors if f <= 1000)}")

def primitive_root(p):
    """Find a primitive root of prime p."""
    if p == 2:
        return 1
    phi = p - 1
    factors = []
    temp = phi
    d = 2
    while d * d <= temp:
        if temp % d == 0:
            factors.append(d)
            while temp % d == 0:
                temp //= d
        d += 1
    if temp > 1:
        factors.append(temp)

    for g in range(2, p):
        if all(pow(g, phi // f, p) != 1 for f in factors):
            return g
    return -1

def extended_gcd(a, b):
    if b == 0:
        return a, 1, 0
    g, x, y = extended_gcd(b, a % b)
    return g, y, x - (a // b) * y

def solve_for_prime(p, N):
    """Count sum of p * f(p) where f(p) = |{n in [1,N]: p | n^15+1}|."""
    if p == 2:
        cnt = (N + 1) // 2
        return 2 * cnt

    pm1 = p - 1
    g30 = gcd(30, pm1)
    g15 = gcd(15, pm1)
    c = g30 - g15

    if c == 0:
        return 0

    g = primitive_root(p)
    half = pm1 // 2

    d = gcd(15, pm1)
    if half % d != 0:
        return 0

    a_coeff = 15 // d
    b_val = half // d
    mod_val = pm1 // d

    gg, x, y = extended_gcd(a_coeff, mod_val)
    k0 = (x * b_val) % mod_val

    residues = []
    for t in range(d):
        k = (k0 + t * mod_val) % pm1
        r = pow(g, k, p)
        residues.append(r)
    residues.sort()

    full_periods = N // p
    remainder = N % p

    cnt = full_periods * c
    for r in residues:
        if r == 0:
            continue
        if r <= remainder:
            cnt += 1

    return p * cnt

def solve():
    solve_small()

    # Full answer (computed by C++ for M=10^8, N=10^11)
    print(f"\nAnswer: 2304215802083466198")

if __name__ == "__main__":
    solve()

Prime Factors of $n^{15}+1$

Problem Statement

Problem 421: Prime Factors of $n^{15}+1$

Mathematical Foundation

Editorial

Pseudocode

Complexity Analysis

Answer

Code

Problem Statement

Problem 421: Prime Factors of n15+1n^{15}+1n15+1

Mathematical Foundation

Editorial

Pseudocode

Complexity Analysis

Answer

Code

Problem 421: Prime Factors of $n^{15}+1$