Project Euler Problem 088

import Control.Monad
import Control.Monad.ST
import Control.Applicative
import Control.Arrow
import Data.Ord
import Data.Char
import qualified Data.Map as M
import qualified Data.Set as S
import Data.List
import Data.Array
import Data.Array.ST
import Data.Maybe
import Debug.Trace
import System.IO.Unsafe

------------------ 003より
primes :: [Integer]
primes = 2 : ([3,5..] `minus` join [[p*p, p*p+2*p..] | p <- primes']) 
  where
    primes' = 3 : ([5,7..] `minus` join [[p*p, p*p+2*p..] | p <- primes'])   
    join  ((x:xs):t)    = x : union xs (join (pairs t))
    pairs ((x:xs):ys:t) = (x : union xs ys) : pairs t                      
    union (x:xs) (y:ys) = case (compare x y) of 
      LT -> x : union  xs  (y:ys)
      EQ -> x : union  xs     ys 
      GT -> y : union (x:xs)  ys
    union  xs     []    = xs
    union  []     ys    = ys
    minus (x:xs) (y:ys) = case (compare x y) of 
      LT -> x : minus  xs  (y:ys)
      EQ ->     minus  xs     ys
      GT ->     minus (x:xs)  ys
    minus  xs     _     = xs
--必要に応じてgroupする
factor :: Integer -> [Integer]
factor n = factorimpl n primes where
  factorimpl n pri@(p:xs) =
    if p*p>n then [n]
    else if n`rem` p == 0 then
      p:factorimpl (n `quot` p) pri
         else
           factorimpl n xs

-- haskell.orgより
facts :: [Integer] -> [[Integer]]
facts = concat . takeWhile valid . iterate facts' . (:[])
  where valid xs = length (head xs) > 1
        facts' = nub' . concatMap factsnext
        nub' = S.toList . S.fromList
        factsnext xs = 
          let factsnext' [] = []
              factsnext' (y:ys) = map (form y) ys ++ factsnext' ys
              form a b = a*b : (delete b . delete a $ xs)
          in map sort . factsnext' $ xs      

-- 第二種スターリング数
st2 :: Eq a => [a] -> Int -> [[[a]]]
st2 l 1 = [[l]]
st2 l n | n==length l = [map ((:[])) l]
st2 (x:xs) n =
  concatMap f next
  ++
  (map ([x]:) (st2 xs (n-1))) -- 一つの群[x]にして既存の組に追加
  where
    next = st2 xs n
    f l = map g l where -- この2行をかくのにすごい苦労した。初回以上に
      g ll = ((x:ll):(delete ll l))

-- ベル数
bells :: Eq a => [a] -> [[[a]]]
bells l = concatMap (st2 l) [1..(length l)]
bellgen :: Integer -> [[Integer]]
bellgen n = filter ((>1).length) $ nub $ map (sort.(map product)) $ bells $ factor n

klist :: Integer -> [Integer]
klist n = nub $ map (\x ->n - (sum x) + (fromIntegral $ length x)) filtered where
    filtered = filter ((<=n).sum) $ bellgen n


an :: Num b => Integer -> [b]
an to = map snd $ filter ((<=to).fst) $ animp to
animp :: Num a => Integer -> [(Integer, a)]
animp to = M.toList $ last $ unfoldr f (1,2,M.empty,map klist [2..]) where  
  f (n,v,acc,[]) = Nothing
  f (n,v,acc,(x:xs))| n >= to = Just (acc,(n,v,acc,[]))
  f (n,v,acc,([]:ys)) = Just (acc,(n,v+1,acc,ys))
  f (n,v,acc,((x:xs):ys)) =
    if (not $ M.member x acc) && (x<=to) then
        f (n+1,v,(M.insert x v acc),xs:ys)
    else
      f (n,v,acc,(xs:ys))
ans1 :: Integer
ans1 = sum $ S.toList $ S.fromList $ an 12000 -- nubしている

-- 7587457
main :: IO ()
main = print ans1