DEV Community

Peter Kim Frank
Peter Kim Frank

Posted on

Write a script to find "Perfect Numbers"

Inspired by this tweet:

Write a script to find "Perfect Numbers." Here's the list of numbers up to 19 digits ๐Ÿ˜ˆ:

6, 28, 496, 8138, 33550336, 8589869056, 137438691328, 2305843008139952128

Look forward to seeing everyone's solutions!

Top comments (40)

Collapse
 
nektro profile image
Meghan (she/her) • Edited
const perfect = [6, 28, 496, 8128, 33550336, 8589869056, 137438691328, 
2305843008139952128];

perfect.map(v => v.toString(2));
/*
[ "110",
  "11100",
  "111110000",
  "1111111000000",
  "1111111111111000000000000",
  "111111111111111110000000000000000",
  "1111111111111111111000000000000000000",
  "1111111111111111111111111111111000000000000000000000000000000" ] */

perfect.map(v => v.toString(3));
/*
[ "20",
  "1001",
  "200101",
  "102011001",
  "2100010112102001",
  "211011122100120010101",
  "111010202100010220210001",
  "120100210022221112202020222210100212001" ] */

perfect.map(v => v.toString(4));
/*
[ "12",
  "130",
  "13300",
  "1333000",
  "1333333000000",
  "13333333300000000",
  "1333333333000000000",
  "1333333333333333000000000000000" ] */

I find base 2 and 4 particularly interesting because for each perfect number, the amount of 1s/3s and 0s is the same.

edit: fixed my sequence because 8138 is not perfect, its 8128.

Collapse
 
shiling profile image
Shi Ling

Yeah, I had a suspicion that there's a geometric solution to finding perfect numbers and was trying to formulate it in my head. Putting it in base2 makes it quite clear. Hm... ๐Ÿค”

Collapse
 
kspeakman profile image
Kasey Speakman • Edited

Another interesting observation. In the binary form of each perfect number, the count of leading 1s is a prime number. In fact, they are all Mersenne Exponents.

Collapse
 
itachiuchiha profile image
Itachi Uchiha

I wrote this script 3 years ago.

def perfect_number(num):
    su = 0
    for i in range(1, num):
        if num % i == 0:
            su += i
    return su == num

def perfect_number_main(n):
    if perfect_number(n) == True:
        print("%s is perfect" % n)
    else:
        pass

for i in range(1,10000):
    perfect_number_main(i)
Collapse
 
chenge profile image
chenge • Edited

Translate yours in Ruby.

def sum_divisions(n)
  (1...n).select{|i| n % i == 0}.sum
end

def perfect_number?(n)
  n == sum_divisions(n)
end

(1..10_000).each do |i|
  puts "#{i} is perfect" if perfect_number?(i)
end

Collapse
 
chenge profile image
chenge • Edited

In JS, not so good cause I'm new to JS. If any good writing comment please.

const sum = (accumulator, currentValue) => accumulator + currentValue;

function range(size, startAt = 0) {
  return [...Array(size).keys()].map(i => i + startAt);
}

function sumDivisions(n){
  arr = range(n-1, startAt=1).filter(x => n % x == 0)
  return arr.reduce(sum, 0);
}

function is_perfect_number(n){
  return sumDivisions(n) == n
}

range(10000, 1).forEach(function(x){
  if(is_perfect_number(x)){
    console.log(x);
  }
})
Collapse
 
chenge profile image
chenge

Do you have a better version for Go?

package main

import (
  "fmt"
)

func sum_divisions(n int) int {
  sum := 0
  for i := 1; i < n; i++ {
    if n % i == 0 {
      sum += i
    }
  }
  return sum
}

func is_perfect(n int) bool {
  return n == sum_divisions(n)
}

func main() {
  for i := 1; i <= 10000; i++ {
    if is_perfect(i) {
      fmt.Println(i)
    }
  }
}

Collapse
 
kaelscion profile image
kaelscion

How efficient do you remember this article being?

Collapse
 
kspeakman profile image
Kasey Speakman • Edited

F#

let getDivisors n =
    let limit = int64 (sqrt (float n))
    [
        for i = 1L to limit do
            if n % i = 0L then
                yield i
                let result = n / i
                if result <> i && result <> n then
                    yield result
    ]


let isPerfect n =
    n = (n |> getDivisors |> List.sum)

Usage

(*
   per @nektro's observation on perfect numbers in binary
   * perfect numbers have an odd number of bits
   * the first half (rounded up) are 1s
   * the second half (rounded down) are 0s
*)

// Creates a perfect number candidate of the given bit length.
// Not all candidates are perfect numbers.
// Example: createCandidate 3 returns 110b or 6
//          createCandidate 5 returns 11100b or 28
let createCandidate bitLength =
    let zeroStart = bitLength / 2 + 1
    let update num i =
        let shifted = num * 2L
        if i < zeroStart then shifted + 1L
        else                  shifted
    List.init bitLength id |> List.fold update 0L

let oddBitLengths = // 1, 3, 5, 7, 9, ... , 61
    Seq.init 31 (fun i -> i * 2 + 1)
    // should be 32, but last one overflows on sum step

oddBitLengths
|> Seq.map createCandidate
|> Seq.filter isPerfect
|> Seq.iter (fun i ->
    printfn "%20i %s" i (Convert.ToString(i, 2))
)
//                   1 1
//                   6 110
//                  28 11100
//                 496 111110000
//                8128 1111111000000
//            33550336 1111111111111000000000000
//          8589869056 111111111111111110000000000000000
//        137438691328 1111111111111111111000000000000000000
// 2305843008139952128 1111111111111111111111111111111000000000000000000000000000000

To filter the numbers down to just likely candidates, I am using @nektro 's observation. I first generate a limited list of possible candidates that fit in a 64 bit integer, then filter down to only perfect ones. It takes 40s or so to run, but all of that time is spent between the last and next-to-last result.

Collapse
 
gmartigny profile image
Guillaume Martigny

How much time did this script took to run ? You approach is really nice, but I'm worried by |> Seq.filter isPerfect perf.

Collapse
 
kspeakman profile image
Kasey Speakman • Edited

Runtime is about 40 seconds in Release mode. The largest cost is finding divisors. I used the least expensive algorithm I could find. But it is still expensive. So to get the runtime down that low, I optimize down the number of inputs using nektroโ€™s observation.

Thread Thread
 
gmartigny profile image
Guillaume Martigny • Edited

You could gain a little by only generating prime numbers for your candidates. Eratosthenes' sieve is not that hard to implement.

Thread Thread
 
kspeakman profile image
Kasey Speakman

I will post another solution using primes. Thanks for the tip!

Thread Thread
 
kspeakman profile image
Kasey Speakman • Edited

I posted the prime-based solution here. The actual perfect number calculation time is now negligible, and the code can be a one-liner. But it requires being fed Mersenne Exponents. Fortunately, very few/small Mersenne Exponents quickly turn into perfect numbers which are too large for traditional memory structures.

Collapse
 
gmartigny profile image
Guillaume Martigny • Edited

First dumb "look at every integer" solution:

const range = n => (new Array(+n)).fill().map((_, i) => i + 1);
return range(limit).filter((n) => {
    return range(Math.ceil(n ** 0.5))
        .filter(d => n % d === 0)
        .reduce((acc, val) => acc + val + (n / val), 0) / 2 === n;
});

This algo has trouble going further than 1e6.

Then, I dig your hypothesis that all Perfect number are in the form (n+1)"1"..(n)"0" in base 2. I needed to find the logical sequence of valid n.

[6, 28, 496, 8128, 33550336, 8589869056, 137438691328].map(n => n.toString(2).match(/0/g).length);
// => [1, 2, 4, 6, 12, 16, 18]

So I looked up this sequence on EOIS and luckily it found a match: "Numbers n such that 2n+1 - 1 is prime".

Which lead to this way more powerful script:

const range = n => (new Array(+n)).fill().map((_, i) => i + 1);
const isPrime = (n) => {
    if (n <= 1) return false;
    if (n <= 3) return true;
    if (n % 2 === 0 || n % 3 === 0) return false;
    for (let i = 5, l = n ** 0.5; i <= l; i += 6) {
        if (n % i === 0 || n % (i + 2) === 0) return false;
    }
    return true;
};
return range(limit).filter(n => isPrime(2 ** (n + 1) - 1)).map((n) => {
    return parseInt("1".repeat(n + 1) + "0".repeat(n), 2);
});

This time, it can't go any higher than 1e3 before crashing, but yield more results. (However, JS precision let me down on the last result)

Collapse
 
shiling profile image
Shi Ling • Edited

Came to a same thesis and worked out a similar solution. Also ran into the limits at the 9th perfect number. But JavaScript now supports native BigInt (only in Chrome and I believe newer versions of NodeJs though!). Now we just need a way to generate large prime numbers fast to generate perfect numbers. ๐Ÿ˜‰

Collapse
 
gmartigny profile image
Guillaume Martigny

Thanks to your input on BigInt, I made a new solution.
Thanks a lot !

Thread Thread
 
shiling profile image
Shi Ling

๐Ÿ‘๐Ÿ‘๐Ÿ‘๐Ÿ‘๐Ÿ‘๐Ÿ‘

I feel a heart is not enough. I'm so happy to see the new version.

Collapse
 
gmartigny profile image
Guillaume Martigny

Yeah, I liked that we came with two solution from the same hypothesis.

I'm looking into BigInt (thanks a lot I didn't know !) and it looks promising.

Collapse
 
modaf profile image
Modaf • Edited

Inspired from programming-algorithms.net/article...
Takes less than 1sec to get to 137438691328 but can't find easily 2305843008139952128

from math import sqrt
def perfect(n) :
    s = 1
    rac = int(sqrt(n))
    for k in range(2, rac+2) :
        if n%k == 0 :
            s += k + n//k
    return s == n
def f(x) :
    return (2**(x-1))*((2**x)-1)
n = 137438691329
m = 1
res = []
while True :
    val = f(m)
    print(val)
    if perfect(val) :
        res.append(val)
    m+=1
    if val > n :
        break
print(res)
Collapse
 
chenge profile image
chenge

site is great for algo, thanks.

Collapse
 
rapidnerd profile image
George • Edited

Here's my solution in C

#include <stdio.h>

void  main(){
  int n, i, sum;
  int mn, mx;

  printf("Starting range: ");
  scanf("%d",&mn);

  printf("Ending range: ");
  scanf("%d",&mx);

  printf("Output: ");
  for(n = mn; n <= mx ; n++){
    i = 1;
    sum = 0;
    while(i < n){
      if(n % i == 0)
           sum = sum + i;
          i++;
    }
    if(sum== n)
      printf("%d ",n);
  }
      printf("\n");
}
Collapse
 
anduser96 profile image
Andrei Gatej

I think you can reduce the number of iterations by using โ€œ while (i <=n/2)โ€.

Collapse
 
stefanrendevski profile image
Stefan Rendevski

You can reduce it even further by using "while (i <= sqrt(n))" and noticing that divisors come in pairs, for example: if 16 / 2 = 8, both 2 and 8 are divisors of 16. There is no need to go up to 8 and check all integers up until that point, you only need to perform the division.

Collapse
 
Sloan, the sloth mascot
Comment deleted
Collapse
 
florimondmanca profile image
Florimond Manca • Edited

A (lazy but slow) Python solution which hasn't even yielded 33550336 in 5 minutes, despite a subtle optimization allowing to find divisors in O(โˆšn)โ€ฆ Probably could use some optimization to prevent trying all possible integers. :-)

def find_divisors(n: int):
    divisors = []
    m = 1
    while m * m < n:
        if n % m == 0:
            yield m
            if m != 1:
                yield n / m
        m += 1


def is_perfect(n: int) -> bool:
    return sum(find_divisors(n)) == n


if __name__ == '__main__':
    from itertools import cycle
    all_integers = cycle()
    for n in filter(is_perfect, all_integers):
        print(n)

Collapse
 
phallstrom profile image
Philip Hallstrom • Edited

Working on getting better at Clojure so thought I'd give it a shot.

(ns philip.perfect-numbers)

(defn brute-force-perfect-number? [num]
  (== num (apply + (filter #(zero? (mod num %)) (range 1 (inc (/ num 2)))))))

(time (println (take 4 (filter brute-force-perfect-number? (iterate inc' 2)))))
; (6 28 496 8128)
; "Elapsed time: 2649.69715 msecs"

(defn perfect-number? [num]
  (loop [divisor 2
         subtotal 1]
    (let [
          q (quot num divisor)
          evenly-divisible? (zero? (mod num divisor))]
    (if (>= divisor q)
      (== num subtotal)
      (recur (inc divisor) (if evenly-divisible? (+ subtotal divisor q) subtotal))))))

(time (println (take 4 (filter perfect-number? (iterate inc' 2)))))
; (6 28 496 8128)
; "Elapsed time: 66.942315 msecs"
Collapse
 
kspeakman profile image
Kasey Speakman • Edited

F#, ridiculously short/fast solution

// perfect num calculation is 2^(p-1)*(2^p - 1)
// source: https://oeis.org/A000396
let getPerfect p =
    pown 2L (p - 1) * (pown 2L p - 1L)

Usage

// prime calculation is out of scope for this problem
// (and this few can be looked up easily)
// so I simply list them here.
let mersenneExponents =
    [1; 2; 3; 5; 7; 13; 17; 19; 31]

List.map getPerfect mersenneExponents
|> printfn "%A"
// [1L; 6L; 28L; 496L; 8128L; 33550336L; 8589869056L; 137438691328L; 2305843008139952128L]

This runs in 0.35ms or ~3500 ticks in Release mode (after a JIT prerun). It basically turns the problem around to be just a matter of inputting Mersenne Exponents. Not many are needed before overflowing traditional memory structures like int64.

Collapse
 
albxncom profile image
albxn • Edited

Well, here goes my Haskell version. Enjoy!

perfectNumber :: Int -> Bool
perfectNumber i = 
    sum (divisorsList i) == i

perfectNumbers :: Int -> [Int]
perfectNumbers n =
    take n [ x | x <- [1..], perfectNumber x ]

divisorsList :: Int -> [Int]
divisorsList i = 
    [ x | x <- [1..(i-1)], i `mod` x == 0]

main = do
    print(perfectNumbers 5) -- prints the first 5 perfect numbers..
Collapse
 
antonrich profile image
Anton

A very minor tweak. i-1 doesn't require parentheses in this case.

divisorsList i = 
    [ x | x <- [1..i-1], i `mod` x == 0]
Collapse
 
gmartigny profile image
Guillaume Martigny

I kind of got crazy with this challenge. After a day at refining my solution, I end up with a pretty efficient solution. The script is able to find the 18th perfect number under 14s on my computer.

Since the solution is a bit long and require some explanations, I wrote everything on a gist.
If anyone found a way to reduce computations, I would greatly appreciate.

Collapse
 
peter profile image
Peter Kim Frank

This is incredible. You should consider converting your gist into a DEV article!

Collapse
 
gmartigny profile image
Guillaume Martigny

Yes, I'm so proud of the result. But in the end, this is nothing more than small bits gathered from the Internet. I don't see an article out of this.

If you (or any one else reading this) want to use this resource for an article, feel free to do so. ;)

Collapse
 
shiling profile image
Shi Ling • Edited

My solution in javascript.

Time taken to get the 7th number which is 8 digits long is 8 milliseconds.
Time taken to get the 8th number which is 19 digits long is 2.6 seconds.
Still waiting for my results for the 9th number...

function isPerfect(num){
    let sum = 1 
    let max = num / 2
    for(var i = 2; i < max; i++){
        if(num % i === 0){
            sum += i // add lower divisor
            max = Math.floor(num / i)
            sum += max //add higher divisor
        }
    }
    return sum === num;
}

function getPerfectSeries(n = 5){ 

    let series = []

    // Create number in binary string
    let l = 1
    while(series.length < n){
        let s = "1"
        for(let i=0; i < l; i++){
            s += "1"
        }
        for(let i=0; i < l; i++){
            s += "0"
        }
        let num = parseInt(s,2) // convert to base10
        if(isPerfect(num)){
            series.push(num)
        }
        l++
    }

    console.log(series)

}

function timer(fn){
    let start = new Date();
    fn()
    let time = new Date().getTime() - start.getTime()
    console.log("Time taken: " + time + "ms")
}

timer(()=>{getPerfectSeries(1)}) // < 1ms
timer(()=>{getPerfectSeries(2)}) // < 1ms
timer(()=>{getPerfectSeries(3)}) // < 1ms
timer(()=>{getPerfectSeries(4)}) // < 1ms
timer(()=>{getPerfectSeries(5)}) // < 1ms
timer(()=>{getPerfectSeries(6)}) // 4ms
timer(()=>{getPerfectSeries(7)}) // 8ms
timer(()=>{getPerfectSeries(8)}) // 26.9 seconds

Anyway, this solution is inspired by @nektro 's observation of the perfect numbers - Had a suspicion that there's a geometric pattern to this which is obvious when you look like the perfect numbers is Base2.

const perfect = [6, 28, 496, 8128, 33550336, 8589869056, 137438691328, 
2305843008139952128];

perfect.map(v => v.toString(2));
/*
[ "110",
  "11100",
  "111110000",
  "1111111000000",
  "1111111111111000000000000",
  "111111111111111110000000000000000",
  "1111111111111111111000000000000000000",
  "1111111111111111111111111111111000000000000000000000000000000" ] */

perfect.map(v => v.toString(3));
/*
[ "20",
  "1001",
  "200101",
  "102011001",
  "2100010112102001",
  "211011122100120010101",
  "111010202100010220210001",
  "120100210022221112202020222210100212001" ] */

perfect.map(v => v.toString(4));
/*
[ "12",
  "130",
  "13300",
  "1333000",
  "1333333000000",
  "13333333300000000",
  "1333333333000000000",
  "1333333333333333000000000000000" ] */

I find base 2 and 4 particularly interesting because for each perfect number, the amount of 1s/3s and 0s is the same.

edit: fixed my sequence because 8138 is not perfect, its 8128.

Perfect numbers in binary always start with "1" followed by an equal number of "1" and "0"s.

I first generate plausible perfect numbers from binary strings first using this pattern, then converted them to Base10, before testing if it is a perfect number.

The isPerfect function can be optimised by reducing the loop to find the divisors by 1) adding both the lower and the higher divisor and 2) also setting the higher divisor as a limiter for iteration.

Collapse
 
shiling profile image
Shi Ling • Edited

Results

Still waiting.

@picocreator - Wanna use GPU.js to speed this up?

Collapse
 
gmartigny profile image
Guillaume Martigny

First dumb "look at every integer" solution:

const range = n => (new Array(+n)).fill().map((_, i) => i + 1);
return range(limit).filter((n) => {
    return range(Math.ceil(n ** 0.5))
        .filter(d => n % d === 0)
        .reduce((acc, val) => acc + val + (n / val), 0) / 2 === n;
});

This algo has trouble going further than 1e6.

Then, I dig your hypothesis that all Perfect number are in the form (n+1)"1"..(n)"0" in base 2. I needed to find the logical sequence of valid n.

[6, 28, 496, 8128, 33550336, 8589869056, 137438691328].map(n => n.toString(2).match(/0/g).length);
// => [1, 2, 4, 6, 12, 16, 18]

So I looked up this sequence on EOIS and luckily it found a match: "Numbers n such that 2n+1 - 1 is prime".

Which lead to this way more powerful script:

const range = n => (new Array(+n)).fill().map((_, i) => i + 1);
const isPrime = (n) => {
    if (n <= 1) return false;
    if (n <= 3) return true;
    if (n % 2 === 0 || n % 3 === 0) return false;
    for (let i = 5, l = n ** 0.5; i <= l; i += 6) {
        if (n % i === 0 || n % (i + 2) === 0) return false;
    }
    return true;
};
return range(limit).filter(n => isPrime(2 ** (n + 1) - 1)).map((n) => {
    return parseInt("1".repeat(n + 1) + "0".repeat(n), 2);
});

This time, it can't go any higher than 1e3 before crashing, but yield more results. (However, JS precision let me down on the last result)

Collapse
 
maxart2501 profile image
Massimo Artizzu

You're... not expecting fast implementations, are you? ๐Ÿคจ

Unless you consider scraping a web page to get them. There are just 50 known perfect numbers...

Collapse
 
peter profile image
Peter Kim Frank • Edited

Going beyond the set of currently known perfect numbers is extra credit. But no - any old method will do, fast or slow or limited to just the first few in the list. :)

It will be interesting to see what approach is most efficient for grabbing the first few.

 
chenge profile image
chenge

Local can be short, may be r is better.

Clear and ASAP.