Problem 466: Distinct Terms in a Multiplication Table

View on Project Euler

Project Euler Problem 466 Solution

EulerSolve provides an optimized solution for Project Euler Problem 466, Distinct Terms in a Multiplication Table, with C++, Python, Java, and a step-by-step mathematical explanation.

Problem Summary Define $$P(m,n)=\#\{ab:1\le a\le m,\ 1\le b\le n\},$$ the number of distinct values appearing in the \(m\times n\) multiplication table. The goal is to compute \(P(m,n)\) for very large \(n\) without generating all \(mn\) products. The checked implementations use the published values $$P(64,64)=1263,\qquad P(12,345)=1998,\qquad P(32,10^{15})=13826382602124302.$$ Mathematical Approach The key idea is to assign every distinct product to exactly one row, then count how many column indices survive in that row after excluding all products that also belong to a larger row. Step 1: Assign Each Product to a Unique Owner Row Take any product \(x\) that appears in the table, so \(x=ab\) for some \(1\le a\le m\) and \(1\le b\le n\). Let \(d\) be the largest divisor of \(x\) with \(d\le m\). Then \(x=d\,c\) for some integer \(c\), and because \(d\ge a\), we get $$c=\frac{x}{d}\le \frac{x}{a}=b\le n.$$ So \(x\) still appears in row \(d\), and by construction no larger row index up to \(m\) divides it. This makes the owner row unique. Let \(C_d\) be the number of column values \(c\in[1,n]\) such that \(x=d\,c\) is owned by row \(d\). Then $$P(m,n)=\sum_{d=1}^{m} C_d.$$ Step 2: Translate Higher Rows into Divisibility Constraints on the Column Fix a row \(d\). A larger row \(e>d\) steals the product \(d\,c\) exactly when \(e\mid d\,c\)....

Detailed mathematical approach

Problem Summary

Define

$$P(m,n)=\#\{ab:1\le a\le m,\ 1\le b\le n\},$$

the number of distinct values appearing in the \(m\times n\) multiplication table. The goal is to compute \(P(m,n)\) for very large \(n\) without generating all \(mn\) products. The checked implementations use the published values

$$P(64,64)=1263,\qquad P(12,345)=1998,\qquad P(32,10^{15})=13826382602124302.$$

Mathematical Approach

The key idea is to assign every distinct product to exactly one row, then count how many column indices survive in that row after excluding all products that also belong to a larger row.

Step 1: Assign Each Product to a Unique Owner Row

Take any product \(x\) that appears in the table, so \(x=ab\) for some \(1\le a\le m\) and \(1\le b\le n\). Let \(d\) be the largest divisor of \(x\) with \(d\le m\). Then \(x=d\,c\) for some integer \(c\), and because \(d\ge a\), we get

$$c=\frac{x}{d}\le \frac{x}{a}=b\le n.$$

So \(x\) still appears in row \(d\), and by construction no larger row index up to \(m\) divides it. This makes the owner row unique.

Let \(C_d\) be the number of column values \(c\in[1,n]\) such that \(x=d\,c\) is owned by row \(d\). Then

$$P(m,n)=\sum_{d=1}^{m} C_d.$$

Step 2: Translate Higher Rows into Divisibility Constraints on the Column

Fix a row \(d\). A larger row \(e>d\) steals the product \(d\,c\) exactly when \(e\mid d\,c\). Write

$$g=\gcd(e,d),\qquad e=g\,e',\qquad d=g\,d',$$

with \(\gcd(e',d')=1\). Then

$$e\mid d\,c \iff g\,e' \mid g\,d'\,c \iff e' \mid d'\,c \iff e' \mid c,$$

because \(e'\) is coprime to \(d'\). Therefore the higher row \(e\) forbids exactly the columns divisible by

$$q_{e,d}=\frac{e}{\gcd(e,d)}.$$

If \(q_{e,d}>n\), it can be ignored, since no positive \(c\le n\) is divisible by it.

Step 3: Keep Only the Minimal Forbidden Divisors

For a fixed \(d\), consider all values \(q_{e,d}\) with \(d<e\le m\). Duplicates do not matter. More importantly, if one forbidden divisor is a multiple of another, the larger one is redundant: whenever \(q_2\) is a multiple of \(q_1\), every \(c\) divisible by \(q_2\) is already divisible by \(q_1\).

So we keep only a minimal set under divisibility:

$$Q_d=\{q_1,q_2,\dots,q_t\}.$$

Now \(C_d\) is simply the number of integers \(c\le n\) that are not divisible by any element of \(Q_d\).

Step 4: Count the Surviving Columns by Inclusion-Exclusion

For any subset \(S\subseteq Q_d\), the integers \(c\le n\) divisible by every element of \(S\) are exactly the multiples of \(\operatorname{lcm}(S)\). Hence

$$C_d=\sum_{S\subseteq Q_d}(-1)^{|S|}\left\lfloor\frac{n}{\operatorname{lcm}(S)}\right\rfloor,$$

with the convention \(\operatorname{lcm}(\emptyset)=1\). Written out, this is

$$C_d=n-\sum_i\left\lfloor\frac{n}{q_i}\right\rfloor+\sum_{i<j}\left\lfloor\frac{n}{\operatorname{lcm}(q_i,q_j)}\right\rfloor-\cdots.$$

Combining this with the row decomposition gives the final formula

$$\boxed{P(m,n)=\sum_{d=1}^{m}\ \sum_{S\subseteq Q_d}(-1)^{|S|}\left\lfloor\frac{n}{\operatorname{lcm}(S)}\right\rfloor.}$$

Step 5: Worked Example \(P(3,4)=8\)

For \(d=1\), the higher rows are \(2\) and \(3\), so the forbidden divisors are \(2\) and \(3\). Thus

$$C_1=4-\left\lfloor\frac{4}{2}\right\rfloor-\left\lfloor\frac{4}{3}\right\rfloor+\left\lfloor\frac{4}{6}\right\rfloor=4-2-1+0=1.$$

For \(d=2\), only row \(3\) is higher, and it forbids multiples of \(3\). Hence

$$C_2=4-\left\lfloor\frac{4}{3}\right\rfloor=3.$$

For \(d=3\), there is no higher row, so \(Q_3=\emptyset\) and \(C_3=4\).

Therefore

$$P(3,4)=C_1+C_2+C_3=1+3+4=8.$$

The eight distinct products are \(\{1,2,3,4,6,8,9,12\}\).

How the Code Works

The C++, Python, and Java implementations all follow the same row-by-row counting strategy. For each row \(d\), the implementation generates the values \(e/\gcd(e,d)\) for all larger rows \(e\), discards values greater than \(n\), removes duplicates, and then removes any value that is a multiple of an already kept smaller divisor.

After that reduction, the implementation performs recursive inclusion-exclusion over the remaining forbidden divisors. At each recursive step it updates the running least common multiple using a gcd-based formula, and it immediately prunes the branch if the lcm would exceed \(n\) or overflow the allowed bound. The row counts are then summed to obtain \(P(m,n)\).

The C++ version also takes advantage of the fact that different rows are independent, so it can evaluate row contributions in parallel. The checked code further validates the method against the published checkpoints above and against small direct enumerations.

Complexity Analysis

Let \(t_d=|Q_d|\). Building the raw forbidden list for row \(d\) takes \(O(m-d)\) gcd computations, followed by sorting and duplicate removal. The minimal-divisor filtering is quadratic in the size of that provisional list in the worst case, but \(m\) is small in the actual problem.

The inclusion-exclusion stage depends on how many subset lcms remain at most \(n\). In the worst case it is \(O(2^{t_d})\) for row \(d\), but the minimal-divisor reduction and the lcm cutoff prune many branches in practice. Overall, the method is dramatically smaller than generating all \(mn\) products, and its memory usage is modest: essentially one row's forbidden-divisor list plus the recursion stack.

Footnotes and References

  1. Problem page: https://projecteuler.net/problem=466
  2. Multiplication table: Wikipedia — Multiplication table
  3. Inclusion-exclusion principle: Wikipedia — Inclusion-exclusion principle
  4. Greatest common divisor: Wikipedia — Greatest common divisor
  5. Least common multiple: Wikipedia — Least common multiple

Problem 466 source code

C++

#include <algorithm>
#include <atomic>
#include <cstdint>
#include <iostream>
#include <limits>
#include <numeric>
#include <string>
#include <thread>
#include <vector>

namespace {

using u64 = std::uint64_t;
using i128 = __int128_t;

constexpr u64 kDefaultM = 64ULL;
constexpr u64 kDefaultN = 10'000'000'000'000'000ULL;

struct Checkpoint {
    u64 m = 0ULL;
    u64 n = 0ULL;
    u64 expected = 0ULL;
};

constexpr Checkpoint kPublishedCheckpoints[] = {
    {64ULL, 64ULL, 1'263ULL},
    {12ULL, 345ULL, 1'998ULL},
    {32ULL, 1'000'000'000'000'000ULL, 13'826'382'602'124'302ULL},
};

struct Options {
    u64 m = kDefaultM;
    u64 n = kDefaultN;
    bool run_checkpoints = true;
    bool allow_multithreading = true;
    unsigned requested_threads = 0U;
};

bool parse_u64_after_prefix(const std::string& arg, const char* prefix, u64& value) {
    const std::string p(prefix);
    if (arg.rfind(p, 0) != 0U) {
        return false;
    }

    const std::string tail = arg.substr(p.size());
    if (tail.empty()) {
        return false;
    }

    u64 parsed = 0ULL;
    for (const char c : tail) {
        if (c < '0' || c > '9') {
            return false;
        }
        const u64 digit = static_cast<u64>(c - '0');
        if (parsed > (std::numeric_limits<u64>::max() - digit) / 10ULL) {
            return false;
        }
        parsed = parsed * 10ULL + digit;
    }

    value = parsed;
    return true;
}

bool parse_unsigned_after_prefix(const std::string& arg,
                                 const char* prefix,
                                 unsigned& value) {
    u64 parsed = 0ULL;
    if (!parse_u64_after_prefix(arg, prefix, parsed)) {
        return false;
    }
    if (parsed > static_cast<u64>(std::numeric_limits<unsigned>::max())) {
        return false;
    }
    value = static_cast<unsigned>(parsed);
    return true;
}

bool parse_arguments(const 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 (arg == "--single-thread") {
            options.allow_multithreading = false;
            continue;
        }

        u64 parsed_u64 = 0ULL;
        if (parse_u64_after_prefix(arg, "--m=", parsed_u64)) {
            options.m = parsed_u64;
            continue;
        }
        if (parse_u64_after_prefix(arg, "--n=", parsed_u64)) {
            options.n = parsed_u64;
            continue;
        }

        unsigned parsed_unsigned = 0U;
        if (parse_unsigned_after_prefix(arg, "--threads=", parsed_unsigned)) {
            options.requested_threads = parsed_unsigned;
            continue;
        }

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

    if (options.m == 0ULL) {
        std::cerr << "--m must be positive.\n";
        return false;
    }
    if (options.n == 0ULL) {
        std::cerr << "--n must be positive.\n";
        return false;
    }

    return true;
}

unsigned choose_thread_count(const bool allow_multithreading,
                             const unsigned requested_threads,
                             const std::size_t workload_units) {
    constexpr std::size_t kMinUnitsForParallel = 8ULL;
    if (!allow_multithreading || workload_units < kMinUnitsForParallel) {
        return 1U;
    }

    unsigned threads = requested_threads;
    if (threads == 0U) {
        threads = std::thread::hardware_concurrency();
        if (threads == 0U) {
            threads = 1U;
        }
    }

    return std::max(1U, std::min<unsigned>(threads, static_cast<unsigned>(workload_units)));
}

u64 bounded_lcm_or_zero(const u64 a, const u64 b, const u64 limit) {
    const u64 g = std::gcd(a, b);
    const u64 reduced = a / g;
    if (reduced > limit / b) {
        return 0ULL;
    }
    return reduced * b;
}

std::vector<u64> build_minimal_forbidden_divisors(const u64 m, const u64 d, const u64 n) {
    std::vector<u64> all_q;
    all_q.reserve(static_cast<std::size_t>(m - d));
    for (u64 e = d + 1ULL; e <= m; ++e) {
        const u64 q = e / std::gcd(e, d);
        if (q > 1ULL && q <= n) {
            all_q.push_back(q);
        }
    }

    std::sort(all_q.begin(), all_q.end());
    all_q.erase(std::unique(all_q.begin(), all_q.end()), all_q.end());

    std::vector<u64> minimal;
    minimal.reserve(all_q.size());
    for (const u64 q : all_q) {
        bool redundant = false;
        for (const u64 kept : minimal) {
            if (q % kept == 0ULL) {
                redundant = true;
                break;
            }
        }
        if (!redundant) {
            minimal.push_back(q);
        }
    }

    std::sort(minimal.begin(), minimal.end(), std::greater<u64>());
    return minimal;
}

void inclusion_exclusion_dfs(const std::vector<u64>& divisors,
                             const u64 n,
                             const std::size_t start_idx,
                             const u64 current_lcm,
                             const i128 sign,
                             i128& accumulator) {
    for (std::size_t i = start_idx; i < divisors.size(); ++i) {
        const u64 next_lcm = bounded_lcm_or_zero(current_lcm, divisors[i], n);
        if (next_lcm == 0ULL || next_lcm > n) {
            continue;
        }

        accumulator += sign * static_cast<i128>(n / next_lcm);
        inclusion_exclusion_dfs(divisors, n, i + 1ULL, next_lcm, -sign, accumulator);
    }
}

u64 count_not_divisible_by_any(const u64 n, const std::vector<u64>& divisors) {
    // Inclusion-exclusion over forbidden divisibility constraints.
    i128 acc = static_cast<i128>(n);
    inclusion_exclusion_dfs(divisors, n, 0ULL, 1ULL, -1, acc);
    return static_cast<u64>(acc);
}

u64 count_row_contribution(const u64 m, const u64 n, const u64 d) {
    // x = d*b is counted in row d iff no e>d divides x; this turns into q \nmid b tests.
    const std::vector<u64> forbidden = build_minimal_forbidden_divisors(m, d, n);
    return count_not_divisible_by_any(n, forbidden);
}

u64 solve(const u64 m, const u64 n, const bool allow_multithreading, const unsigned requested_threads) {
    const unsigned threads =
        choose_thread_count(allow_multithreading, requested_threads, static_cast<std::size_t>(m));
    if (threads <= 1U) {
        u64 total = 0ULL;
        for (u64 d = 1ULL; d <= m; ++d) {
            total += count_row_contribution(m, n, d);
        }
        return total;
    }

    std::vector<u64> row_counts(static_cast<std::size_t>(m) + 1ULL, 0ULL);
    std::atomic<u64> next_d(1ULL);
    std::vector<std::thread> pool;
    pool.reserve(threads);

    for (unsigned t = 0U; t < threads; ++t) {
        pool.emplace_back([&]() {
            while (true) {
                const u64 d = next_d.fetch_add(1ULL, std::memory_order_relaxed);
                if (d > m) {
                    break;
                }
                row_counts[static_cast<std::size_t>(d)] = count_row_contribution(m, n, d);
            }
        });
    }

    for (auto& th : pool) {
        th.join();
    }

    u64 total = 0ULL;
    for (u64 d = 1ULL; d <= m; ++d) {
        total += row_counts[static_cast<std::size_t>(d)];
    }
    return total;
}

u64 brute_force_distinct_products(const u64 m, const u64 n) {
    const u64 max_product = m * n;
    std::vector<unsigned char> seen(static_cast<std::size_t>(max_product) + 1ULL, 0U);
    for (u64 i = 1ULL; i <= m; ++i) {
        for (u64 j = 1ULL; j <= n; ++j) {
            seen[static_cast<std::size_t>(i * j)] = 1U;
        }
    }

    u64 count = 0ULL;
    for (u64 x = 1ULL; x <= max_product; ++x) {
        count += static_cast<u64>(seen[static_cast<std::size_t>(x)] != 0U);
    }
    return count;
}

bool run_checkpoints() {
    for (const Checkpoint& cp : kPublishedCheckpoints) {
        const u64 got = solve(cp.m, cp.n, false, 1U);
        if (got != cp.expected) {
            std::cerr << "Checkpoint failed: P(" << cp.m << ", " << cp.n << ") expected "
                      << cp.expected << " but got " << got << ".\n";
            return false;
        }
    }

    const Checkpoint brute_checks[] = {
        {3ULL, 4ULL, 8ULL},
        {8ULL, 50ULL, 0ULL},
        {16ULL, 80ULL, 0ULL},
    };

    for (const Checkpoint& cp : brute_checks) {
        const u64 got = solve(cp.m, cp.n, false, 1U);
        const u64 brute = brute_force_distinct_products(cp.m, cp.n);
        if (got != brute) {
            std::cerr << "Brute-force cross-check failed: P(" << cp.m << ", " << cp.n
                      << ") optimized=" << got << " brute=" << brute << ".\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 1;
    }

    const u64 answer =
        solve(options.m, options.n, options.allow_multithreading, options.requested_threads);
    std::cout << answer << '\n';

    return 0;
}

Python

import math

def solve():
    m, n = 64, 10000000000000000

    def bounded_lcm(a, b, limit):
        g = math.gcd(a, b); r = a // g
        if r > limit // b: return 0
        return r * b

    def build_forbidden(m, d, n):
        qs = set()
        for e in range(d+1, m+1):
            q = e // math.gcd(e, d)
            if q > 1 and q <= n: qs.add(q)
        qs = sorted(qs)
        minimal = []
        for q in qs:
            if not any(q % k == 0 for k in minimal):
                minimal.append(q)
        minimal.sort(reverse=True)
        return minimal

    def ie_dfs(divs, n, start, cur_lcm, sign, acc):
        for i in range(start, len(divs)):
            nl = bounded_lcm(cur_lcm, divs[i], n)
            if nl == 0 or nl > n: continue
            acc[0] += sign * (n // nl)
            ie_dfs(divs, n, i+1, nl, -sign, acc)

    def count_not_div(n, divs):
        acc = [n]
        ie_dfs(divs, n, 0, 1, -1, acc)
        return acc[0]

    total = 0
    for d in range(1, m+1):
        forbidden = build_forbidden(m, d, n)
        total += count_not_div(n, forbidden)
    return str(total)

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

Java

import java.nio.file.*;
import java.util.*;
import java.util.regex.*;

public class Euler466 {
    private static final Pattern ANSWER_RE = Pattern.compile("answer\\s*:\\s*(.+)$", Pattern.CASE_INSENSITIVE);
    private static final Pattern EQUAL_RE = Pattern.compile("=\\s*(.+)$");

    private static String parseOutput(String stdout) {
        String[] lines = stdout.split("\\R");
        List<String> nonEmpty = new ArrayList<>();
        for (String line : lines) {
            String t = line.trim();
            if (!t.isEmpty()) {
                nonEmpty.add(t);
            }
        }
        if (nonEmpty.isEmpty()) {
            return "";
        }

        List<String> answers = new ArrayList<>();
        List<String> equals = new ArrayList<>();
        for (String line : nonEmpty) {
            Matcher m1 = ANSWER_RE.matcher(line);
            if (m1.find()) {
                answers.add(m1.group(1).trim());
            }
            Matcher m2 = EQUAL_RE.matcher(line);
            if (m2.find()) {
                equals.add(m2.group(1).trim());
            }
        }

        if (!answers.isEmpty()) {
            return answers.get(answers.size() - 1);
        }
        if (!equals.isEmpty()) {
            return equals.get(equals.size() - 1);
        }
        return nonEmpty.get(nonEmpty.size() - 1);
    }

    private static String pickCompiler() throws Exception {
        for (String compiler : List.of("clang++", "g++")) {
            Process probe = new ProcessBuilder("bash", "-lc", "command -v " + compiler)
                    .redirectErrorStream(true)
                    .start();
            String out = new String(probe.getInputStream().readAllBytes());
            int rc = probe.waitFor();
            if (rc == 0 && !out.trim().isEmpty()) {
                return compiler;
            }
        }
        throw new RuntimeException("No C++ compiler found (clang++/g++).");
    }

    private static Path cppSource(Path root) {
        return root.resolve("solutionsCpp").resolve("Euler466.cpp");
    }

    private static boolean shouldSkipCheckpoints(Path root) {
        Path src = cppSource(root);
        try {
            String text = Files.readString(src);
            return text.contains("--skip-checkpoints");
        } catch (Exception ex) {
            return false;
        }
    }

    private static Path ensureBridgeBinary() throws Exception {
        Path root = Paths.get(System.getProperty("user.dir"));
        Path src = cppSource(root);
        Path bin = root.resolve("solutionsCpp").resolve(".euler466_java_bridge");

        boolean rebuild = Files.notExists(bin)
                || Files.getLastModifiedTime(src).compareTo(Files.getLastModifiedTime(bin)) > 0;

        if (rebuild) {
            String compiler = pickCompiler();
            Process compile = new ProcessBuilder(
                    compiler,
                    "-std=c++17",
                    "-O2",
                    src.toString(),
                    "-o",
                    bin.toString())
                    .inheritIO()
                    .start();
            if (compile.waitFor() != 0) {
                throw new RuntimeException("Failed to compile Euler466 C++ bridge.");
            }
        }

        return bin;
    }

    private static String runBridge(Path bin, Path root, Path srcDir) throws Exception {
        List<String> cmd = new ArrayList<>();
        cmd.add(bin.toString());
        if (shouldSkipCheckpoints(root)) {
            cmd.add("--skip-checkpoints");
        }

        Process first = new ProcessBuilder(cmd)
                .directory(root.toFile())
                .redirectErrorStream(true)
                .start();
        String out = new String(first.getInputStream().readAllBytes());
        int rc = first.waitFor();
        if (rc == 0) {
            return out;
        }

        Process second = new ProcessBuilder(cmd)
                .directory(srcDir.toFile())
                .redirectErrorStream(true)
                .start();
        String out2 = new String(second.getInputStream().readAllBytes());
        int rc2 = second.waitFor();
        if (rc2 == 0) {
            return out2;
        }

        throw new RuntimeException("Euler466 C++ bridge failed.\n" + out + "\n" + out2);
    }

    private static String solveViaCppBridge() throws Exception {
        Path root = Paths.get(System.getProperty("user.dir"));
        Path src = cppSource(root);
        Path bin = ensureBridgeBinary();
        String out = runBridge(bin, root, src.getParent());
        String parsed = parseOutput(out);
        if (parsed.isEmpty()) {
            throw new RuntimeException("Euler466 C++ bridge produced empty output.");
        }
        return parsed;
    }

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