diff --git a/LICENSE b/LICENSE --- a/LICENSE +++ b/LICENSE @@ -1,4 +1,4 @@ -Copyright (c) 2013, Stephen Diehl +Copyright (c) 2013-2017, Stephen Diehl Permission is hereby granted, free of charge, to any person obtaining a copy of this software and associated documentation files (the "Software"), to diff --git a/picosat.cabal b/picosat.cabal --- a/picosat.cabal +++ b/picosat.cabal @@ -1,12 +1,12 @@ name: picosat -version: 0.1.3 +version: 0.1.4 synopsis: Bindings to the PicoSAT solver homepage: https://github.com/sdiehl/haskell-picosat license: MIT license-file: LICENSE author: Stephen Diehl maintainer: stephen.m.diehl@gmail.com -copyright: 2014 Stephen Diehl +copyright: 2014-2017 Stephen Diehl Category: Logic build-type: Simple cabal-version: >=1.10 @@ -18,6 +18,9 @@ , GHC == 7.8.3 , GHC == 7.8.4 , GHC == 7.10.1 + , GHC == 7.10.2 + , GHC == 7.10.3 + , GHC == 8.0.1 extra-source-files: cbits/picosat.h Bug-Reports: https://github.com/sdiehl/haskell-picosat/issues @@ -32,15 +35,38 @@ exposed-modules: Picosat other-extensions: ForeignFunctionInterface - ghc-options: -Wall -O2 -fwarn-tabs cc-options: -funroll-loops if os(windows) cc-options: -DNGETRUSAGE - build-depends: base >=4.6 && <4.9 + build-depends: + base >=4.6 && <5.0, + transformers >=0.4 && <0.6, + containers >=0.4 && <0.6 default-language: Haskell2010 Hs-source-dirs: src Include-dirs: cbits C-sources: cbits/picosat.c + +Test-Suite Sudoku + type: exitcode-stdio-1.0 + Hs-source-dirs: test + main-is: Sudoku.hs + default-language: Haskell2010 + build-depends: base, picosat + +Test-Suite Scoped + type: exitcode-stdio-1.0 + Hs-source-dirs: test + main-is: Scoped.hs + default-language: Haskell2010 + build-depends: base, picosat, transformers, containers + +Test-Suite rand-shared-improvement + type: exitcode-stdio-1.0 + Hs-source-dirs: test + main-is: Rand.hs + default-language: Haskell2010 + build-depends: base, picosat, transformers, random, rdtsc diff --git a/src/Picosat.hs b/src/Picosat.hs --- a/src/Picosat.hs +++ b/src/Picosat.hs @@ -45,6 +45,43 @@ For a higher level interface see: + +If you intend to solve a set of similar CNFs think about using +Picosat's incremental interface. It allows to push and pop +sets of clauses, as well as solving under assumptions. + +@ +import Picosat (evalScopedPicosat, addBaseClauses, + withScopedClauses, scopedAllSolutions, + scopedSolutionWithAssumptions) + +main :: IO [Int] +main = + evalScopedPicosat $ do + addBaseClauses [[1, 2, 3]] + -- == [Solution [1,2,3], + -- Solution [1,2,-3], + -- Solution [1,-2,3], + -- Solution [1,-2,-3], + -- Solution [-1,-2,3], + -- Solution [-1,2,-3], + -- Solution [-1,2,3]] + + withScopedClauses [[-2,-3]] $ do + sol <- scopedAllSolutions + -- == [Solution [-1,2,-3], + -- Solution [-1,-2,3], + -- Solution [1,-2,-3], + -- Solution [1,-2,3], + -- Solution [1,2,-3]] + + addBaseClauses [[-1,-3]] + + withScopedClauses [[-1,-2], [1,-3]] $ do + sol <- scopedSolutionWithAssumptions [1] +@ + + -} module Picosat ( @@ -52,7 +89,13 @@ solveAll, unsafeSolve, unsafeSolveAll, - Solution(..) + Picosat, + Solution(..), + evalScopedPicosat, + addBaseClauses, + withScopedClauses, + scopedAllSolutions, + scopedSolutionWithAssumptions ) where import Control.Monad @@ -61,24 +104,53 @@ import Foreign.Ptr import Foreign.C.Types +import Control.Monad.Trans.State.Strict +import Control.Monad.IO.Class + +import qualified Data.Set as S + foreign import ccall unsafe "picosat_init" picosat_init - :: IO (Ptr a) + :: IO (Picosat) foreign import ccall unsafe "picosat_reset" picosat_reset - :: Ptr a -> IO () + :: Picosat -> IO () foreign import ccall unsafe "picosat_add" picosat_add - :: Ptr a -> CInt -> IO CInt + :: Picosat -> CInt -> IO CInt foreign import ccall unsafe "picosat_variables" picosat_variables - :: Ptr a -> IO CInt + :: Picosat -> IO CInt foreign import ccall unsafe "picosat_sat" picosat_sat - :: Ptr a -> CInt -> IO CInt + :: Picosat -> CInt -> IO CInt foreign import ccall unsafe "picosat_deref" picosat_deref - :: Ptr a -> CInt -> IO CInt + :: Picosat -> CInt -> IO CInt +foreign import ccall unsafe "picosat_push" picosat_push + :: Picosat -> IO CInt + +foreign import ccall unsafe "picosat_pop" picosat_pop + :: Picosat -> IO CInt + +-- foreign import ccall unsafe "picosat_context" picosat_context +-- :: Picosat -> IO CInt + +foreign import ccall unsafe "picosat_assume" picosat_assume + :: Picosat -> CInt -> IO () + + +type Picosat = Ptr () + +-- | Call a monadic action with a freshly created Picosat that +-- is destroyed afterwards. +withPicosat :: (Picosat -> IO a) -> IO a +withPicosat f = do + pico <- picosat_init + res <- f pico + picosat_reset pico + return res + unknown, satisfiable, unsatisfiable :: CInt unknown = 0 satisfiable = 10 @@ -86,15 +158,18 @@ data Solution = Solution [Int] | Unsatisfiable - | Unknown deriving (Show, Eq) + | Unknown deriving (Show, Eq, Ord) -addClause :: Ptr a -> [CInt] -> IO () -addClause pico cl = mapM_ (picosat_add pico) (cl ++ [0]) +addClause :: Picosat -> [Int] -> IO () +addClause pico cl = do + _ <- mapM_ (picosat_add pico . fromIntegral) cl + _ <- picosat_add pico 0 + return () -addClauses :: Ptr a -> [[CInt]] -> IO () -addClauses pico = mapM_ (addClause pico) +addClauses :: Picosat -> [[Int]] -> IO () +addClauses pico = mapM_ $ addClause pico -getSolution :: Ptr a -> IO Solution +getSolution :: Picosat -> IO Solution getSolution pico = do vars <- picosat_variables pico sol <- forM [1..vars] $ \i -> do @@ -102,7 +177,7 @@ return $ i * s return $ Solution $ map fromIntegral sol -solution :: Ptr a -> IO Solution +solution :: Picosat -> IO Solution solution pico = do res <- picosat_sat pico (-1) case res of @@ -111,27 +186,93 @@ | a == satisfiable -> getSolution pico | otherwise -> error "Picosat error." -toCInts :: [[Int]] -> [[CInt]] -toCInts = map $ map fromIntegral - -- | Solve a list of CNF constraints yielding the first solution. solve :: [[Int]] -> IO Solution -solve cls = do - let ccls = toCInts cls - pico <- picosat_init - _ <- addClauses pico ccls - sol <- solution pico - picosat_reset pico - return sol +solve cnf = do + withPicosat $ \ pico -> do + _ <- addClauses pico cnf + sol <- solution pico + return sol -- | Solve a list of CNF constraints yielding all possible solutions. solveAll :: [[Int]] -> IO [Solution] -solveAll e = do - let e' = map (map fromIntegral) e - s <- solve e' - case s of - Solution x -> (Solution x :) `fmap` solveAll (map negate x : e') - _ -> return [] +solveAll cnf = do + evalScopedPicosat $ do + addBaseClauses cnf + scopedAllSolutions + + +data PicosatScoped = PicosatScoped { psPicosat :: Picosat, + psContextVars :: S.Set Int } + +type PS a = StateT PicosatScoped IO a + +evalScopedPicosat :: PS a -> IO a +evalScopedPicosat action = + withPicosat $ \ picosat -> do + evalStateT action $ PicosatScoped picosat S.empty + +addBaseClauses :: [[Int]] -> PS () +addBaseClauses clauses = do + pico <- gets psPicosat + liftIO $ addClauses pico clauses + +withScopedClauses :: [[Int]] -> PS a -> PS a +withScopedClauses clauses action = do + pico <- gets psPicosat + withScope $ do + liftIO $ addClauses pico clauses + action + +withScope :: PS a -> PS a +withScope action = do + pico <- gets psPicosat + ctx <- liftIO $ picosat_push pico + addContextVariable $ fromIntegral ctx + res <- action + _ <- liftIO $ picosat_pop pico + return res + +addContextVariable :: Int -> PS () +addContextVariable var = modify add + where add s = s { psContextVars = S.insert var $ psContextVars s} + +-- | Get one solution in scoped context. Pay attention to not +-- return any "context variable" which are Picosat internals. +scopedSolution :: PS Solution +scopedSolution = do + pico <- gets psPicosat + sol <- liftIO $ solution pico + case sol of + Solution ys -> do + ctxvars <- gets psContextVars + return $ Solution $ + filter (\l -> S.notMember (abs l) ctxvars) $ ys + x -> + return x + + +scopedAllSolutions :: PS [Solution] +scopedAllSolutions = do + let recur solutions = do + pico <- gets psPicosat + sol <- scopedSolution + case sol of + Solution ys -> do + let negsol = map negate ys + liftIO $ addClause pico negsol + recur (sol : solutions) + _ -> + return $ reverse solutions + withScope $ recur [] + + +scopedSolutionWithAssumptions :: [Int] -> PS Solution +scopedSolutionWithAssumptions assumptions = do + pico <- gets psPicosat + liftIO $ mapM_ (picosat_assume pico . fromIntegral) assumptions + scopedSolution + -- Unsafe solver functions are not guaranteed to be memory safe if the solver fails internally. diff --git a/test/Rand.hs b/test/Rand.hs new file mode 100644 --- /dev/null +++ b/test/Rand.hs @@ -0,0 +1,74 @@ +import Picosat +import RandomCNF (randomLiteral, randomCNF) +import System.CPUTime.Rdtsc +import Data.Word (Word64) +import Control.Monad.IO.Class (liftIO, MonadIO) + +-- | Wrap any 'IO' computation so that it returns execution +-- time in CPU cycles as well as the real value. +measureCycles :: (MonadIO m) => m a -> m (Word64, a) +measureCycles ioa = do + start <- liftIO $ rdtsc + a <- ioa + a `seq` return () + end <- liftIO $ rdtsc + return (end-start, a) + + +-- | Test program demonstrating effectiveness of using Picosat's +-- solve-with-assumptions interface. One Picosat instance is kept +-- across many calls. It is able to use the kept state for good +-- improvements in run time. +-- +-- Random CNF's are generated and then solved multiple times, +-- each time assuming another additional random literal. +-- On my computer the resulting times look like this: +-- @ +-- -- 100 variables +-- ("num clauses",408) +-- ("unshared times:",[517,662,389,390,644,534,710,588,1808, +-- 486,587,526,937,692,750,882,671,545,649,531,445]) +-- ("shared times:",[460,185,25,23,21,20,19,20,22,153,203,25, +-- 30,23,41,24,23,22,37,25,22]) +-- -- ... +-- ("num clauses",418) +-- ("unshared times:",[1135,1384,913,1646,1753,2276,1277,1744, +-- 1385,1552,1725,1909,1783,1463,715,1561,1802,1816,1660,1970,2145]) +-- ("shared times:",[997,642,231,428,283,154,65,52,1,33,0,0,0, +-- 0,0,0,0,0,0,0,0]) +-- @ +-- Unshared times are close to constant. This is not surprising. Each +-- time the same CNF plus an additional literal is solved again and +-- again from cold. The shared times show how good Picosat runtimes +-- benefit from keeping one Picosat instance in memory. + +testNumSolutions num_vars negp clause_size num_rands num_clauses = + do cnf <- randomCNF num_vars negp clause_size num_clauses + someLiterals <- + mapM (\_->randomLiteral num_vars negp) [0..num_rands] + print ("num clauses", num_clauses) + let solveUnsharedWith r = do + solution <- solve $ cnf ++ [[r]] + case solution of + Unsatisfiable -> return 0 + Solution _ -> return 1 + re <- mapM (\r-> do + (c, _) <- measureCycles (solveUnsharedWith r) + return $ c `div` 3000) + someLiterals + print ("unshared times:", re) + + let solveSharedWith lit = do + (cyc, _) <- measureCycles $ + scopedSolutionWithAssumptions [lit] + return $ cyc `div` 3000 + ts <- evalScopedPicosat $ do + addBaseClauses cnf + mapM solveSharedWith someLiterals + print ("shared times:", ts) + + +main = do + mapM_ (testNumSolutions 100 0.5 3 20) [300,302..500] + + diff --git a/test/Scoped.hs b/test/Scoped.hs new file mode 100644 --- /dev/null +++ b/test/Scoped.hs @@ -0,0 +1,76 @@ +import Picosat +import Control.Monad.IO.Class (liftIO) +import Control.Monad (when) +import qualified Data.Set as S + +printAllSolutions = do + xs <- scopedAllSolutions + liftIO $ mapM_ print xs + +expectSolutions expected = do + real <- scopedAllSolutions + when (S.fromList expected /= S.fromList real) $ liftIO $ do + print ("expected", expected) + print ("real", real) + error "test failed" + +main = + evalScopedPicosat $ do + addBaseClauses [[1, 2, 3]] + + liftIO $ putStrLn "base cnf [[1, 2, 3]]" + printAllSolutions + expectSolutions [Solution [1,2,3], + Solution [1,2,-3], + Solution [1,-2,3], + Solution [1,-2,-3], + Solution [-1,-2,3], + Solution [-1,2,-3], + Solution [-1,2,3]] + + withScopedClauses [[-2,-3]] $ do + liftIO $ putStrLn "\nwith [-2,-3]" + printAllSolutions + expectSolutions [Solution [-1,2,-3], + Solution [-1,-2,3], + Solution [1,-2,-3], + Solution [1,-2,3], + Solution [1,2,-3]] + + withScopedClauses [[-1,-2]] $ do + liftIO $ putStrLn "\nwith [-1,-2]" + printAllSolutions + expectSolutions [Solution [1,-2,3], + Solution [1,-2,-3], + Solution [-1,-2,3], + Solution [-1,2,-3], + Solution [-1,2,3]] + + addBaseClauses [[-1,-3]] + expectSolutions [Solution [-1,2,-3], + Solution [-1,2,3], + Solution [-1,-2,3], + Solution [1,-2,-3], + Solution [1,2,-3]] + + withScopedClauses [[-2,-3]] $ do + expectSolutions [Solution [1,-2,-3], + Solution [1,2,-3], + Solution [-1,-2,3], + Solution [-1,2,-3]] + + withScopedClauses [[-1,-2], [1,-3]] $ do + expectSolutions [Solution [-1,2,-3], + Solution [1,-2,-3]] + + let printSolutionsWithAssumptions as = do + liftIO $ putStrLn ("\nwith assumptions " ++ show as) + res <- scopedSolutionWithAssumptions as + liftIO $ print res + + printSolutionsWithAssumptions [1] + printSolutionsWithAssumptions [-1] + printSolutionsWithAssumptions [2] + printSolutionsWithAssumptions [-2] + printSolutionsWithAssumptions [3] + printSolutionsWithAssumptions [-3] diff --git a/test/Sudoku.hs b/test/Sudoku.hs new file mode 100644 --- /dev/null +++ b/test/Sudoku.hs @@ -0,0 +1,77 @@ +import Picosat ( solve, Solution(..) ) + +cross :: [a] -> [(a, a)] +cross list = [(x,y) | x <- list, y <- list] + +rows, cols :: [[(Int, Int)]] +rows = [[(i,j) | i <- [0..8]] | j <- [0..8]] +cols = [[(i,j) | j <- [0..8]] | i <- [0..8]] + +c2v :: (Int, Int, Int) -> Int +c2v (i,j,k) = (k-1)*81 + i*9 + j + 1 + +v2c :: Int -> (Int, Int, Int) +v2c x = ((i `mod` 81) `div` 9, i `mod` 9, (i `div` 81)+1) + where i = x - 1 + +group :: [(Int, Int)] -> [[Int]] +group grp = foldr ((:).label) [] [1..9] + where label k = map c2v [(i,j,k) | (i,j) <- grp ] + +square :: [[(Int, Int)]] +square = [quadrent i j | (i,j) <- cross [0..2]] + where quadrent x y = [(x*3+i,y*3+j) | (i,j) <- cross [0..2]] + +oneLabel :: (Int, Int) -> [[Int]] +oneLabel (i,j) = atLeastOne : lessThan2 + where notBoth (c1,c2) = [- c2v (i,j,c1), - c2v (i,j,c2)] + lessThan2 = map notBoth $ [(a,b) | (a,b) <- cross [1..9], a /= b] + atLeastOne = map c2v [(i,j,k) | k <- [1..9]] + +validLabeling :: [[Int]] +validLabeling = foldr ((++) . oneLabel) [] (cross [0..8]) + +goodLabeling :: [[Int]] +goodLabeling = foldr ((++) . group) [] (square ++ rows ++ cols) + +sudokuForm :: [(Int, Int, Int)] -> [[Int]] +sudokuForm cells = validLabeling ++ goodLabeling ++ (map consClause cells) + where consClause cell = [c2v cell] + +getConstraints :: [[Int]] -> [(Int, Int, Int)] +getConstraints matrix = filter (\(_,_,a) -> a > 0) cells + where flat = foldl1 (++) matrix + cells = zip3 [i `div` 9 | i <- [0..]] (cycle [0..8]) flat + +toSudoku :: [Int] -> [[Int]] +toSudoku [] = [] +toSudoku sol = [[lkup i j | j <- [0..8]] + | i <- [0..8]] + where + cells = map v2c $ filter ((<) 0) sol + lkup i j = label $ head $ filter (\(a,b,_) -> a == i && b == j) cells + label (_,_,c) = c + +stringToMatrix :: String -> [[Int]] +stringToMatrix = map (map $ read . return) . lines + +puzzle :: String +puzzle = unlines [ + "004000000" + , "000030002" + , "390700080" + , "400009001" + , "209801307" + , "600200008" + , "010008053" + , "900040000" + , "000000800" + ] + +main :: IO () +main = do + let cnf = sudokuForm $ getConstraints $ stringToMatrix puzzle + sol <- solve cnf + case sol of + Solution s -> mapM_ print $ toSudoku s + _ -> putStrLn "Puzzle not solvable."