Problem 23: Non-abundant Sums

View on Project Euler

Project Euler Problem 23 Solution

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

Problem Summary For a positive integer \(n\), let \(s(n)\) denote the sum of its proper divisors. A number is called abundant when \(s(n) \gt n\). The task is to add all positive integers that cannot be written as \(a+b\) with both \(a\) and \(b\) abundant. The crucial reduction is given in the problem statement itself: every integer greater than \(28123\) can be written as the sum of two abundant numbers. Therefore the computation only needs the finite interval \(1 \le n \le L\) with \(L=28123\). Mathematical Approach The solution has two mathematical ingredients. First, classify every integer up to \(L\) by its proper-divisor sum. Second, build the set of all sums of two abundant numbers up to the same limit. Once those reachable sums are known, the answer is just the sum of the numbers that were never reached. Proper-Divisor Sums and the Definition of Abundance Define the aliquot sum by $$s(n)=\sum_{\substack{d \mid n \\ 1 \le d \lt n}} d.$$ Then \(n\) is deficient if \(s(n) \lt n\), perfect if \(s(n)=n\), and abundant if \(s(n) \gt n\). The first abundant number is \(12\), because $$s(12)=1+2+3+4+6=16 \gt 12.$$ This already gives one useful observation: no integer below \(24\) can be a sum of two abundant numbers, since the smallest possible such sum is \(12+12\)....

Detailed mathematical approach

Problem Summary

For a positive integer \(n\), let \(s(n)\) denote the sum of its proper divisors. A number is called abundant when \(s(n) \gt n\). The task is to add all positive integers that cannot be written as \(a+b\) with both \(a\) and \(b\) abundant.

The crucial reduction is given in the problem statement itself: every integer greater than \(28123\) can be written as the sum of two abundant numbers. Therefore the computation only needs the finite interval \(1 \le n \le L\) with \(L=28123\).

Mathematical Approach

The solution has two mathematical ingredients. First, classify every integer up to \(L\) by its proper-divisor sum. Second, build the set of all sums of two abundant numbers up to the same limit. Once those reachable sums are known, the answer is just the sum of the numbers that were never reached.

Proper-Divisor Sums and the Definition of Abundance

Define the aliquot sum by

$$s(n)=\sum_{\substack{d \mid n \\ 1 \le d \lt n}} d.$$

Then \(n\) is deficient if \(s(n) \lt n\), perfect if \(s(n)=n\), and abundant if \(s(n) \gt n\). The first abundant number is \(12\), because

$$s(12)=1+2+3+4+6=16 \gt 12.$$

This already gives one useful observation: no integer below \(24\) can be a sum of two abundant numbers, since the smallest possible such sum is \(12+12\).

Reducing the Problem to a Finite Sumset

Let

$$A=\{a\in\{1,2,\dots,L\}: s(a)\gt a\}$$

be the set of abundant numbers up to the search limit. We only care whether a value \(n\le L\) belongs to

$$E=(A+A)\cap\{1,2,\dots,L\}=\{a_i+a_j:\ a_i,a_j\in A,\ a_i+a_j\le L\}.$$

The requested total is therefore

$$\sum_{\substack{1 \le n \le L \\ n\notin E}} n.$$

That formula is the whole problem in compact form: build \(A\), mark every element of \(E\), and add the integers that remain outside \(E\).

Computing Every \(s(n)\) at Once

The implementations do not factor each number independently. Instead, they use a divisor sieve. For each possible proper divisor \(d\), add \(d\) to every multiple \(2d,3d,4d,\dots\le L\). When this process finishes, the entry at position \(n\) has received exactly the contribution of every proper divisor of \(n\), so it equals \(s(n)\).

The running time of this phase is controlled by the harmonic series:

$$\sum_{d=1}^{\lfloor L/2 \rfloor}\left\lfloor \frac{L}{d}\right\rfloor = O(L\log L).$$

For \(L=28123\), this preprocessing is small enough that a simple array-based sieve is entirely sufficient.

Why the Pair Scan Can Stop Early

After the sieve, the abundant numbers are collected in increasing order. Suppose we fix one abundant number \(a_i\) and then try partners \(a_j\) with \(j\ge i\). Because the list is sorted, the sums \(a_i+a_j\) grow as \(j\) grows. Therefore, once \(a_i+a_j\gt L\), every later partner also produces a sum above \(L\), so the inner scan can stop immediately.

This monotonicity is the key invariant behind the second phase. It means the algorithm never examines useless pairs beyond the search limit, even though it still conceptually explores the sumset \(A+A\).

Worked Example

Up to \(25\), the abundant numbers are \(12,18,20,24\). The first reachable sum is

$$12+12=24.$$

The next partner already gives

$$12+18=30 \gt 25,$$

so there is no reason to test \(12+20\) or \(12+24\) when the limit is \(25\). Thus \(24\) is marked as expressible, while \(25\) is never marked and must be included in the final answer.

How the Code Works

Phase 1: Build the Divisor-Sum Table

The C++, Python, and Java implementations allocate an integer table of length \(L+1\) and fill it by the divisor sieve. A single left-to-right pass over that table then identifies every abundant number from \(12\) to \(L\). Because this pass is increasing, the abundant list is automatically sorted.

Phase 2: Mark Every Reachable Sum

A Boolean table indexed by \(0,1,\dots,L\) records whether a number can be written as the sum of two abundant numbers. The implementations loop over abundant pairs with the second index starting at the first, so each unordered pair is handled once. If the sum stays within the limit, the Boolean entry is marked; if the sum exceeds the limit, the inner loop breaks immediately.

Phase 3: Accumulate the Missing Integers

The final pass simply adds all \(n\) with \(1 \le n \le L\) whose Boolean entry is still false. The implementations also verify two checkpoint cases before using the default bound, confirming that the partial answers for \(L=50\) and \(L=100\) are \(891\) and \(2766\).

Complexity Analysis

The divisor sieve costs \(O(L\log L)\). If \(k=|A|\) is the number of abundant numbers up to \(L\), then the marking phase is \(O(k^2)\) in the worst case, with practical savings from the early break in the sorted list. The final accumulation is \(O(L)\).

Space usage is \(O(L)\): one integer array for proper-divisor sums, one Boolean array for representability, and the list of abundant numbers. For \(L=28123\), these structures are easily small enough for a straightforward implementation in all three languages.

Footnotes and References

  1. Project Euler problem page: https://projecteuler.net/problem=23
  2. Abundant number: Wikipedia - Abundant number
  3. Aliquot sum: Wikipedia - Aliquot sum
  4. Additional background: MathWorld - Abundant Number

Problem 23 source code

C++

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

namespace {

using i64 = std::int64_t;

struct Options {
    int limit = 28123;
    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;
    }

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

    value = static_cast<int>(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, "--limit=", options.limit)) {
            continue;
        }

        std::cerr << "Unknown argument: " << arg << '\n';
        return false;
    }
    return options.limit >= 1;
}

std::vector<int> proper_divisor_sums(const int limit) {
    std::vector<int> sums(static_cast<std::size_t>(limit + 1), 0);
    for (int d = 1; d <= limit / 2; ++d) {
        for (int m = d * 2; m <= limit; m += d) {
            sums[static_cast<std::size_t>(m)] += d;
        }
    }
    return sums;
}

i64 solve(const int limit) {
    const std::vector<int> sums = proper_divisor_sums(limit);

    std::vector<int> abundant;
    abundant.reserve(static_cast<std::size_t>(limit / 2));
    for (int n = 12; n <= limit; ++n) {
        if (sums[static_cast<std::size_t>(n)] > n) {
            abundant.push_back(n);
        }
    }

    std::vector<bool> is_abundant_sum(static_cast<std::size_t>(limit + 1), false);
    for (std::size_t i = 0; i < abundant.size(); ++i) {
        for (std::size_t j = i; j < abundant.size(); ++j) {
            const int value = abundant[i] + abundant[j];
            if (value > limit) {
                break;
            }
            is_abundant_sum[static_cast<std::size_t>(value)] = true;
        }
    }

    i64 total = 0;
    for (int n = 1; n <= limit; ++n) {
        if (!is_abundant_sum[static_cast<std::size_t>(n)]) {
            total += n;
        }
    }

    return total;
}

bool run_checkpoints() {
    if (solve(50) != 891) {
        std::cerr << "Checkpoint failed for limit=50" << '\n';
        return false;
    }
    if (solve(100) != 2766) {
        std::cerr << "Checkpoint failed for limit=100" << '\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.limit) << '\n';
    return 0;
}

Python

import sys

def proper_divisor_sums(limit):
    sums = [0] * (limit + 1)
    for d in range(1, limit // 2 + 1):
        for m in range(d * 2, limit + 1, d):
            sums[m] += d
    return sums

def solve(limit):
    sums = proper_divisor_sums(limit)
    
    abundant = []
    for n in range(12, limit + 1):
        if sums[n] > n:
            abundant.append(n)
    
    is_abundant_sum = [False] * (limit + 1)
    n_abundant = len(abundant)
    for i in range(n_abundant):
        for j in range(i, n_abundant):
            value = abundant[i] + abundant[j]
            if value > limit:
                break
            is_abundant_sum[value] = True
    
    total = 0
    for n in range(1, limit + 1):
        if not is_abundant_sum[n]:
            total += n
    
    return total

def run_checkpoints():
    if solve(50) != 891:
        sys.stderr.write("Checkpoint failed for limit=50\n")
        return False
    if solve(100) != 2766:
        sys.stderr.write("Checkpoint failed for limit=100\n")
        return False
    return True

def parse_arguments(args):
    limit = 28123
    run_checkpoints_flag = True
    
    i = 1
    while i < len(args):
        arg = args[i]
        if arg == "--skip-checkpoints":
            run_checkpoints_flag = False
        elif arg.startswith("--limit="):
            try:
                value_str = arg[8:]
                if not value_str:
                    return None, None
                value = int(value_str)
                if value < 1:
                    return None, None
                limit = value
            except ValueError:
                return None, None
        else:
            sys.stderr.write(f"Unknown argument: {arg}\n")
            return None, None
        i += 1
    
    if limit < 1:
        return None, None
    
    return limit, run_checkpoints_flag

def main():
    args = sys.argv
    limit, run_checkpoints_flag = parse_arguments(args)
    
    if limit is None:
        sys.exit(1)
    
    if run_checkpoints_flag and not run_checkpoints():
        sys.exit(2)
    
    print(solve(limit))

if __name__ == "__main__":
    main()

Java

class Euler23 {

    private static class Options {
        int limit = 28123;
        boolean runCheckpoints = true;
    }

    private static boolean parseArguments(String[] args, Options options) {
        for (String arg : args) {
            if ("--skip-checkpoints".equals(arg)) {
                options.runCheckpoints = false;
                continue;
            }
            if (arg.startsWith("--limit=")) {
                try {
                    options.limit = Integer.parseInt(arg.substring(8));
                    if (options.limit < 1) {
                        System.err.println("Limit must be positive");
                        return false;
                    }
                    continue;
                } catch (NumberFormatException e) {
                    System.err.println("Invalid limit value: " + arg);
                    return false;
                }
            }
            System.err.println("Unknown argument: " + arg);
            return false;
        }
        return true;
    }

    private static int[] properDivisorSums(int limit) {
        int[] sums = new int[limit + 1];
        for (int d = 1; d <= limit / 2; d++) {
            for (int m = d * 2; m <= limit; m += d) {
                sums[m] += d;
            }
        }
        return sums;
    }

    private static long solve(int limit) {
        int[] sums = properDivisorSums(limit);

        java.util.ArrayList<Integer> abundant = new java.util.ArrayList<>();
        for (int n = 12; n <= limit; n++) {
            if (sums[n] > n) {
                abundant.add(n);
            }
        }

        boolean[] isAbundantSum = new boolean[limit + 1];
        for (int i = 0; i < abundant.size(); i++) {
            for (int j = i; j < abundant.size(); j++) {
                int value = abundant.get(i) + abundant.get(j);
                if (value > limit) {
                    break;
                }
                isAbundantSum[value] = true;
            }
        }

        long total = 0;
        for (int n = 1; n <= limit; n++) {
            if (!isAbundantSum[n]) {
                total += n;
            }
        }

        return total;
    }

    private static boolean runCheckpoints() {
        if (solve(50) != 891) {
            System.err.println("Checkpoint failed for limit=50");
            return false;
        }
        if (solve(100) != 2766) {
            System.err.println("Checkpoint failed for limit=100");
            return false;
        }
        return true;
    }

    public static void main(String[] args) {
        Options options = new Options();
        if (!parseArguments(args, options)) {
            System.exit(1);
        }

        if (options.runCheckpoints && !runCheckpoints()) {
            System.exit(2);
        }

        System.out.println(solve(options.limit));
    }
}