Problem23が解けない

久しぶりの30分プログラム、その291。Problem23がなかなか解けない。

現在の方針は、

  • 素因数分解を利用して約数の和を計算
  • 過剰数のリストを生成する
  • それで表わせない数のリストを生成する
  • 和を計算する

という非常に素直なもの。どうも"それで表わせない数のリストを生成する"が重そう。いったんじっくり考えてみよう。

使い方

$ python problem23.py
....

ソースコード

#! /usr/bin/python
# -*- mode:python; coding:utf-8 -*-
#
# problem23.py -
#
# Copyright(C) 2008 by mzp
# Author: MIZUNO Hiroki / mzpppp at gmail dot com
# http://howdyworld.org
#
# Timestamp: 2008/04/19 09:26:46
#
# This program is free software; you can redistribute it and/or
# modify it under MIT Lincence.
#
## 素数リストを作る
def sieve(max):
    (n,xs) = (2,[])
    while n < max:
        if all(map(lambda x: n % x != 0,xs)):
            xs.append(n)
        n += 1
    return xs
primes = sieve(100)

## 素因数分解を利用した約数の計算
def factor(n,primes):
    def f(x,y,n):
        if x % y == 0:
            return f(x/y,y,n+1)
        else:
            return n
    xs = []
    for prime in primes:
        if n < prime:
            return xs
        elif n % prime == 0:
            xs.append((prime,f(n,prime,0)))
    return xs

def product(xs):
    return reduce(lambda x,y: x*y,xs,1)

def divisor_sum(n):
    xs = factor(n,primes)
    return product(map(lambda (p,a):(p**(a+1)-1)/(p-1),xs))

## 過剰数の計算
def is_abundant(n):
    return divisor_sum(n) > n

def abudant_list(n):
    return filter(is_abundant,xrange(12,n))

## 組み合せで表現可能かの計算
def is_combination(n,l):
    """is n expressed by some combination of l"""
    for x in l:
        if n - x in l:
            return True
    return False

def main(n):
    xs = abudant_list(n)
    ys = map(lambda x: not is_combination(x,xs),
            xrange(1,n))
    return sum(ys)

if __name__ == '__main__':
    print main(28123)