Problem 84

http://projecteuler.net/index.php?section=problems&id=84
モノポリーのどのマスに止まりやすいか、という問題。
マルコフ連鎖だとおもう。ただ、状態数が120ぐらいある。
遷移行列を求めるのが面倒だ、というか面倒そう。

import Data.List
import Data.Array.IArray
import Control.Arrow
import Numeric.LinearAlgebra

-- dice = 6
dice = 4
len = 40
jail = 10
go = 0
g2j = 30
cc = [2,17,33]
ch@[ch1,ch2,ch3]=[7,22,36]
c1 = 11
e3 = 24
h2 = 39
[r1,r2,r3] = [5,15,25]
[u1,u2] = [12,28]

xi (i,j) = i*3+j
ix =flip divMod 3

next (i,j) (d1,d2) 
         | j'==3   = (jail,0)  -- 3 consective double
         | i'==g2j = (jail,j')  -- go to jail
         | i'/=g2j = (i',j')   
         where i' = mod (i+d1+d2) len
               j' = if d1==d2 then j+1 else 0

card ((i,(j,n)),p)
    | elem j cc = zipWith pack [go,jail,j][p/16,p/16,7*p/8]
    | elem j ch = zipWith pack 
                  [go,jail,c1,e3,h2,r1,nextU j,mod (j-3) len,nextR j,j]
                  (replicate 8 (p/16)++[p/8,3*p/8])
    | otherwise = [((i,(j,n)),p)]
    where pack j' p' = ((i,(j',n)),p')
          nextR j | j==ch1 = r2
                  | j==ch2 = r3
                  | j==ch3 = r1
          nextU j | j==ch1 = u1
                  | j==ch2 = u2
                  | j==ch3 = u1

move ix = concatMap (card.f).group.sort.map (const ix&&&next ix) $range ((1,1),(dice,dice))
    where f = head&&&(/(fromIntegral$dice^2)).genericLength

transM =  fromArray2D.accumArray (+) 0 ((0,0),(len*3-1,len*3-1))$moves::Matrix Double
    where moves = concatMap (map toXi.move) .map ix $[0..len*3-1]
          toXi = first(xi***xi)

toProb ::Vector Double -> [Double]
toProb = f.toList
    where f [] = []
          f (a:b:c:ds) = (sum.sort$[a,b,c]):f ds
steady = toProb.fst$until enough next (fromList.(1:).replicate (len*3-1)$0,constant 0 $len*3)
    where next (sM,s) = (sM<>transM,sM)
          enough (sM,s) = pnorm Infinity (sM-s) < 0.00000001
main = print.take 3.reverse.sort.zip steady $[0..]

やはり予想どおり、めんどくさかった。
条件を勘違いしていたり、間違えていたりと大変だった。
行列計算にhmatrixを利用。
ちなみにdice=6で実行すると

[(6.233011585917959e-2,10),(3.183021123222021e-2,24),(3.0887948067887767e-2,0)]

だいたいあっているでしょう。