Problem 249: Prime Subset Sums

Project Euler Problem 249 Solution

EulerSolve provides an optimized solution for Project Euler Problem 249, Prime Subset Sums, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary Let $\mathcal{P}$ be the set of all primes below 5000. There are 669 such primes, so a brute-force search would have to reason about $2^{669}$ different subsets. For each subset $A \subseteq \mathcal{P}$, form the sum $\sum_{p \in A} p$. The task is to count how many subsets have a prime-valued sum and report that count modulo $10^{16}$. The count is over subsets, not over distinct sums. If two different subsets produce the same prime total, both must be counted. The empty subset contributes the sum 0, which is not prime, so it never enters the final answer, but it is still the correct base case for the recurrence. Mathematical Approach The natural state space is the set of all achievable subset sums. Since every input value is itself prime and may be chosen at most once, the problem is a counting version of 0/1 subset sum rather than an unbounded partition problem. The generating function for subset sums Order the input primes increasingly as $$p_1 < p_2 < \cdots < p_m,\qquad m=669.$$ If $$G(x)=\prod_{i=1}^{m}(1+x^{p_i}),$$ then choosing the term $1$ means “do not use $p_i$”, while choosing $x^{p_i}$ means “use $p_i$ once”....

Detailed mathematical approach

Problem Summary

Let $S=\{2,3,5,\dots,4999\}$ be the set of the $669$ primes below $5000$. We must count the subsets of $S$ whose elements add up to a prime number. This count is enormous, so only its rightmost $16$ digits are asked for, that is, its remainder modulo $10^{16}$. ^[1] Subsets with the same sum are counted separately: $\{5\}$ and $\{2,3\}$ both have the prime sum $5$ and count as two. The empty subset has the sum $0$, which is not prime. The statement gives no smaller example; for the four primes below $10$ the count is $7$, and for the eight primes below $20$ it is $65$, the two values that the C++ program checks first.

The whole solution in one sentence: we add the primes one at a time and keep, for every possible sum $s$ from $0$ to $1\,548\,136$, the number of subsets with sum $s$, because a subset either avoids the new prime or contains it, and at the end we add up these numbers over the prime values of $s$, keeping only the last $16$ digits throughout.

Mathematical Approach

Listing the subsets is hopeless, but their sums take only about $1.5$ million values, so we count the subsets of each sum and ask which sums are prime only at the end. All sums must be kept, not only the prime ones: $\{2,3,5\}$ has the composite sum $10$, and adding $7$ gives the prime sum $17$.

1. The four primes below 10

The $16$ subsets of $\{2,3,5,7\}$ have these sums: $0$ for the empty subset; $2$, $3$, $5$ and $7$ for the single primes; $5$, $7$, $9$, $8$, $10$ and $12$ for the pairs $\{2,3\}$, $\{2,5\}$, $\{2,7\}$, $\{3,5\}$, $\{3,7\}$ and $\{5,7\}$; $10$, $12$, $14$ and $15$ for the triples without $7$, $5$, $3$ and $2$ respectively; and $17$ for the whole set. Seven of these sums are prime: those of the four single primes, of $\{2,3\}$ with sum $5$, of $\{2,5\}$ with sum $7$ and of the whole set with sum $17$. The $16$ subsets have only $12$ different sums, and all we need is how many subsets have each sum.

2. Counting the subsets of each sum

For the real set, $2^{669}$ is a number with $202$ digits, but every sum lies between $0$ and $\sigma=1\,548\,136$, the sum of all primes in $S$. Write the primes in increasing order as $p_1=2$, $p_2=3$, …, $p_{669}=4999$, put $\sigma_0=0$ and $\sigma_k=p_1+\dots+p_k$, and let $c_k(s)$ be the number of subsets of $\{p_1,\dots,p_k\}$ with sum $s$. With no primes there is only the empty subset, so $c_0(0)=1$ and $c_0(s)=0$ for $s\ne0$.

Claim. For $k\ge1$ and every integer $s$,

$$c_k(s)=c_{k-1}(s)+c_{k-1}(s-p_k),$$

where counts at negative sums are $0$.

Proof. Split the subsets of $\{p_1,\dots,p_k\}$ with sum $s$ by whether they contain $p_k$. Those without $p_k$ are exactly the subsets of $\{p_1,\dots,p_{k-1}\}$ with sum $s$. Removing $p_k$ from one that contains it leaves a subset of $\{p_1,\dots,p_{k-1}\}$ with sum $s-p_k$, and adding $p_k$ back reverses this, so there are $c_{k-1}(s-p_k)$ of them. $\square$

The quantity asked for is therefore the sum of $c_{669}(q)$ over the primes $q\le\sigma$, modulo $10^{16}$. Filling a table of counts from smaller counts of the same kind, each computed once, is dynamic programming; for subset sums its cost is the number of elements times the number of possible sums. ^[2]

3. Worked example: the four primes below 10

Starting from $c_0(0)=1$, the claim gives the following counts; all counts not listed are $0$.

After $p_1=2$: $c_1(0)=c_1(2)=1$.
After $p_2=3$: $c_2(s)=1$ for $s=0,2,3,5$.
After $p_3=5$: $c_3(5)=2$, and $c_3(s)=1$ for $s=0,2,3,7,8,10$.
After $p_4=7$: $c_4(s)=2$ for $s=5,7,10,12$, and $c_4(s)=1$ for $s=0,2,3,8,9,14,15,17$.

For example, $c_3(5)=c_2(5)+c_2(0)$${}=1+1=2$ counts $\{2,3\}$ and $\{5\}$, and $c_4(12)=c_3(12)+c_3(5)$${}=0+2=2$ counts $\{5,7\}$ and $\{2,3,7\}$, both made by adding $7$ to a subset with sum $5$. The prime sums $2,3,5,7,17$ carry $1+1+2+2+1=7$ subsets, as found in Section 1.

4. The same counts as a polynomial

Multiply out $(1+x^2)(1+x^3)(1+x^5)(1+x^7)$. Every term of the expansion takes from each factor either $1$, leaving the prime out, or $x^p$, taking it in; so the terms match the subsets one to one, and the term of a subset is $x^s$ for its sum $s$. Collecting equal powers gives

$$\begin{aligned}&1+x^2+x^3+2x^5+2x^7+x^8+x^9\\&\quad+2x^{10}+2x^{12}+x^{14}+x^{15}+x^{17},\end{aligned}$$

with the counts $c_4(s)$ of Section 3 as coefficients. In the same way the coefficient of $x^s$ in $\prod_{p\in S}(1+x^p)$ is the number of subsets of $S$ with sum $s$, and multiplying by one more factor $1+x^{p_k}$ adds to every coefficient the one $p_k$ places lower, which is the claim of Section 2. A power series whose coefficients count objects is a generating function, and multiplying two of them convolves their coefficient sequences. ^[3] The product over all primes is the generating function of OEIS A000586, the number of partitions of $n$ into distinct primes, which begins $1,0,1,1,0,2,0,2,1,1,2$. ^[4] The first prime after $4999$ is $5003$, so every prime up to $5002$ lies in $S$, and $c_{669}(s)$ equals this sequence for all $s\le5002$.

5. One array, swept downwards

The programs keep one array $w$, called ways, with room for the sums $0$ to $\sigma$, instead of $669$ separate tables. Before the prime $p=p_k$ is added it holds $w(s)=c_{k-1}(s)$, and adding $p$ means

$$\begin{aligned}w(s+p)&\leftarrow w(s+p)+w(s),\\s&=\sigma_{k-1},\dots,1,0,\end{aligned}$$

in this order, from the largest sum down to $0$.

Claim. After this sweep, $w(t)=c_k(t)$ for every $t$.

Proof. The step for $s$ changes only the entry $s+p$, which lies above $s$. The entry $s$ can therefore be changed only by the step for $s-p$, which comes later in the decreasing order, so the step for $s$ reads the old value $c_{k-1}(s)$. Each entry $t$ is changed at most once, by the step for $t-p$, which adds $c_{k-1}(t-p)$; where this step does not exist, because $t-p\lt0$ or $t-p\gt\sigma_{k-1}$, that count is $0$ anyway. So $w(t)=c_{k-1}(t)+c_{k-1}(t-p)$, which is $c_k(t)$ by Section 2. $\square$

The order matters. After $2$ and $3$ the array holds $1$ at $0,2,3,5$. Sweeping upwards for $p=5$, the step for $s=0$ raises $w(5)$ to $2$, and the later step for $s=5$ adds this new value to $w(10)$, which ends at $2$ instead of $c_3(10)=1$: the extra "subset" is $\{5\}$ with $5$ added a second time. Run upwards over the whole array, the update would count sums in which each prime may be used any number of times; for items used at most once the loop must run downwards. ^[5] Starting the sweep at $\sigma_{k-1}$ rather than at the top of the array, $s=\sigma-p$, is safe because $c_{k-1}(s)=0$ for $s\gt\sigma_{k-1}$, and it reduces the work from $1\,034\,155\,517$ to $325\,967\,779$ steps.

6. Which sums occur

Replacing a subset $A$ by its complement $S\setminus A$ turns the sum $s$ into $\sigma-s$, so $c_{669}(s)=c_{669}(\sigma-s)$; in Section 3, with the total $17$, for example $c_4(5)=c_4(12)=2$. Every positive integer except $1$, $4$ and $6$ is a sum of distinct primes, a consequence of Richert's theorem, ^[6] and for our $S$ the exact counts show that the only sums from $0$ to $\sigma$ that do not occur are $1$, $4$, $6$ and, by the symmetry, $\sigma-6$, $\sigma-4$ and $\sigma-1$. A check of the reachable sums finds the same pattern, with $\sigma_k$ in place of $\sigma$, for every $k$ from $4$ to $669$, so the programs' zero test skips the addition in only $3\,997$ of the $325\,967\,779$ loop steps. The test looks at the stored remainders, and a run of the sweeps modulo $10^{16}$ gives the same $3\,997$, so no positive count visited there has the remainder $0$.

7. Only the last 16 digits

The counts are far too large for machine integers: the largest, $c_{669}(774\,068)$ at the middle sum $\sigma/2$, has $197$ digits. The programs, however, only add, and congruence modulo $m$ is compatible with addition: if $x\equiv x'$ and $y\equiv y'\pmod m$, then $x+y\equiv x'+y'\pmod m$. ^[7] So each count may be replaced by its remainder modulo $m=10^{16}$ after every addition. By induction over the sweeps, every entry of $w$ is then $c_k(s)\bmod m$, and the sum of these remainders over the prime sums, reduced once more, is the requested count modulo $m$. Two remainders are below $m$, so their sum is below $2m$, and subtracting $m$ once whenever the sum reaches $m$ brings it back into the range from $0$ to $m-1$. The largest intermediate value, below $2\cdot10^{16}$, is about $461$ times smaller than $2^{63}-1$, the largest signed $64$-bit integer.

8. The prime sums

The sieve of Eratosthenes finds all primes up to $\sigma$ at once. For $i=2,3,\dots$ while $i^2\le\sigma$, if $i$ has not been crossed out, it crosses out $i^2$, $i^2+i$, $i^2+2i$, …; the numbers from $2$ on that are never crossed out are the primes. Starting at $i^2$ suffices because a smaller multiple $ji$ with $2\le j\lt i$ has a prime factor $r\le j$ and was crossed out by $r$ already; stopping after $i=\lfloor\sqrt\sigma\rfloor=1244$ suffices because every composite number up to $\sigma$ has a prime factor at most $\sqrt\sigma$. ^[8] The sieve finds $117\,531$ primes up to $\sigma$, and the same sieve run below $5000$ produces $S$ itself.

For the eight primes below $20$, with the total $77$, the $21$ primes $q=2,3,5,\dots,73$ up to $77$ carry the counts $1,1,2,2,1,2,2,3,4,6,6,$$7,5,5,4,5,4,3,2,0,0$, whose sum is the checkpoint $65$. The last two are $0$ by Section 6, because $77-71=6$ and $77-73=4$ are not sums of distinct primes.

9. How large the count is

Choose each prime independently with probability $1/2$, which makes all $2^{669}$ subsets equally likely. Each prime then contributes $0$ or $p$ with equal chances, with the mean $p/2$ and the variance $p^2/4$; means add, and for independent choices variances add too. So the sum has the mean $\sigma/2=774\,068$ and the standard deviation $\tfrac12\sqrt{V}\approx35\,531$, where $V=\sum_{p\in S}p^2$, only $4.6$ per cent of the mean. Exactly half of the subsets have an odd sum, because adding or removing $3$ changes the parity and pairs the subsets up, and an even sum is prime only for $\{2\}$. By the prime number theorem the density of primes near $774\,068$ is about $1/\ln774\,068$, ^[9] so about $2/\ln774\,068\approx0.1475$ of the odd numbers there are prime. If the odd sums meet primes at this rate, about $1/\ln774\,068\approx0.07375$ of all subsets have a prime sum; this is a heuristic estimate, not a proven bound. The actual proportion, rounded to five decimals, is $0.07385$, so the count has $201$ digits.

How the Code Works

The C++ program ^[10] starts in main. parse_arguments fills an Options record whose defaults are prime_limit $=5000$ and modulo $=10^{16}$; it accepts --prime-limit=, --mod= and --skip-checkpoints, and an unknown argument, a limit below $2$ or a zero modulus ends the program with exit code $1$. Unless the checkpoints are skipped, run_checkpoints requires solve(10, …) to return $7$ and solve(20, …) to return $65$, and otherwise ends the program with exit code $2$. Then the result of solve(prime_limit, modulo) is printed.

primes_below(limit) is the sieve of Section 8 on a byte array and returns the primes in increasing order. solve adds them up to max_sum $=\sigma$, allocates ways as $\sigma+1$ unsigned 64-bit counts with ways[0] = 1, and for each prime runs the downward sweep of Section 5 from current_sum $=\sigma_{k-1}$ to $0$: it skips zero entries, adds ways[s] to ways[s + p] and subtracts modulo once if the result reaches it (Section 7), and then increases current_sum by $p$. A second byte sieve, is_prime_sum, marks the primes up to max_sum, and a last loop adds ways[s] for every prime $s$, reducing with % after each addition. The single subtraction is exact for every modulus up to $2^{63}$, as two remainders then add up to less than $2^{64}$; a larger --mod can wrap around. With --prime-limit=2 there are no primes, and is_prime_sum[1] = 0 writes past the end of a one-element array. Neither case affects the default run.

The Python and Java programs implement the same steps in one function solve(), with the limit $5000$ and the modulus $10^{16}$ fixed and without options or checkpoints. Python's sieves are bytearrays, ways is a list of Python integers, which have unlimited precision, and every update reduces with % instead of the conditional subtraction. Java uses boolean arrays, collects the primes in a List<Integer> while adding them to maxSum, and stores the counts in a long[] with the conditional subtraction of the C++ version; its signed 64-bit long reaches $2^{63}-1$, far above $2\cdot10^{16}$. ^[11]

Complexity and Verification

The sweeps make $\sum_{k=1}^{669}(\sigma_{k-1}+1)=325\,967\,779$ steps of constant cost, about $0.31$ times $669\,\sigma$, because the partial sums $\sigma_k$ grow roughly quadratically in $k$. For the $K$ primes below a general limit $L$, with the total $\sigma$, the time is $O(K\sigma)$: polynomial in $L$, but not in its number of digits, which is what pseudo-polynomial means. ^[2] The sieve up to $\sigma$ costs $O(\sigma\log\log\sigma)$ ^[8] and the final loop $O(\sigma)$. The memory is $O(\sigma)$: in C++ the $\sigma+1$ counts take $8$ bytes each, $12.4$ MB, and the second sieve $\sigma+1$ bytes, with a measured peak of $15.8$ MB; Python, whose list refers to separate integer objects, peaks at about $91$ MB. On the test machine (Apple M4 Max; clang++ -O2, OpenJDK 17, Python 3.11) the C++, Java and Python programs take about $0.08$, $1.2$ and $33$ seconds.

The C++ program checks itself with the two checkpoints $7$ and $65$. For this article, a brute force over all $2^k$ subsets of the $k$ primes below $10$, $20$, $30$, $40$, $50$ and $60$ gave $7$, $65$, $237$, $846$, $6159$ and $23484$, in agreement with the recurrence and with the C++ program run with --prime-limit. An exact computation of the whole table with big integers and no modulus confirmed that the counts add up to $2^{669}$, that $2^{668}$ subsets have an odd sum, the symmetry and the six missing sums of Section 6, the first $60$ terms of OEIS A000586 and the numbers of Sections 7 and 9. Reduced modulo $10^{16}$, the exact number of subsets with a prime sum equals the output of all three programs, which has $16$ digits, so no leading zeros are lost. These checks are not part of the published solutions.

Footnotes and References

The references below supply the background facts used above; the recurrence, the proof of the downward sweep and all counts and timings are derived or computed in this article.

Project Euler 249 — Prime Subset Sums. The official statement lets S = {2, 3, 5, …, 4999} be the set of primes less than 5000 and asks for the rightmost 16 digits of the number of subsets of S whose elements have a prime sum; it gives no smaller example.
Wikipedia — Subset sum problem. The problem is NP-complete, but dynamic programming solves it in pseudo-polynomial time, proportional to the number of elements times the range of possible sums, which is not polynomial in the input size in bits.
Wikipedia — Generating function. An ordinary generating function encodes a sequence as the coefficients of a power series, and multiplying two of them yields the discrete convolution, the Cauchy product, of their sequences.
OEIS A000586 — Number of partitions of n into distinct primes. The entry gives the generating function as the product of (1 + x^p) over all primes p and lists the terms 1, 0, 1, 1, 0, 2, 0, 2, 1, 1, 2, … of Section 4.
cp-algorithms — Knapsack Problem. In a one-dimensional array the 0-1 knapsack update must run in decreasing order of the capacity, so that each item is used at most once; the increasing order gives the complete knapsack, with unlimited repetitions.
ProofWiki — Number as Sum of Distinct Primes. Every positive integer n other than 1, 4 and 6 is a sum of distinct primes; the proof checks 7 to 19 by hand and then applies Richert's theorem and the Bertrand–Chebyshev theorem.
Wikipedia — Modular arithmetic. Congruence modulo m is compatible with addition, and a mod m denotes the unique integer r with 0 ≤ r < m and r ≡ a (mod m).
Wikipedia — Sieve of Eratosthenes. The multiples of each prime p are crossed out starting from p², the algorithm may stop once p² exceeds n, and it runs in O(n log log n) time.
Wikipedia — Prime number theorem. The theorem states π(x) ~ x / log x, and the article notes that the logarithmic integral suggests a density 1 / log t of the primes around t.
C++ — Euler249.cpp. The linked immutable revision contains parse_arguments, primes_below, solve and run_checkpoints; it is the implementation record for this article.
Java — Primitive Data Types gives the range −2⁶³ to 2⁶³ − 1 of long, and Python — numeric types states that integers have unlimited precision.

Problem 249 source code

C++

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

namespace {

using u64 = std::uint64_t;

struct Options {
    int prime_limit = 5000;
    u64 modulo = 10000000000000000ULL;
    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_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 = 0;
    for (char c : tail) {
        if (c < '0' || c > '9') {
            return false;
        }
        parsed = parsed * 10 + 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_int_after_prefix(arg, "--prime-limit=", options.prime_limit) ||
            parse_u64_after_prefix(arg, "--mod=", options.modulo)) {
            continue;
        }

        std::cerr << "Unknown argument: " << arg << '\n';
        return false;
    }

    return options.prime_limit >= 2 && options.modulo > 0;
}

std::vector<int> primes_below(const int limit) {
    std::vector<std::uint8_t> sieve(static_cast<std::size_t>(limit), 1);
    sieve[0] = 0;
    sieve[1] = 0;

    for (int i = 2; i * i < limit; ++i) {
        if (!sieve[static_cast<std::size_t>(i)]) {
            continue;
        }
        for (int j = i * i; j < limit; j += i) {
            sieve[static_cast<std::size_t>(j)] = 0;
        }
    }

    std::vector<int> primes;
    for (int i = 2; i < limit; ++i) {
        if (sieve[static_cast<std::size_t>(i)]) {
            primes.push_back(i);
        }
    }
    return primes;
}

u64 solve(const int prime_limit, const u64 modulo) {
    const std::vector<int> primes = primes_below(prime_limit);

    int max_sum = 0;
    for (int p : primes) {
        max_sum += p;
    }

    std::vector<u64> ways(static_cast<std::size_t>(max_sum + 1), 0);
    ways[0] = 1;

    int current_sum = 0;
    for (int p : primes) {
        for (int s = current_sum; s >= 0; --s) {
            const u64 add = ways[static_cast<std::size_t>(s)];
            if (add == 0) {
                continue;
            }
            u64& dst = ways[static_cast<std::size_t>(s + p)];
            dst += add;
            if (dst >= modulo) {
                dst -= modulo;
            }
        }
        current_sum += p;
    }

    std::vector<std::uint8_t> is_prime_sum(static_cast<std::size_t>(max_sum + 1), 1);
    is_prime_sum[0] = 0;
    is_prime_sum[1] = 0;
    for (int i = 2; i * i <= max_sum; ++i) {
        if (!is_prime_sum[static_cast<std::size_t>(i)]) {
            continue;
        }
        for (int j = i * i; j <= max_sum; j += i) {
            is_prime_sum[static_cast<std::size_t>(j)] = 0;
        }
    }

    u64 answer = 0;
    for (int s = 2; s <= max_sum; ++s) {
        if (!is_prime_sum[static_cast<std::size_t>(s)]) {
            continue;
        }
        answer += ways[static_cast<std::size_t>(s)];
        answer %= modulo;
    }

    return answer;
}

bool run_checkpoints() {
    if (solve(10, 10000000000000000ULL) != 7ULL) {
        std::cerr << "Checkpoint failed for primes below 10" << '\n';
        return false;
    }
    if (solve(20, 10000000000000000ULL) != 65ULL) {
        std::cerr << "Checkpoint failed for primes below 20" << '\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.prime_limit, options.modulo) << '\n';
    return 0;
}

Python

def solve():
    prime_limit = 5000
    modulo = 10**16

    sieve = bytearray(b'\x01' * prime_limit)
    sieve[0] = 0
    sieve[1] = 0
    for i in range(2, prime_limit):
        if i * i >= prime_limit:
            break
        if sieve[i]:
            for j in range(i * i, prime_limit, i):
                sieve[j] = 0
    primes = [i for i in range(2, prime_limit) if sieve[i]]

    max_sum = sum(primes)
    ways = [0] * (max_sum + 1)
    ways[0] = 1
    current_sum = 0
    for p in primes:
        for s in range(current_sum, -1, -1):
            add = ways[s]
            if add == 0:
                continue
            ways[s + p] = (ways[s + p] + add) % modulo
        current_sum += p

    is_prime_sum = bytearray(b'\x01' * (max_sum + 1))
    is_prime_sum[0] = 0
    is_prime_sum[1] = 0
    for i in range(2, max_sum + 1):
        if i * i > max_sum:
            break
        if is_prime_sum[i]:
            for j in range(i * i, max_sum + 1, i):
                is_prime_sum[j] = 0

    answer = 0
    for s in range(2, max_sum + 1):
        if is_prime_sum[s]:
            answer = (answer + ways[s]) % modulo
    return str(answer)

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

Java

import java.util.ArrayList;
import java.util.List;

public class Euler249 {
    public static String solve() {
        int primeLimit = 5000;
        long modulo = 10000000000000000L;

        boolean[] sieve = new boolean[primeLimit];
        for (int i = 2; i < primeLimit; i++)
            sieve[i] = true;

        for (int i = 2; i * i < primeLimit; ++i) {
            if (sieve[i]) {
                for (int j = i * i; j < primeLimit; j += i) {
                    sieve[j] = false;
                }
            }
        }

        List<Integer> primes = new ArrayList<>();
        int maxSum = 0;
        for (int i = 2; i < primeLimit; ++i) {
            if (sieve[i]) {
                primes.add(i);
                maxSum += i;
            }
        }

        long[] ways = new long[maxSum + 1];
        ways[0] = 1;

        int currentSum = 0;
        for (int p : primes) {
            for (int s = currentSum; s >= 0; --s) {
                long add = ways[s];
                if (add == 0)
                    continue;

                int dst = s + p;
                ways[dst] += add;
                if (ways[dst] >= modulo) {
                    ways[dst] -= modulo;
                }
            }
            currentSum += p;
        }

        boolean[] isPrimeSum = new boolean[maxSum + 1];
        for (int i = 2; i <= maxSum; i++)
            isPrimeSum[i] = true;

        for (int i = 2; i * i <= maxSum; ++i) {
            if (isPrimeSum[i]) {
                for (int j = i * i; j <= maxSum; j += i) {
                    isPrimeSum[j] = false;
                }
            }
        }

        long answer = 0;
        for (int s = 2; s <= maxSum; ++s) {
            if (isPrimeSum[s]) {
                answer += ways[s];
                answer %= modulo;
            }
        }

        return String.valueOf(answer);
    }

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