Non-Abundant Sums

January 28, 2020

It’s been a long time since we did an exercise from Project Euler; here is number 23:

A perfect number is a number for which the sum of its proper divisors is exactly equal to the number. For example, the sum of the proper divisors of 28 would be 1 + 2 + 4 + 7 + 14 = 28, which means that 28 is a perfect number.

A number n is called deficient if the sum of its proper divisors is less than n and it is called abundant if this sum exceeds n.

As 12 is the smallest abundant number, 1 + 2 + 3 + 4 + 6 = 16, the smallest number that can be written as the sum of two abundant numbers is 24. By mathematical analysis, it can be shown that all integers greater than 28123 can be written as the sum of two abundant numbers. However, this upper limit cannot be reduced any further by analysis even though it is known that the greatest number that cannot be expressed as the sum of two abundant numbers is less than this limit.

Find the sum of all the positive integers which cannot be written as the sum of two abundant numbers.

Your task is to solve Project Euler 23; in the spirit of Project Euler, please do not display the solution. When you are finished, you are welcome to read or run a suggested solution, or to post your own solution or discuss the exercise in the comments below.

Pages: 1 2

5 Responses to “Non-Abundant Sums”

  1. Steve said

    Link to read does not work, and there is no run link.

  2. matthew said

    Here’s a Python solution that makes use of the multiplicativity of divisor sum:

    import sympy
    
    def dsum(N):
        primes = sympy.primerange(2,N)
        a = [0]+[1]*(N-1)
        b = [1]*N
        for p in primes:
            pp = p
            while pp < N:
                for i in range(pp,N,pp):
                    b[i] += pp
                pp *= p
            for i in range(p,N,p):
                a[i] *= b[i]
                b[i] = 1
        return a
    
    def abundant(N):
        s = dsum(N)
        return [n for n in range(N) if s[n] > 2*n]
    
    N = 28124
    a = abundant(N)
    b = [1]*N
    for p in a:
        lim = min(p,N-p-1)
        for q in a:
            if q > lim: break
            b[p+q] = 0
    print(sum(b)) # 1457 (counting from zero)
    
  3. matthew said

    That prints the total number of non-abundant sum numbers (like the Praxis solution at https://programmingpraxis.com/2020/01/28/non-abundant-sums/2/, except mine includes 0).

  4. Zack said

    Here is my take on this in Julia: https://pastebin.com/9XL9fWs5

    Total run time (including the time-consuming printing of the numbers): 44.5 sec.
    Memory usage: 59.04 MB

    Naturally, this code could be improved, something I may venture another time.

    Cheers for the intriguing problem! Definitely saved me a lot of time going through the PE archive and skimming through all the overly mathematical problems there!

  5. matthew said

    Reusing some code from https://programmingpraxis.com/2016/12/20/highly-abundant-numbers/, here’s a C++ version of that Python solution. Once again, just the count of non-sum numbers is printed:

    #include <stdio.h>
    #include <stdlib.h>
    #include <numeric>
    #include <vector>
     
    int main(int argc, char *argv[]) {
      int N = (argc == 1) ? 1000 : strtoul(argv[1],NULL,0);
      // a contains products of prime power sums
      // b contains prime power sum for current prime
      // initialize divisor sums to 1
      std::vector<int> a(N,1), b(N,1);
      a[0] = 0; // dsum(0) = 0
      for (int p = 2; p < N; p++) {
        // for each prime < N
        if (a[p] == 1) {
          // for each prime power: beware overflow here
          for (long pp = p; pp < N; pp *= p) {
            // add to sum for all multiples
            for (int i = pp; i < N; i += pp) {
              b[i] += pp;
            }
          }
          // for each multiple of p
          for (int i = p; i < N; i += p) {
            a[i] *= b[i]; // add in the prime power sum
            b[i] = 1;     // reset for future use
          }
        }
      }
    
      // Cross off sums in b, which is initially all ones.
      for (int i = 0; i < N; i++) {
        if (a[i] > 2*i) {
          for (int j = 0; j <= i && i+j < N; j++) {
            if (a[j] > 2*j) {
              b[i+j] = 0;
            }
          }
        }
      }
      
      printf("%d\n", std::accumulate(b.begin(),b.end(),0)); // 1457 for N = 30000
    }
    

Leave a Reply

Fill in your details below or click an icon to log in:

WordPress.com Logo

You are commenting using your WordPress.com account. Log Out /  Change )

Google photo

You are commenting using your Google account. Log Out /  Change )

Twitter picture

You are commenting using your Twitter account. Log Out /  Change )

Facebook photo

You are commenting using your Facebook account. Log Out /  Change )

Connecting to %s

%d bloggers like this: