Haskell, n=30 (57s)
With a lot of invaluable contributions by @Angs: use Vector
, use short circuit products, look at odd n.
import Control.Parallel.Strategies
import qualified Data.Vector.Unboxed as V
import Data.Maybe
import Data.Int
type Row = V.Vector Int8
x :: Row -> [Row] -> Integer -> Int -> Integer
x p (v:vs) m c = let c' = c - 1
r = if c>0 then parTuple2 rseq rseq else r0
(a,b) = ( x p vs m c' ,
x (V.zipWith(-) p v) vs (-m) c' )
`using` r
in a+b
x p _ m _ = m * prod p
prod :: Row -> Integer
prod p = fromMaybe 0 $ V.foldM' (\a b -> if b==0 then Nothing else
Just $ a * fromIntegral b) 1 p
p, pt :: [Row] -> Integer
p (v:vs) = x (foldl (V.zipWith (+)) v vs) (map (V.map (2*)) vs) 1 11
`div` 2^(length vs)
p [] = 1 -- handle 0x0 matrices too :-)
pt (v:vs) | even (length vs) = p ((V.map (2*) v) : vs ) `div` 2
pt mat = p mat
main = getContents >>= print . pt . map V.fromList . read
My first attempts at parallelism in Haskell. You can see a lot of optimization steps through the revision history. Amazingly, it were mostly very small changes. The code is based on the formula in the section "Balasubramanian-Bax/Franklin-Glynn formula" in the Wikipedia article on computing the permanent.
p
computes the permanent. It is called via pt
which transforms the matrix in a way that is always valid, but especially useful for the matrices that we get here.
Compile with ghc -O2 -threaded -o <name> <name>.hs
. To run with parallelisation, give it runtime parameters like this: ./<name> +RTS -N
. Input is from stdin with nested comma separated lists in brackets like [[1,2],[3,4]]
as in the last example (newlines allowed everywhere).