primes = (2 : filter (check primes) [3..]) where check (p:ps) n | p*p > n = True | n `mod` p == 0 = False | otherwise = check ps n primeDivs n = divs n primes where divs 1 _ = [] divs n (p:ps) = if p*p > n then [(n, 1)] else let (c, r) = cnt n p 0 in if c == 0 then divs n ps else (p, c) : divs r ps cnt n p c | mod n p == 0 = cnt (div n p) p (c+1) | otherwise = (c, n) combine divs = comb [1] divs where comb list [] = list comb list ((p, c) : ps) = comb [ x*y | x <- take (c+1) $ iterate (*p) 1, y <- list ] ps sumDiv n = (sum $ combine $ primeDivs n) - n friends :: [(Integer, Integer)] friends = filter (\(x, s) -> x < s && sumDiv s == x) $ map (\x -> (x, sumDiv x)) [1..] main = mapM print $ take 30 friends