n倍完全数

• ref:ダブル完全数

自分自身も含めた約数の和が自分自身の n 倍になる数を n 倍完全数というらしい。ということで、n 倍完全数を列挙するプログラムとか。
愚直な解き方。

def divisors(n):
    return (i for i in xrange(1, n + 1) if n % i == 0)

def perfect_numbers(n, limit):
    def pred(x):
        return x * n == sum(divisors(x))

    return (i for i in xrange(1, 1 + limit) if pred(i))

if __name__ == '__main__':
    for i in perfect_numbers(3, 10000):
        print i

i が n の約数なら n / i も n の約数。これなら、自身の平方根までで探索を打ち切れる。

from itertools import count, takewhile

def divisors(n):
    for i in takewhile(lambda i: i * i <= n, count(1)):
        if n % i == 0:
            yield i
            if i * i != n:
                yield n / i

def perfect_numbers(n, limit):
    def pred(x):
        return x * n == sum(divisors(x))

    return (i for i in xrange(1, 1 + limit) if pred(i))

if __name__ == '__main__':
    for i in perfect_numbers(3, 10000):
        print i

i * j == n で i が素数なら n の約数は j の約数 or j の約数 * i になるはず。で、 2 から探して最初に見付かった約数は素数。

from itertools import takewhile, count
import operator

def divisors(n):
    for i in takewhile(lambda i: i * i <= n, count(2)):
        if n % i == 0:
            div = divisors(n / i)
            return set(j * i for j in div) | div
    else:
        return set((1, n))

def perfect_numbers(n, limit):
    def pred(x):
        return x * n == sum(divisors(x))

    return (i for i in xrange(2, 1 + limit) if pred(i))
    

if __name__ == '__main__':
    for n in perfect_numbers(3, 10000):
        print n

約数を見付けるより倍数を列挙する方が簡単なので、表を作って、自身の倍数を表の要素に加えていく。

from itertools import takewhile, count

def perfect_numbers(n, limit):
    table = [1] * (limit + 1)

    for i in xrange(2, limit / 2 + 1):
        for j in takewhile(lambda j: i * j <= limit, count(1)):
            table[i * j] += i

    return [i for i, s in enumerate(table) if i * n == s]

if __name__ == '__main__':
    for n in perfect_numbers(3, 10000):
        print n

速度比較。番号はここに列挙している順。

% for i in 1 2 3 4; do time python perfect_numbers$i.py > /dev/null; done
python perfect_numbers$i.py > /dev/null  12.21s user 0.02s system 99% cpu 12.328 total
python perfect_numbers$i.py > /dev/null  0.51s user 0.02s system 97% cpu 0.542 total
python perfect_numbers$i.py > /dev/null  0.41s user 0.02s system 94% cpu 0.456 total
python perfect_numbers$i.py > /dev/null  0.10s user 0.00s system 80% cpu 0.119 total

まぁ、こんなものか。ちなみに、3 番目のは memoize すれば速くなるかと思って、下のような memoize decorator を作ってみたけど全然効果がなかった。

def memoize(func):
    memo = dict()

    def memoized_func(*args):
        if args not in memo:
            memo[args] = func(*args)
        return memo[args]
    return memoized_func