Problem 268: Counting Numbers with at Least Four Distinct Prime Factors Less Than 100

Mathematical Analysis

Setup

The primes less than 100 are: 2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, 53, 59, 61, 67, 71, 73, 79, 83, 89, 97

There are 25 such primes.

Inclusion-Exclusion Principle

Let A_i be the set of numbers below N = 10^16 that are divisible by the i-th prime p_i.

We want |A_{i1} ∩ A_{i2} ∩ A_{i3} ∩ A_{i4} ∩ …| for at least 4 primes.

Using inclusion-exclusion on “at least 4”:

Count = sum_{k=4}^{25} (-1)^(k-4) * C(k-1, 3) * S_k

where S_k = sum over all k-element subsets {p_{i1}, …, p_{ik}} of floor(N / (p_{i1} * … * p_{ik}))

Derivation

Let f(k) = sum over all k-subsets of primes, of floor(N / product).

By inclusion-exclusion, the count of numbers divisible by at least 4 distinct primes from our set is:

Count = sum_{k=4}^{25} (-1)^(k-4) * C(k-1, 3) * f(k)

This follows from the standard inclusion-exclusion identity: to count elements in at least m sets out of possible sets, we use:

|at least m| = sum_{k=m}^{n} (-1)^(k-m) * C(k, m) …

Wait, let’s be more precise. Let S_k be the k-th elementary symmetric sum (sum over k-subsets of floor(N/product)).

The number of integers below N divisible by at least r distinct primes from the set equals:

sum_{k=r}^{25} (-1)^{k-r} * C(k-1, r-1) * S_k

For r = 4:

Count = sum_{k=4}^{25} (-1)^{k-4} * C(k-1, 3) * S_k

Computation

We enumerate all subsets of size k (for k = 4 to 25) of the 25 primes. For each subset, if the product exceeds N, the floor term is 0, so we skip it. Since the primes grow, many large subsets will have product > N.

The product of the 4 smallest primes is 235*7 = 210, and 10^16 / 210 is large. But the product of many primes quickly exceeds 10^16. For instance, the product of all 25 primes is enormous.

We enumerate subsets using DFS/recursion, pruning when the product exceeds N.

Complexity

• Worst case: all C(25, k) subsets for each k, but pruning makes this feasible
• The number of subsets with product <= 10^16 is manageable (on the order of millions)

Correctness

Theorem. The method described above computes exactly the quantity requested in the problem statement.

Proof. The preceding analysis identifies the admissible objects and derives the formula, recurrence, or exhaustive search carried out by the algorithm. The computation evaluates exactly that specification, so every valid contribution is included once and no invalid contribution is counted. Therefore the returned value is the required answer. $\square$

Answer

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

int main() {
    long long N = 10000000000000000LL; // 10^16
    // but the problem says "below 10^16", so we count numbers in [1, 10^16 - 1]
    // floor((N-1)/product) for divisibility count... actually floor((10^16-1)/p) = floor(10^16/p) - (p divides 10^16 ? 1 : 0)
    // Simpler: count of multiples of p in [1, N-1] = floor((N-1)/p)
    // But for large N, floor((10^16 - 1)/p) = floor(10^16/p) - 1 if p | 10^16, else floor(10^16/p)
    // Actually it's simpler: numbers below N = 10^16 means [1, 10^16 - 1]
    // count of multiples of d in [1, M] = floor(M/d)
    // So we use M = 10^16 - 1

    long long M = N - 1;

    vector<int> primes = {2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47,
                          53, 59, 61, 67, 71, 73, 79, 83, 89, 97};
    int np = primes.size(); // 25

    // Inclusion-exclusion: count of numbers with at least 4 distinct prime factors from the list
    // = sum_{k=4}^{25} (-1)^{k-4} * C(k-1, 3) * S_k
    // where S_k = sum over k-subsets of floor(M / product)

    // We'll enumerate subsets via recursion with pruning
    // For each subset of size k with product p, contribute (-1)^{k-4} * C(k-1,3) * floor(M/p)

    long long answer = 0;

    // Recursive enumeration
    // We process primes in order, building subsets
    function<void(int, int, long long)> solve = [&](int idx, int cnt, long long prod) {
        if (cnt >= 4) {
            // Contribute to answer
            long long term = M / prod;
            int sign = ((cnt - 4) % 2 == 0) ? 1 : -1;
            // C(cnt-1, 3)
            long long binom = 1;
            for (int i = 0; i < 3; i++) {
                binom *= (cnt - 1 - i);
                binom /= (i + 1);
            }
            answer += sign * binom * term;
        }

        for (int i = idx; i < np; i++) {
            // Check if product would overflow or exceed M
            if (prod > M / primes[i]) break; // primes are sorted, so all further products are larger
            solve(i + 1, cnt + 1, prod * primes[i]);
        }
    };

    solve(0, 0, 1);

    cout << answer << endl;
    return 0;
}

"""
Problem 268: Counting numbers with at least four distinct prime factors < 100

Count numbers below 10^16 divisible by at least 4 distinct primes less than 100.

Approach: Inclusion-exclusion.
For r = 4, the count of numbers with at least r distinct prime factors from a set is:
  sum_{k=r}^{25} (-1)^{k-r} * C(k-1, r-1) * S_k
where S_k = sum over all k-subsets of primes of floor(M / product).
"""

from math import comb
from itertools import combinations

def solve():
    M = 10**16 - 1  # numbers below 10^16 means [1, 10^16 - 1]

    primes = [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47,
              53, 59, 61, 67, 71, 73, 79, 83, 89, 97]

    # Recursive enumeration with pruning
    answer = 0

    def dfs(idx, cnt, prod):
        nonlocal answer
        if cnt >= 4:
            term = M // prod
            sign = 1 if (cnt - 4) % 2 == 0 else -1
            binom = comb(cnt - 1, 3)
            answer += sign * binom * term

        for i in range(idx, len(primes)):
            if prod * primes[i] > M:
                break
            dfs(i + 1, cnt + 1, prod * primes[i])

    dfs(0, 0, 1)
    print(answer)

solve()