Probelm51

30分プログラム、その327。Problem51 - Project Euler。
ちょっと方針を変えて、適当な素数を1つ選んで、それを置き換えて素数かどうか確認する方法にしてみた。とりあえず、計算は進んでるみたいだけれども、まだ答えはでてない。

追記: 2時間ほどかけて、答えがでた。まあ、これでいいや。

使い方

$ time ./problem51
.....
121313
[[1,2,1,3,1,3],[2,2,2,3,2,3],[3,2,3,3,3,3],[4,2,4,3,4,3],[5,2,5,3,5,3],[6,2,6,3,6,3],[8,2,8,3,8,3],[9,2,9,3,9,3]]
./problem51  9413.65s user 81.08s system 99% cpu 2:39:30.06 total

ソースコード

import Control.Monad.List
import Data.Char
import Data.List
import Debug.Trace

sieve (x:xs) = x:sieve [y | y <- xs, y `mod` x /= 0]

primes = sieve [2..]

isPrime n = n == (head $ dropWhile (< n) primes)

splitInt :: Int -> [Int]
splitInt n = unfoldr f n
    where f 0 = Nothing
          f n = Just (n `mod` 10,n `div` 10)

uniqC :: Num a => [a] -> [(a,Int)]
uniqC []  = []
uniqC [x] = [(x,1)]
uniqC (x:xs) = let yss@((y,n):ys) = uniqC xs
               in if x == y then
                      (y,n+1):ys
                  else 
                      (x,1):yss

replace :: [Int] -> Int -> Int -> [Bool] -> [Int]
replace [] _ _ [] = []
replace [] _ _ _  = error "invalid mask"
replace x _ _ [] = x
replace (x:xs) old new mask@(y:ys)= case (x == old,y) of
                                      (True,True)  ->
                                          new:(replace xs old new ys)
                                      (True,False) ->
                                          old:(replace xs old new ys)
                                      (_,_) ->
                                          x:(replace xs old new mask)

perm [] = [[]]
perm xs = nub $ do x <- xs
                   ys <- perm $ delete x xs
                   return (x:ys)

primesN n = let min = 10^(n-1)
                max = 10^n
            in takeWhile (< max) $ dropWhile (<min) primes

makeInt xs = foldl1 (\x y->10*x+y) xs 

replacedInt p = do let p' = reverse $ splitInt p
                   (old,count) <- uniqC $ sort $ p'
                   count'    <- [1..count]
                   mask <- perm $ take count $ (replicate count' True) ++ (cycle [False])
                   return $ filter check $ map (\new -> replace p' old new mask)  [0..9]
    where check (0:_) = False
          check xs    = isPrime $ makeInt xs

solveN n = do p <- primes
              sameGroup <- replacedInt (trace (show p) p)
              guard $ length sameGroup == n
              return sameGroup

answer = head $ solveN 8
main = print answer

参考

• 過去の30分プログラム
• Problem 51 - PukiWiki