Probrem 21

今日は Probrem 21 をやってみた。

 cf. http://projecteuler.net/index.php?section=problems&id=21

数 n に対して約数(nを含まない)すべてを足しあわせる関数 d があるとして,d(a) = b かつ d(b) = a となるような a と b (ただしa≠b)を amicable pair,それぞれを amicable number という(んだそうだ)。で,10000未満の amicable number をすべて足しあわせろ,という問題。
内包表記大活躍……なのはいいけど,すっげー時間がかかる。10000未満なんてやったら待ってられないので,下のコードでは1000未満にしてある。いい方法はないだろうか。

module Main (main) where

d :: Int -> Int
d n = sum [x | x <- [1..(n `div` 2)], n `mod` x == 0]

euler021 :: Int -> Int
euler021 n = sum $ map (\(a,b) -> a + b) pairs
  where
    pairs = [(a,b) | a <- [1..n], b <- [1..n], a < b, d a == b, d b == a]

main :: IO ()
main = print $ euler021 1000
^o^ >euler021.exe
504

1000未満では 220 と 284 のペアしかないみたいだ。

追記:

いい方法はないだろうか,と書いたらコメントやらトラックバックやらでいろいろ教えてくれた。ありがとうございます。

 cf. http://haskell.g.hatena.ne.jp/nobsun/20080314/p1

おおっ,速い!

追記2:

ううっ,うっかりとエントリに上書きして消してしまった。できるだけ復元したけど元のエントリとは少し違ってるかも。

「Probrem 21」への5件のフィードバック

  1. 取り敢えず
    d n = (1+) $ sum $ do
     x <- takeWhile (¥x -> x*x <= n) [2..]
     let (y,r) = n `divMod` x
     guard $ r == 0
     if (y == x) then [x] else [x,y]
    で少し速くなります。

    後は
    ・素因数分解してオイラー関数を使う
     http://ja.wikipedia.org/wiki/%E3%82%AA%E3%82%A4%E3%83%A9%E3%83%BC%E3%81%AE%CF%86%E9%96%A2%E6%95%B0
    ・メモ化
    かな。

  2. ああ。pairs が同じ計算を何度もやるからおそいですよ。
    eular21 n = sum $ map (uncurry (+)) pairs
    where pairs = [(a,b) | a <- [1..n], let b = d a, a /= b && a == d b ]
    てな感じでいかが?

コメントを残す

メールアドレスが公開されることはありません。 が付いている欄は必須項目です

このサイトはスパムを低減するために Akismet を使っています。コメントデータの処理方法の詳細はこちらをご覧ください