Problem 135

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

まず、ナイーブな方法。約数を求めて、利用。

import Data.List
divisors n = (1,n):d n 2
where d x m | m^2 > x = []
| mod x m ==  = (m,div x m):d x (m+1)
| otherwise = d x (m+1)
sol = filter((> ).fst).concatMap f.divisors
where f (d1,d2) | mod (d1+d2) 4 /=  = []
| otherwise = let t = div (d1+d2) 4
in nub[(d1-t,t),(d2-t,t)]
main = print.length.filter((==10).length.sol)$ [1..10^6-1]

やはり、予想通り遅い。だって、約数もとめるのに割り算たくさんしてるもん。

ということで、方針転換。探索範囲が決まっているので、

import Data.Array.IO
p135 (lim+1) = do count <- newArray (1,lim)  :: IO (IOUArray Int Int)
let step s = takeWhile(<=lim) [(s+t)*(3*t-s)|t<-[div s 3 +1..]]
increment n = readArray count n>>=writeArray count n.succ
mapM_ increment.concatMap step $ [1..lim]
getElems count
main = p135 (10^6)>>=print.length.filter(==10)

速くなりました。

More Reading
Older// Problem 136