Project Euler

Additive Prime Counting

An additive prime is a prime whose digit sum is also prime. Find the count of additive primes below 10^6.

Source sync Apr 19, 2026

Problem #0969

Level Level 33

Solved By 248

Languages C++, Python

Answer 412543690

Length 215 words

number_theorymodular_arithmeticprobability

Starting at zero, a kangaroo hops along the real number line in the positive direction. Each successive hop takes the kangaroo forward a uniformly random distance between and . Let be the expected number of hops needed for the kangaroo to pass on the real line.

If we write , then for all positive integers , can be expressed as a polynomial function of with rational coefficients. For example . Define to be the sum of all integer coefficients in this polynomial form of . Therefore and .
You are also given .
Find . Give your answer modulo .

Problem 969: Additive Prime Counting

Mathematical Analysis

Digit Sum Properties

Definition. The digit sum $d(n) = \sum_{i} d_i$ where $d_i$ are the base-10 digits of $n$ .

Theorem. For $n < 10^6$ , the digit sum satisfies $1 \le d(n) \le 54$ (max is $999999$ with digit sum $54$ ).

Proposition. $n \equiv d(n) \pmod{9}$ , so for prime $p > 3$ , $d(p) \not\equiv 0 \pmod{3}$ (otherwise $3 \mid p$ , contradiction). Thus $d(p)$ avoids multiples of 3 except when $p = 3$ .

Density Estimate

Among primes below $N$ , the digit sum $d(p)$ is roughly normally distributed with mean $\approx 4.5 \cdot \log_{10} N$ and standard deviation $\approx \sqrt{4.5 \cdot \log_{10} N \cdot 0.83}$ . The probability that a random integer in the digit-sum range is prime is approximately $1/\ln(\text{digit sum})$ .

For $N = 10^6$ : average digit sum $\approx 27$ , and primes up to 54 include 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53 (16 primes out of 54). Roughly 30% of digit sums are prime, so we expect $\approx 0.3 \times 78498 \approx 23500$ additive primes. The actual count is slightly higher.

Derivation

Editorial

These are called “additive primes.”. We sieve primes below $10^6$ . We then precompute the set of primes up to 54. Finally, iterate over each prime $p$ , compute digit sum and check if it is prime.

Pseudocode

Sieve primes below $10^6$
Precompute the set of primes up to 54
For each prime $p$, compute digit sum and check if it is prime
Count matches

Proof of Correctness

The sieve is exact. Digit sum computation by repeated modular division is exact. Primality check for numbers up to 54 is trivial.

Complexity Analysis

$O(N \log \log N)$ for sieve, $O(\pi(N) \cdot \log_{10} N)$ for digit sum checks.

Answer

\boxed{412543690}

C++ project_euler/problem_969/solution.cpp

#include <bits/stdc++.h>
using namespace std;
int main() {
    const int N = 1000000;
    vector<bool> is_p(N, true); is_p[0]=is_p[1]=false;
    for (int i=2;(long long)i*i<N;i++) if(is_p[i]) for(int j=i*i;j<N;j+=i) is_p[j]=false;
    set<int> small_p;
    for (int i=2;i<55;i++) if(is_p[i]) small_p.insert(i);
    int cnt=0;
    for(int p=2;p<N;p++){
        if(!is_p[p]) continue;
        int ds=0,t=p; while(t){ds+=t%10;t/=10;}
        if(small_p.count(ds)) cnt++;
    }
    cout<<cnt<<endl;
    return 0;
}

Python project_euler/problem_969/solution.py

"""
Problem 969: Additive Prime Counting

Find the number of primes p < 10^6 such that the digit sum of p is also prime.
These are called "additive primes."

Key results:
    - A prime p is "additive" if digit_sum(p) is also prime.
    - For p < 10^6, digit sums range from 2 to 54, so we only need
      small primes up to 54.
    - answer = count of additive primes below 10^6

Methods:
    1. sieve_of_eratosthenes  — standard boolean sieve
    2. digit_sum              — iterative digit extraction
    3. count_additive_primes  — main scan with digit-sum primality check
    4. additive_ratio_by_magnitude — fraction of primes that are additive per decade
    5. digit_sum_histogram    — full digit-sum frequency among primes
"""

from math import isqrt
from collections import Counter

def sieve_of_eratosthenes(n):
    """Return bytearray where s[i] is 1 iff i is prime, for 0 <= i <= n."""
    s = bytearray(b'\x01') * (n + 1)
    s[0] = s[1] = 0
    for i in range(2, isqrt(n) + 1):
        if s[i]:
            s[i * i :: i] = bytearray(len(s[i * i :: i]))
    return s

def digit_sum(n):
    """Sum of decimal digits of n."""
    s = 0
    while n:
        s += n % 10
        n //= 10
    return s

def count_additive_primes(N):
    """Return (count, is_prime_array, digit_sum_counter, additive_ds_counter)."""
    is_p = sieve_of_eratosthenes(N - 1)
    small_primes = set(p for p in range(2, 55) if is_p[p])
    count = 0
    ds_all = Counter()       # digit sums of ALL primes
    ds_additive = Counter()  # digit sums of additive primes only
    for p in range(2, N):
        if is_p[p]:
            d = digit_sum(p)
            ds_all[d] += 1
            if d in small_primes:
                count += 1
                ds_additive[d] += 1
    return count, is_p, ds_all, ds_additive

def additive_ratio_by_magnitude(is_p, max_exp=6):
    """Fraction of primes that are additive in [2, 10^k) for k=1..max_exp."""
    small_primes = set(p for p in range(2, 55) if is_p[p])
    exponents, ratios = [], []
    for k in range(1, max_exp + 1):
        bound = 10 ** k
        total = 0
        additive = 0
        for p in range(2, min(bound, len(is_p))):
            if is_p[p]:
                total += 1
                if digit_sum(p) in small_primes:
                    additive += 1
        exponents.append(k)
        ratios.append(additive / total if total else 0)
    return exponents, ratios

def digit_sum_histogram(ds_counter):
    """Sorted keys and values from a digit-sum counter."""
    keys = sorted(ds_counter.keys())
    vals = [ds_counter[k] for k in keys]
    return keys, vals

#  Verification
# Small known values: additive primes below 100
# 2,3,5,7,11,23,29,41,43,47,61,67,83,89 => 14
is_p_small = sieve_of_eratosthenes(99)
sp = set(p for p in range(2, 55) if is_p_small[p])
small_count = sum(1 for p in range(2, 100) if is_p_small[p] and digit_sum(p) in sp)
assert small_count == 14, f"Expected 14 additive primes < 100, got {small_count}"

#  Main computation
N = 10 ** 6
answer, is_p, ds_all, ds_additive = count_additive_primes(N)
print(answer)

#  Visualization — 4-panel figure